|
|
1.1 root 1: subroutine dglsbp(nu, order, bc, e)
2: integer nu
3: integer order(nu, nu, 2), bc(nu, 2, 2), e(nu, 2, 2)
4: common /cstak/ ds
5: double precision ds(500)
6: integer ice, istkgt, max0, i, j, l
7: integer ipps, inow, imaord, inoold, is(1000)
8: real rs(1000)
9: logical allero, ls(1000)
10: complex cs(500)
11: double precision ws(500)
12: integer temp, temp1
13: equivalence (ds(1), cs(1), ws(1), rs(1), is(1), ls(1))
14: c To determine which ODE should use which boundary condition.
15: c Mnemonic - Double precision Galerkin's method for Linear Systems,
16: c Boundary condition Placement.
17: c Scratch Space Allocated -
18: c S(DglsBP) <= Nu*(4*Nu+15)
19: c Integer words.
20: c Define Node = Is(inow) -> Nodei
21: c Check the input for errors.
22: if (nu .lt. 1) call seterr(16hDglsBP - Nu.lt.1, 16, 1, 2)
23: do 3 l = 1, 2
24: do 2 i = 1, nu
25: c Is Order(i,.,l) = (-1,...,-1)?
26: allero = .true.
27: do 1 j = 1, nu
28: allero = allero .and. order(i, j, l) .eq. (-1)
29: if (order(i, j, l) .lt. (-1) .or. order(i, j, l) .gt. 2)
30: 1 call seterr(
31: 2 41hDglsBP - Order(i,j,l) not one of -1,0,1,2, 41, 2, 2
32: 3 )
33: 1 continue
34: if (bc(i, 1, l) .ne. (-2) .and. bc(i, 1, l) .ne. 0 .or. bc(i
35: 1 , 2, l) .ne. (-2) .and. bc(i, 2, l) .ne. 1) call seterr(
36: 2 36hDglsBP - BC(i,.,l) not one of -2,0,1, 36, 3, 2)
37: if (allero) call seterr(
38: 1 33hDglsBP - Order(i,.,l)=(-1,...,-1), 33, 4, 2)
39: 2 continue
40: 3 continue
41: call enter(1)
42: c Complement of E.
43: ice = istkgt(nu, 2)
44: c Maxord(i,l) = Max over j=1,...,Nu Order(i,j,l).
45: imaord = istkgt(2*nu, 2)
46: call seti(2*nu, -1, is(imaord))
47: do 6 l = 1, 2
48: do 5 i = 1, nu
49: do 4 j = 1, nu
50: temp1 = imaord+i-1+(l-1)*nu
51: temp = imaord+i-1+(l-1)*nu
52: is(temp1) = max0(is(temp), order(i, j, l))
53: 4 continue
54: 5 continue
55: 6 continue
56: i = 0
57: ipps = 1
58: 7 if (i .ge. 4*nu .and. ipps .eq. 1) goto 15
59: goto 12
60: c Make a node.
61: 8 inoold = inow
62: i = i+1
63: inow = istkgt(nu+3, 2)
64: is(inow) = inoold
65: c Get the candidates for E(i).
66: call d6lsbp(i, nu, order, bc, e, is(imaord), is(ice), is(
67: 1 inow+3), is(inow+1))
68: is(inow+2) = 0
69: ipps = 0
70: goto 13
71: goto 13
72: c Searching a node.
73: 9 is(inow+2) = is(inow+2)+1
74: if (is(inow+2) .le. is(inow+1)) goto 10
75: ipps = -1
76: c Back-up.
77: goto 7
78: 10 temp = inow+2+is(inow+2)-1
79: e(i, 1, 1) = is(temp+1)
80: ipps = 1
81: goto 13
82: goto 13
83: c Backing up a Node.
84: 11 inow = is(inow)
85: call istkrl(1)
86: i = i-1
87: ipps = 0
88: goto 13
89: goto 13
90: 12 temp = ipps+2
91: if (temp .gt. 0 .and. temp .le. 3) goto ( 11, 9, 8), temp
92: c End Switch.
93: 13 if (i .ne. 0) goto 14
94: call seterr(37hDglsBP - Improper Boundary Conditions, 37, 5,
95: 1 1)
96: goto 15
97: 14 continue
98: goto 7
99: c End While.
100: 15 call leave
101: return
102: end
103: subroutine d6lsbp(i, nu, order, bc, e, maxord, ce, r, n)
104: integer nu
105: integer i, order(nu, nu, 2), bc(1), e(1), maxord(nu, 2), ce(nu)
106: integer r(nu), n
107: integer mod, max0, j, l, nbcs, dm
108: integer ii, lr
109: integer temp
110: c BC(Nu,2,2),E(Nu,2,2),
111: c E(i-1),R(N).
112: if (bc(i) .ge. 0) goto 1
113: n = 1
114: r(n) = 0
115: return
116: c LR = 1 for left, LR = 2 for right.
117: 1 lr = (i-1)/(2*nu)+1
118: c DM = 1 for Dirichlet, DM = 2 for Mixed boundary conditions.
119: dm = mod((i-1)/nu, 2)+1
120: ii = mod(i, nu)
121: if (ii .eq. 0) ii = nu
122: c B(i) = B(ii,DM,LR).
123: n = 0
124: do 2 j = 1, nu
125: ce(j) = j
126: 2 continue
127: c CE = Complement of E.
128: if (i .gt. 2*nu) goto 7
129: j = 1
130: goto 4
131: 3 j = j+1
132: 4 if (j .ge. i) goto 6
133: if (bc(j) .lt. 0) goto 5
134: temp = e(j)
135: ce(temp) = 0
136: 5 continue
137: goto 3
138: 6 continue
139: goto 12
140: 7 j = 2*nu+1
141: goto 9
142: 8 j = j+1
143: 9 if (j .ge. i) goto 11
144: if (bc(j) .lt. 0) goto 10
145: temp = e(j)
146: ce(temp) = 0
147: 10 continue
148: goto 8
149: 11 continue
150: 12 do 18 j = 1, nu
151: if (ce(j) .eq. 0) goto 18
152: nbcs = 0
153: l = 1
154: goto 14
155: 13 l = l+1
156: 14 if (l .ge. i) goto 15
157: if (e(l) .eq. j .and. bc(l) .ge. 0) nbcs = nbcs+1
158: goto 13
159: 15 if ((dm .ne. 1 .or. maxord(j, lr) .le. bc(i)) .and. (dm .ne. 2
160: 1 .or. order(j, ii, lr) .le. bc(i))) goto 17
161: if (nbcs .ge. max0(maxord(j, 1), maxord(j, 2))) goto 16
162: n = n+1
163: r(n) = j
164: 16 continue
165: 17 continue
166: 18 continue
167: return
168: end
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.