/[MITgcm]/MITgcm/model/src/impldiff.F
ViewVC logotype

Contents of /MITgcm/model/src/impldiff.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.16 - (show annotations) (download)
Mon May 14 21:46:17 2001 UTC (23 years, 1 month ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint40pre3, checkpoint40pre1, checkpoint40pre7, checkpoint40pre6, checkpoint40pre9, checkpoint40pre8, checkpoint40pre2, checkpoint40pre4, checkpoint39, checkpoint40pre5, checkpoint40
Changes since 1.15: +9 -2 lines
Modifications/fixes to support TAMC differentiability
(mostly missing or wrong directives).

1 C $Header: /u/gcmpack/models/MITgcmUV/model/src/impldiff.F,v 1.15 2001/02/04 14:38:47 cnh Exp $
2 C $Name: checkpoint38 $
3
4 #include "CPP_OPTIONS.h"
5
6 C |==========================================================|
7 C | S/R IMPLDIFF |
8 C | o Solve implicit diffusion equation for vertical |
9 C | diffusivity. |
10 C | o Recoded from 2d intermediate fields to 3d to reduce |
11 C | TAMC storage |
12 C | o Fixed missing masks for fields a(), c() |
13 C |==========================================================|
14 SUBROUTINE IMPLDIFF( bi, bj, iMin, iMax, jMin, jMax,
15 I deltaTX,KappaRX,recip_hFac,
16 U gXNm1,
17 I myThid )
18 IMPLICIT NONE
19 C == Global data ==
20 #include "SIZE.h"
21 #include "DYNVARS.h"
22 #include "EEPARAMS.h"
23 #include "PARAMS.h"
24 #include "GRID.h"
25
26 #ifdef ALLOW_AUTODIFF_TAMC
27 #include "tamc_keys.h"
28 #endif
29
30 C == Routine Arguments ==
31 INTEGER bi,bj,iMin,iMax,jMin,jMax
32 _RL deltaTX
33 _RL KappaRX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
34 _RS recip_hFac(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
35 _RL gXnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
36 _RL gYnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
37 INTEGER myThid
38
39 C == Local variables ==
40 INTEGER i,j,k
41 _RL a(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
42 _RL b(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
43 _RL c(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
44 _RL bet(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
45 _RL gam(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
46
47 C-- Initialise
48 DO k=1,Nr
49 DO j=jMin,jMax
50 DO i=iMin,iMax
51 gYNm1(i,j,k,bi,bj) = 0. _d 0
52 ENDDO
53 ENDDO
54 ENDDO
55
56 C-- Old aLower
57 DO j=jMin,jMax
58 DO i=iMin,iMax
59 a(i,j,1) = 0. _d 0
60 ENDDO
61 ENDDO
62 DO k=2,Nr
63 DO j=jMin,jMax
64 DO i=iMin,iMax
65 a(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
66 & *KappaRX(i,j, k )*recip_drC( k )
67 ENDDO
68 ENDDO
69 ENDDO
70
71 C-- Old aUpper
72 DO k=1,Nr-1
73 DO j=jMin,jMax
74 DO i=iMin,iMax
75 c(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
76 & *KappaRX(i,j,k+1)*recip_drC(k+1)
77 IF (recip_hFac(i,j,k+1,bi,bj).EQ.0.) c(i,j,k)=0.
78 ENDDO
79 ENDDO
80 ENDDO
81 DO j=jMin,jMax
82 DO i=iMin,iMax
83 c(i,j,Nr) = 0. _d 0
84 ENDDO
85 ENDDO
86
87 C-- Old aCenter
88 DO k=1,Nr
89 DO j=jMin,jMax
90 DO i=iMin,iMax
91 b(i,j,k) = 1. _d 0 - c(i,j,k) - a(i,j,k)
92 ENDDO
93 ENDDO
94 ENDDO
95
96 C-- Old and new gam, bet are the same
97 DO k=1,Nr
98 DO j=jMin,jMax
99 DO i=iMin,iMax
100 bet(i,j,k) = 0. _d 0
101 gam(i,j,k) = 0. _d 0
102 ENDDO
103 ENDDO
104 ENDDO
105
106 C-- Only need do anything if Nr>1
107 IF (Nr.GT.1) THEN
108
109 k = 1
110 C-- Beginning of forward sweep (top level)
111 DO j=jMin,jMax
112 DO i=iMin,iMax
113 IF (b(i,j,1).NE.0.) bet(i,j,1) = 1. _d 0 / b(i,j,1)
114 ENDDO
115 ENDDO
116
117 ENDIF
118
119 C-- Middle of forward sweep
120 IF (Nr.GT.2) THEN
121
122 CADJ loop = sequential
123 DO k=2,Nr
124
125 DO j=jMin,jMax
126 DO i=iMin,iMax
127 gam(i,j,k) = c(i,j,k-1)*bet(i,j,k-1)
128 IF ( ( b(i,j,k) - a(i,j,k)*gam(i,j,k) ) .NE. 0.)
129 & bet(i,j,k) = 1. _d 0 / ( b(i,j,k) - a(i,j,k)*gam(i,j,k) )
130 ENDDO
131 ENDDO
132
133 ENDDO
134
135 ENDIF
136
137
138 DO j=jMin,jMax
139 DO i=iMin,iMax
140 gYNm1(i,j,1,bi,bj) = gXNm1(i,j,1,bi,bj)*bet(i,j,1)
141 ENDDO
142 ENDDO
143 DO k=2,Nr
144 DO j=jMin,jMax
145 DO i=iMin,iMax
146 gYnm1(i,j,k,bi,bj) = bet(i,j,k)*
147 & (gXnm1(i,j,k,bi,bj) - a(i,j,k)*gYnm1(i,j,k-1,bi,bj))
148 ENDDO
149 ENDDO
150 ENDDO
151
152
153 C-- Backward sweep
154 CADJ loop = sequential
155 DO k=Nr-1,1,-1
156 DO j=jMin,jMax
157 DO i=iMin,iMax
158 gYnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
159 & -gam(i,j,k+1)*gYnm1(i,j,k+1,bi,bj)
160 ENDDO
161 ENDDO
162 ENDDO
163
164 DO k=1,Nr
165 DO j=jMin,jMax
166 DO i=iMin,iMax
167 gXnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
168 ENDDO
169 ENDDO
170 ENDDO
171
172 RETURN
173 END

  ViewVC Help
Powered by ViewVC 1.1.22