/[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.13 - (hide annotations) (download)
Wed Nov 29 22:38:31 2000 UTC (23 years, 6 months ago) by adcroft
Branch: MAIN
CVS Tags: branch-atmos-merge-freeze, branch-atmos-merge-start, branch-atmos-merge-shapiro, checkpoint33, branch-atmos-merge-zonalfilt, branch-atmos-merge-phase5, branch-atmos-merge-phase4, branch-atmos-merge-phase7, branch-atmos-merge-phase6, branch-atmos-merge-phase1, branch-atmos-merge-phase3, branch-atmos-merge-phase2
Branch point for: branch-atmos-merge
Changes since 1.12: +5 -1 lines 
Bug fix for lower boundaries over topography of non-zero
height. Thanks to CAE for spotting it.
1 adcroft 1.13 C $Header: /u/gcmpack/models/MITgcmUV/model/src/impldiff.F,v 1.12 2000/11/13 16:30:32 heimbach Exp $
2 adcroft 1.1
3 cnh 1.7 #include "CPP_OPTIONS.h"
4 adcroft 1.1
5     C /==========================================================\
6     C | S/R IMPLDIFF |
7 cnh 1.5 C | o Solve implicit diffusion equation for vertical |
8     C | diffusivity. |
9 adcroft 1.13 C | o Recoded from 2d intermediate fields to 3d to reduce |
10     C | TAMC storage |
11     C | o Fixed missing masks for fields a(), c() |
12 adcroft 1.1 C \==========================================================/
13     SUBROUTINE IMPLDIFF( bi, bj, iMin, iMax, jMin, jMax,
14 adcroft 1.8 I deltaTX,KappaRX,recip_hFac,
15     U gXNm1,
16 adcroft 1.1 I myThid )
17 cnh 1.5 IMPLICIT NONE
18     C == Global data ==
19 adcroft 1.1 #include "SIZE.h"
20     #include "DYNVARS.h"
21 cnh 1.2 #include "EEPARAMS.h"
22 adcroft 1.1 #include "PARAMS.h"
23     #include "GRID.h"
24 cnh 1.5
25 heimbach 1.9 #ifdef ALLOW_AUTODIFF_TAMC
26     #include "tamc_keys.h"
27     #endif
28    
29 adcroft 1.1 C == Routine Arguments ==
30     INTEGER bi,bj,iMin,iMax,jMin,jMax
31 adcroft 1.8 _RL deltaTX
32     _RL KappaRX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
33     _RS recip_hFac(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
34     _RL gXnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
35 heimbach 1.12 _RL gYnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
36 adcroft 1.1 INTEGER myThid
37 cnh 1.5
38 adcroft 1.1 C == Local variables ==
39     INTEGER i,j,k
40 heimbach 1.12 _RL a(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
41     _RL b(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
42     _RL c(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
43     _RL bet(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
44 cnh 1.6 _RL gam(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
45 adcroft 1.1
46 heimbach 1.9 #ifdef ALLOW_AUTODIFF_TAMC
47     INTEGER kkey
48     #endif
49    
50 heimbach 1.12 C-- Initialise
51    
52     C-- Old aLower
53     DO j=1-Oly,sNy+Oly
54     DO i=1-Olx,sNx+Olx
55     a(i,j,1) = 0. _d 0
56     ENDDO
57     ENDDO
58     DO k=2,Nr
59     DO j=1-Oly,sNy+Oly
60     DO i=1-Olx,sNx+Olx
61     a(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
62     & *KappaRX(i,j, k )*recip_drC( k )
63     ENDDO
64     ENDDO
65     ENDDO
66    
67     C-- Old aUpper
68     DO k=1,Nr-1
69     DO j=1-Oly,sNy+Oly
70     DO i=1-Olx,sNx+Olx
71     c(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
72     & *KappaRX(i,j,k+1)*recip_drC(k+1)
73 adcroft 1.13 IF (recip_hFac(i,j,k+1,bi,bj).EQ.0.) c(i,j,k)=0.
74 heimbach 1.12 ENDDO
75     ENDDO
76     ENDDO
77     DO j=1-Oly,sNy+Oly
78     DO i=1-Olx,sNx+Olx
79     c(i,j,Nr) = 0. _d 0
80     ENDDO
81     ENDDO
82    
83     C-- Old aCenter
84     DO k=1,Nr
85     DO j=1-Oly,sNy+Oly
86     DO i=1-Olx,sNx+Olx
87     b(i,j,k) = 1. _d 0 - c(i,j,k) - a(i,j,k)
88     ENDDO
89     ENDDO
90     ENDDO
91    
92     C-- Old and new gam, bet are the same
93     DO k=1,Nr
94     DO j=1-Oly,sNy+Oly
95     DO i=1-Olx,sNx+Olx
96     bet(i,j,k) = 0. _d 0
97     gam(i,j,k) = 0. _d 0
98     ENDDO
99     ENDDO
100     ENDDO
101    
102 heimbach 1.10 C-- Only need do anything if Nr>1
103     IF (Nr.GT.1) THEN
104    
105 heimbach 1.12 k = 1
106 cnh 1.5 C-- Beginning of forward sweep (top level)
107 adcroft 1.1 DO j=jMin,jMax
108     DO i=iMin,iMax
109 heimbach 1.12 IF (b(i,j,1).NE.0.) bet(i,j,1) = 1. _d 0 / b(i,j,1)
110 adcroft 1.1 ENDDO
111     ENDDO
112 heimbach 1.10
113 adcroft 1.1 ENDIF
114 heimbach 1.9
115 cnh 1.5 C-- Middle of forward sweep
116 cnh 1.6 IF (Nr.GT.2) THEN
117 heimbach 1.10
118 heimbach 1.12 CADJ loop = sequential
119     DO k=2,Nr
120 heimbach 1.9
121     #ifdef ALLOW_AUTODIFF_TAMC
122     kkey = (idkey-1)*(Nr-2) + k-1
123     #endif
124    
125 adcroft 1.1 DO j=jMin,jMax
126     DO i=iMin,iMax
127 heimbach 1.12 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 adcroft 1.1 ENDDO
131     ENDDO
132 heimbach 1.9
133 adcroft 1.1 ENDDO
134 heimbach 1.10
135 adcroft 1.1 ENDIF
136 heimbach 1.10
137 heimbach 1.11
138 heimbach 1.12 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 heimbach 1.10 ENDDO
142 heimbach 1.12 ENDDO
143     DO k=2,Nr
144 heimbach 1.10 DO j=jMin,jMax
145     DO i=iMin,iMax
146 heimbach 1.12 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 heimbach 1.9 ENDDO
149     ENDDO
150 heimbach 1.12 ENDDO
151 heimbach 1.9
152    
153 heimbach 1.12 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 adcroft 1.1 ENDDO
162     ENDDO
163 heimbach 1.9
164 heimbach 1.12 DO k=1,Nr
165 adcroft 1.1 DO j=jMin,jMax
166     DO i=iMin,iMax
167 heimbach 1.12 gXnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
168 adcroft 1.1 ENDDO
169     ENDDO
170     ENDDO
171    
172     RETURN
173     END
  ViewVC Help
Powered by ViewVC 1.1.22