1 |
SUBROUTINE qckxyz(q, qtmp, imr, jnp, nlay, j1, j2, cosp, acosp, cross, ic, & |
2 |
nstep) |
3 |
|
4 |
PARAMETER (tiny=1.E-60) |
5 |
DIMENSION q(imr, jnp, nlay), qtmp(imr, jnp), cosp(*), acosp(*) |
6 |
LOGICAL cross |
7 |
|
8 |
nlaym1 = nlay - 1 |
9 |
len = imr*(j2-j1+1) |
10 |
ip = 0 |
11 |
|
12 |
! Top layer |
13 |
l = 1 |
14 |
icr = 1 |
15 |
CALL filns(q(1,1,l), imr, jnp, j1, j2, cosp, acosp, ipy, tiny) |
16 |
IF (ipy==0) GO TO 50 |
17 |
CALL filew(q(1,1,l), qtmp, imr, jnp, j1, j2, ipx, tiny) |
18 |
IF (ipx==0) GO TO 50 |
19 |
|
20 |
IF (cross) THEN |
21 |
CALL filcr(q(1,1,l), imr, jnp, j1, j2, cosp, acosp, icr, tiny) |
22 |
END IF |
23 |
IF (icr==0) GO TO 50 |
24 |
|
25 |
! Vertical filling... |
26 |
DO i = 1, len |
27 |
IF (q(i,j1,1)<0.) THEN |
28 |
ip = ip + 1 |
29 |
q(i, j1, 2) = q(i, j1, 2) + q(i, j1, 1) |
30 |
q(i, j1, 1) = 0. |
31 |
END IF |
32 |
END DO |
33 |
|
34 |
50 CONTINUE |
35 |
DO l = 2, nlaym1 |
36 |
icr = 1 |
37 |
|
38 |
CALL filns(q(1,1,l), imr, jnp, j1, j2, cosp, acosp, ipy, tiny) |
39 |
IF (ipy==0) GO TO 225 |
40 |
CALL filew(q(1,1,l), qtmp, imr, jnp, j1, j2, ipx, tiny) |
41 |
IF (ipx==0) GO TO 225 |
42 |
IF (cross) THEN |
43 |
CALL filcr(q(1,1,l), imr, jnp, j1, j2, cosp, acosp, icr, tiny) |
44 |
END IF |
45 |
IF (icr==0) GO TO 225 |
46 |
|
47 |
DO i = 1, len |
48 |
IF (q(i,j1,l)<0.) THEN |
49 |
|
50 |
ip = ip + 1 |
51 |
! From above |
52 |
qup = q(i, j1, l-1) |
53 |
qly = -q(i, j1, l) |
54 |
dup = min(qly, qup) |
55 |
q(i, j1, l-1) = qup - dup |
56 |
q(i, j1, l) = dup - qly |
57 |
! Below |
58 |
q(i, j1, l+1) = q(i, j1, l+1) + q(i, j1, l) |
59 |
q(i, j1, l) = 0. |
60 |
END IF |
61 |
END DO |
62 |
225 END DO |
63 |
|
64 |
! BOTTOM LAYER |
65 |
sum = 0. |
66 |
l = nlay |
67 |
|
68 |
CALL filns(q(1,1,l), imr, jnp, j1, j2, cosp, acosp, ipy, tiny) |
69 |
IF (ipy==0) GO TO 911 |
70 |
CALL filew(q(1,1,l), qtmp, imr, jnp, j1, j2, ipx, tiny) |
71 |
IF (ipx==0) GO TO 911 |
72 |
|
73 |
CALL filcr(q(1,1,l), imr, jnp, j1, j2, cosp, acosp, icr, tiny) |
74 |
IF (icr==0) GO TO 911 |
75 |
|
76 |
DO i = 1, len |
77 |
IF (q(i,j1,l)<0.) THEN |
78 |
ip = ip + 1 |
79 |
|
80 |
! From above |
81 |
|
82 |
qup = q(i, j1, nlaym1) |
83 |
qly = -q(i, j1, l) |
84 |
dup = min(qly, qup) |
85 |
q(i, j1, nlaym1) = qup - dup |
86 |
! From "below" the surface. |
87 |
sum = sum + qly - dup |
88 |
q(i, j1, l) = 0. |
89 |
END IF |
90 |
END DO |
91 |
|
92 |
911 CONTINUE |
93 |
|
94 |
IF (ip>imr) THEN |
95 |
WRITE (6, *) 'IC=', ic, ' STEP=', nstep, ' Vertical filling pts=', ip |
96 |
END IF |
97 |
|
98 |
IF (sum>1.E-25) THEN |
99 |
WRITE (6, *) ic, nstep, ' Mass source from the ground=', sum |
100 |
END IF |
101 |
RETURN |
102 |
END SUBROUTINE qckxyz |
103 |
|