|
|
BSD 4.0
subroutine dglsbp(nu, order, bc, e)
integer nu
integer order(nu, nu, 2), bc(nu, 2, 2), e(nu, 2, 2)
common /cstak/ ds
double precision ds(500)
integer ice, istkgt, max0, i, j, l
integer ipps, inow, imaord, inoold, is(1000)
real rs(1000)
logical allero, ls(1000)
complex cs(500)
double precision ws(500)
integer temp, temp1
equivalence (ds(1), cs(1), ws(1), rs(1), is(1), ls(1))
c To determine which ODE should use which boundary condition.
c Mnemonic - Double precision Galerkin's method for Linear Systems,
c Boundary condition Placement.
c Scratch Space Allocated -
c S(DglsBP) <= Nu*(4*Nu+15)
c Integer words.
c Define Node = Is(inow) -> Nodei
c Check the input for errors.
if (nu .lt. 1) call seterr(16hDglsBP - Nu.lt.1, 16, 1, 2)
do 3 l = 1, 2
do 2 i = 1, nu
c Is Order(i,.,l) = (-1,...,-1)?
allero = .true.
do 1 j = 1, nu
allero = allero .and. order(i, j, l) .eq. (-1)
if (order(i, j, l) .lt. (-1) .or. order(i, j, l) .gt. 2)
1 call seterr(
2 41hDglsBP - Order(i,j,l) not one of -1,0,1,2, 41, 2, 2
3 )
1 continue
if (bc(i, 1, l) .ne. (-2) .and. bc(i, 1, l) .ne. 0 .or. bc(i
1 , 2, l) .ne. (-2) .and. bc(i, 2, l) .ne. 1) call seterr(
2 36hDglsBP - BC(i,.,l) not one of -2,0,1, 36, 3, 2)
if (allero) call seterr(
1 33hDglsBP - Order(i,.,l)=(-1,...,-1), 33, 4, 2)
2 continue
3 continue
call enter(1)
c Complement of E.
ice = istkgt(nu, 2)
c Maxord(i,l) = Max over j=1,...,Nu Order(i,j,l).
imaord = istkgt(2*nu, 2)
call seti(2*nu, -1, is(imaord))
do 6 l = 1, 2
do 5 i = 1, nu
do 4 j = 1, nu
temp1 = imaord+i-1+(l-1)*nu
temp = imaord+i-1+(l-1)*nu
is(temp1) = max0(is(temp), order(i, j, l))
4 continue
5 continue
6 continue
i = 0
ipps = 1
7 if (i .ge. 4*nu .and. ipps .eq. 1) goto 15
goto 12
c Make a node.
8 inoold = inow
i = i+1
inow = istkgt(nu+3, 2)
is(inow) = inoold
c Get the candidates for E(i).
call d6lsbp(i, nu, order, bc, e, is(imaord), is(ice), is(
1 inow+3), is(inow+1))
is(inow+2) = 0
ipps = 0
goto 13
goto 13
c Searching a node.
9 is(inow+2) = is(inow+2)+1
if (is(inow+2) .le. is(inow+1)) goto 10
ipps = -1
c Back-up.
goto 7
10 temp = inow+2+is(inow+2)-1
e(i, 1, 1) = is(temp+1)
ipps = 1
goto 13
goto 13
c Backing up a Node.
11 inow = is(inow)
call istkrl(1)
i = i-1
ipps = 0
goto 13
goto 13
12 temp = ipps+2
if (temp .gt. 0 .and. temp .le. 3) goto ( 11, 9, 8), temp
c End Switch.
13 if (i .ne. 0) goto 14
call seterr(37hDglsBP - Improper Boundary Conditions, 37, 5,
1 1)
goto 15
14 continue
goto 7
c End While.
15 call leave
return
end
subroutine d6lsbp(i, nu, order, bc, e, maxord, ce, r, n)
integer nu
integer i, order(nu, nu, 2), bc(1), e(1), maxord(nu, 2), ce(nu)
integer r(nu), n
integer mod, max0, j, l, nbcs, dm
integer ii, lr
integer temp
c BC(Nu,2,2),E(Nu,2,2),
c E(i-1),R(N).
if (bc(i) .ge. 0) goto 1
n = 1
r(n) = 0
return
c LR = 1 for left, LR = 2 for right.
1 lr = (i-1)/(2*nu)+1
c DM = 1 for Dirichlet, DM = 2 for Mixed boundary conditions.
dm = mod((i-1)/nu, 2)+1
ii = mod(i, nu)
if (ii .eq. 0) ii = nu
c B(i) = B(ii,DM,LR).
n = 0
do 2 j = 1, nu
ce(j) = j
2 continue
c CE = Complement of E.
if (i .gt. 2*nu) goto 7
j = 1
goto 4
3 j = j+1
4 if (j .ge. i) goto 6
if (bc(j) .lt. 0) goto 5
temp = e(j)
ce(temp) = 0
5 continue
goto 3
6 continue
goto 12
7 j = 2*nu+1
goto 9
8 j = j+1
9 if (j .ge. i) goto 11
if (bc(j) .lt. 0) goto 10
temp = e(j)
ce(temp) = 0
10 continue
goto 8
11 continue
12 do 18 j = 1, nu
if (ce(j) .eq. 0) goto 18
nbcs = 0
l = 1
goto 14
13 l = l+1
14 if (l .ge. i) goto 15
if (e(l) .eq. j .and. bc(l) .ge. 0) nbcs = nbcs+1
goto 13
15 if ((dm .ne. 1 .or. maxord(j, lr) .le. bc(i)) .and. (dm .ne. 2
1 .or. order(j, ii, lr) .le. bc(i))) goto 17
if (nbcs .ge. max0(maxord(j, 1), maxord(j, 2))) goto 16
n = n+1
r(n) = j
16 continue
17 continue
18 continue
return
end
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.