1 |
dgoldberg |
1.1 |
#include "STREAMICE_OPTIONS.h" |
2 |
|
|
|
3 |
|
|
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| |
4 |
|
|
|
5 |
|
|
CBOP |
6 |
|
|
|
7 |
|
|
SUBROUTINE STREAMICE_TRIDIAG_SOLVE( |
8 |
|
|
U cg_Uin, ! x-velocities |
9 |
|
|
U cg_Vin, ! x-velocities |
10 |
|
|
U cg_Bu, ! force in x dir |
11 |
|
|
I A_uu, ! section of matrix that multiplies u and projects on u |
12 |
|
|
I umask, |
13 |
|
|
I myThid ) |
14 |
|
|
C /============================================================\ |
15 |
|
|
C | SUBROUTINE | |
16 |
|
|
C | o | |
17 |
|
|
C |============================================================| |
18 |
|
|
C | | |
19 |
|
|
C \============================================================/ |
20 |
|
|
IMPLICIT NONE |
21 |
|
|
|
22 |
|
|
#include "SIZE.h" |
23 |
|
|
#include "EEPARAMS.h" |
24 |
|
|
#include "PARAMS.h" |
25 |
|
|
#include "STREAMICE.h" |
26 |
|
|
#include "STREAMICE_CG.h" |
27 |
|
|
|
28 |
|
|
_RL cg_Uin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
29 |
|
|
_RL cg_Vin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
30 |
|
|
_RL cg_Bu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
31 |
|
|
_RL A_uu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1) |
32 |
|
|
_RS umask (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
33 |
|
|
INTEGER myThid |
34 |
|
|
|
35 |
dgoldberg |
1.2 |
INTEGER iMin,iMax,i,j,k |
36 |
dgoldberg |
1.1 |
_RL aMat(1:Nx) |
37 |
|
|
_RL bMat(1:Nx) |
38 |
|
|
_RL cMat(1:Nx) |
39 |
|
|
_RL yMat(1:Nx) |
40 |
|
|
_RL bet(1:Nx) |
41 |
|
|
_RL tmpvar |
42 |
|
|
INTEGER errCode |
43 |
|
|
|
44 |
dgoldberg |
1.2 |
|
45 |
|
|
! CALL WRITE_FLD_XY_RL ("taud_tri","",cg_Bu,0,mythid) |
46 |
|
|
! CALL WRITE_FLD_XY_RL ("A_m1m1","",A_uu(:,:,:,:,-1,-1),0,mythid) |
47 |
|
|
! CALL WRITE_FLD_XY_RL ("A_m1_0","",A_uu(:,:,:,:,-1,0),0,mythid) |
48 |
|
|
! CALL WRITE_FLD_XY_RL ("A_m1p1","",A_uu(:,:,:,:,-1,1),0,mythid) |
49 |
|
|
! CALL WRITE_FLD_XY_RL ("A_0_m1","",A_uu(:,:,:,:,0,-1),0,mythid) |
50 |
|
|
! CALL WRITE_FLD_XY_RL ("A_0_0","",A_uu(:,:,:,:,0,0),0,mythid) |
51 |
|
|
! CALL WRITE_FLD_XY_RL ("A_0_p1","",A_uu(:,:,:,:,0,1),0,mythid) |
52 |
|
|
! CALL WRITE_FLD_XY_RL ("A_p1m1","",A_uu(:,:,:,:,1,-1),0,mythid) |
53 |
|
|
! CALL WRITE_FLD_XY_RL ("A_p1_0","",A_uu(:,:,:,:,1,0),0,mythid) |
54 |
|
|
! CALL WRITE_FLD_XY_RL ("A_p1p1","",A_uu(:,:,:,:,1,1),0,mythid) |
55 |
|
|
|
56 |
|
|
|
57 |
dgoldberg |
1.1 |
|
58 |
|
|
IF (nPx.gt.1 .or. nSx.gt.1) THEN |
59 |
|
|
STOP 'must be serial for tridiag solve' |
60 |
|
|
ENDIF |
61 |
|
|
|
62 |
|
|
errCode = 0 |
63 |
|
|
|
64 |
|
|
imax = 0 |
65 |
|
|
iMin = 2 |
66 |
|
|
do i=imin,Nx |
67 |
|
|
if (umask(i,1,1,1).eq.1.0) THEN |
68 |
|
|
|
69 |
|
|
aMat(i)=0.0 |
70 |
|
|
bmat(i)=0.0 |
71 |
|
|
cmat(i)=0.0 |
72 |
|
|
ymat(i)=0.0 |
73 |
|
|
do j=-1,1 |
74 |
dgoldberg |
1.2 |
do k=1,3 |
75 |
|
|
aMat(i) = amat(i)+A_uu(i,k,1,1,-1,j) |
76 |
|
|
bMat(i) = bmat(i)+A_uu(i,k,1,1,0,j) |
77 |
|
|
cMat(i) = cmat(i)+A_uu(i,k,1,1,1,j) |
78 |
|
|
enddo |
79 |
|
|
yMat(i) = ymat(i)+cg_Bu(i,j+2,1,1) |
80 |
dgoldberg |
1.1 |
enddo |
81 |
|
|
else |
82 |
|
|
iMax = i-1 |
83 |
|
|
exit |
84 |
|
|
endif |
85 |
|
|
enddo |
86 |
|
|
|
87 |
|
|
IF(imax.eq.0) THEN |
88 |
|
|
imax=Nx |
89 |
|
|
ENDIF |
90 |
|
|
|
91 |
|
|
|
92 |
|
|
IF ( bMat(imin).NE.0. _d 0 ) THEN |
93 |
|
|
bet(imin) = 1. _d 0 / bMat(imin) |
94 |
|
|
ELSE |
95 |
|
|
bet(imin) = 0. _d 0 |
96 |
|
|
errCode = 1 |
97 |
|
|
ENDIF |
98 |
|
|
|
99 |
|
|
DO i=imin+1,imax |
100 |
|
|
tmpvar = bmat(i) - amat(i)*cmat(i-1)*bet(i-1) |
101 |
|
|
IF ( tmpvar .NE. 0. _d 0 ) THEN |
102 |
|
|
bet(i) = 1. _d 0 / tmpvar |
103 |
|
|
ELSE |
104 |
|
|
bet(i) = 0. _d 0 |
105 |
|
|
errCode = 1 |
106 |
|
|
ENDIF |
107 |
|
|
ENDDO |
108 |
|
|
|
109 |
|
|
|
110 |
|
|
ymat(imin) = ymat(imin)*bet(imin) |
111 |
|
|
|
112 |
|
|
DO i=imin+1,imax |
113 |
|
|
ymat(i) = ( ymat(i) |
114 |
|
|
& - amat(i)*ymat(i-1) |
115 |
|
|
& )*bet(i) |
116 |
|
|
ENDDO |
117 |
|
|
|
118 |
|
|
|
119 |
|
|
DO i=imax-1,imin,-1 |
120 |
|
|
ymat(i) = ymat(i) |
121 |
|
|
& - cmat(i)*bet(i)*ymat(i+1) |
122 |
|
|
ENDDO |
123 |
|
|
|
124 |
|
|
DO j=1,sNy |
125 |
|
|
DO i=imin,imax |
126 |
|
|
cg_Uin (i,j,1,1) = ymat(i) |
127 |
|
|
ENDDO |
128 |
|
|
ENDDO |
129 |
|
|
|
130 |
|
|
DO j=1,sNy |
131 |
|
|
DO i=1,sNx |
132 |
|
|
cg_Vin (i,j,1,1) = 0. _d 0 |
133 |
|
|
ENDDO |
134 |
|
|
ENDDO |
135 |
|
|
|
136 |
|
|
print *, "ERRORCODE", errcode |
137 |
|
|
|
138 |
|
|
RETURN |
139 |
|
|
END |