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

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

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

Revision 1.20 - (hide annotations) (download)
Fri Jan 2 23:07:14 2004 UTC (20 years, 5 months ago) by jmc
Branch: MAIN
CVS Tags: hrcube4, checkpoint52j_pre, checkpoint52f_post, checkpoint52i_pre, hrcube_1, hrcube_2, hrcube_3, checkpoint52e_pre, checkpoint52e_post, checkpoint52f_pre, checkpoint52i_post, checkpoint52h_pre
Changes since 1.19: +2 -2 lines 
fix bug when Nr=2
1 jmc 1.20 C $Header: /u/gcmpack/MITgcm/model/src/impldiff.F,v 1.19 2003/12/28 17:25:25 jmc Exp $
2 cnh 1.17 C $Name: $
3 adcroft 1.1
4 cnh 1.7 #include "CPP_OPTIONS.h"
5 adcroft 1.1
6 cnh 1.17 CBOP
7     C !ROUTINE: IMPLDIFF
8     C !INTERFACE:
9 adcroft 1.1 SUBROUTINE IMPLDIFF( bi, bj, iMin, iMax, jMin, jMax,
10 adcroft 1.8 I deltaTX,KappaRX,recip_hFac,
11     U gXNm1,
12 adcroft 1.1 I myThid )
13 cnh 1.17 C !DESCRIPTION: \bv
14     C *==========================================================*
15     C | S/R IMPLDIFF
16     C | o Solve implicit diffusion equation for vertical
17     C | diffusivity.
18     C *==========================================================*
19     C | o Recoded from 2d intermediate fields to 3d to reduce
20     C | TAMC storage
21     C | o Fixed missing masks for fields a(), c()
22     C *==========================================================*
23     C \ev
24    
25     C !USES:
26 cnh 1.5 IMPLICIT NONE
27     C == Global data ==
28 adcroft 1.1 #include "SIZE.h"
29     #include "DYNVARS.h"
30 cnh 1.2 #include "EEPARAMS.h"
31 adcroft 1.1 #include "PARAMS.h"
32     #include "GRID.h"
33 heimbach 1.9 #ifdef ALLOW_AUTODIFF_TAMC
34     #include "tamc_keys.h"
35     #endif
36    
37 cnh 1.17 C !INPUT/OUTPUT PARAMETERS:
38 adcroft 1.1 C == Routine Arguments ==
39     INTEGER bi,bj,iMin,iMax,jMin,jMax
40 adcroft 1.8 _RL deltaTX
41     _RL KappaRX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
42     _RS recip_hFac(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
43     _RL gXnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
44 adcroft 1.1 INTEGER myThid
45 cnh 1.5
46 cnh 1.17 C !LOCAL VARIABLES:
47 adcroft 1.1 C == Local variables ==
48     INTEGER i,j,k
49 cnh 1.17 _RL gYnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
50 heimbach 1.12 _RL a(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
51     _RL b(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
52     _RL c(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
53     _RL bet(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
54 cnh 1.6 _RL gam(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
55 cnh 1.17 CEOP
56 jmc 1.19
57     IF (Nr.LE.1) RETURN
58 adcroft 1.1
59 heimbach 1.12 C-- Initialise
60 heimbach 1.16 DO k=1,Nr
61     DO j=jMin,jMax
62     DO i=iMin,iMax
63     gYNm1(i,j,k,bi,bj) = 0. _d 0
64     ENDDO
65     ENDDO
66     ENDDO
67 heimbach 1.12
68     C-- Old aLower
69 heimbach 1.14 DO j=jMin,jMax
70     DO i=iMin,iMax
71 heimbach 1.12 a(i,j,1) = 0. _d 0
72     ENDDO
73     ENDDO
74     DO k=2,Nr
75 heimbach 1.14 DO j=jMin,jMax
76     DO i=iMin,iMax
77 heimbach 1.12 a(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
78     & *KappaRX(i,j, k )*recip_drC( k )
79 mlosch 1.18 IF (recip_hFac(i,j,k-1,bi,bj).EQ.0.) a(i,j,k)=0.
80 heimbach 1.12 ENDDO
81     ENDDO
82     ENDDO
83    
84     C-- Old aUpper
85     DO k=1,Nr-1
86 heimbach 1.14 DO j=jMin,jMax
87     DO i=iMin,iMax
88 heimbach 1.12 c(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
89     & *KappaRX(i,j,k+1)*recip_drC(k+1)
90 adcroft 1.13 IF (recip_hFac(i,j,k+1,bi,bj).EQ.0.) c(i,j,k)=0.
91 heimbach 1.12 ENDDO
92     ENDDO
93     ENDDO
94 heimbach 1.14 DO j=jMin,jMax
95     DO i=iMin,iMax
96 heimbach 1.12 c(i,j,Nr) = 0. _d 0
97     ENDDO
98     ENDDO
99    
100     C-- Old aCenter
101     DO k=1,Nr
102 heimbach 1.14 DO j=jMin,jMax
103     DO i=iMin,iMax
104 heimbach 1.12 b(i,j,k) = 1. _d 0 - c(i,j,k) - a(i,j,k)
105     ENDDO
106     ENDDO
107     ENDDO
108    
109     C-- Old and new gam, bet are the same
110     DO k=1,Nr
111 heimbach 1.14 DO j=jMin,jMax
112     DO i=iMin,iMax
113 heimbach 1.12 bet(i,j,k) = 0. _d 0
114     gam(i,j,k) = 0. _d 0
115     ENDDO
116     ENDDO
117     ENDDO
118    
119 heimbach 1.10 C-- Only need do anything if Nr>1
120     IF (Nr.GT.1) THEN
121    
122 heimbach 1.12 k = 1
123 cnh 1.5 C-- Beginning of forward sweep (top level)
124 adcroft 1.1 DO j=jMin,jMax
125     DO i=iMin,iMax
126 heimbach 1.12 IF (b(i,j,1).NE.0.) bet(i,j,1) = 1. _d 0 / b(i,j,1)
127 adcroft 1.1 ENDDO
128     ENDDO
129 heimbach 1.10
130 adcroft 1.1 ENDIF
131 heimbach 1.9
132 cnh 1.5 C-- Middle of forward sweep
133 jmc 1.20 IF (Nr.GE.2) THEN
134 heimbach 1.10
135 heimbach 1.12 CADJ loop = sequential
136     DO k=2,Nr
137 heimbach 1.9
138 adcroft 1.1 DO j=jMin,jMax
139     DO i=iMin,iMax
140 heimbach 1.12 gam(i,j,k) = c(i,j,k-1)*bet(i,j,k-1)
141     IF ( ( b(i,j,k) - a(i,j,k)*gam(i,j,k) ) .NE. 0.)
142     & bet(i,j,k) = 1. _d 0 / ( b(i,j,k) - a(i,j,k)*gam(i,j,k) )
143 adcroft 1.1 ENDDO
144     ENDDO
145 heimbach 1.9
146 adcroft 1.1 ENDDO
147 heimbach 1.10
148 adcroft 1.1 ENDIF
149 heimbach 1.10
150 heimbach 1.11
151 heimbach 1.12 DO j=jMin,jMax
152     DO i=iMin,iMax
153     gYNm1(i,j,1,bi,bj) = gXNm1(i,j,1,bi,bj)*bet(i,j,1)
154 heimbach 1.10 ENDDO
155 heimbach 1.12 ENDDO
156     DO k=2,Nr
157 heimbach 1.10 DO j=jMin,jMax
158     DO i=iMin,iMax
159 heimbach 1.12 gYnm1(i,j,k,bi,bj) = bet(i,j,k)*
160     & (gXnm1(i,j,k,bi,bj) - a(i,j,k)*gYnm1(i,j,k-1,bi,bj))
161 heimbach 1.9 ENDDO
162     ENDDO
163 heimbach 1.12 ENDDO
164 heimbach 1.9
165    
166 heimbach 1.12 C-- Backward sweep
167     CADJ loop = sequential
168     DO k=Nr-1,1,-1
169     DO j=jMin,jMax
170     DO i=iMin,iMax
171     gYnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
172     & -gam(i,j,k+1)*gYnm1(i,j,k+1,bi,bj)
173     ENDDO
174 adcroft 1.1 ENDDO
175     ENDDO
176 heimbach 1.9
177 heimbach 1.12 DO k=1,Nr
178 adcroft 1.1 DO j=jMin,jMax
179     DO i=iMin,iMax
180 heimbach 1.12 gXnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
181 adcroft 1.1 ENDDO
182     ENDDO
183     ENDDO
184    
185     RETURN
186     END
  ViewVC Help
Powered by ViewVC 1.1.22