/[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.18 - (hide annotations) (download)
Mon Dec 16 21:40:07 2002 UTC (21 years, 5 months ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint51k_post, checkpoint47e_post, checkpoint50c_post, checkpoint52d_pre, checkpoint48e_post, checkpoint50c_pre, checkpoint51o_pre, checkpoint51l_post, checkpoint48i_post, checkpoint51, checkpoint50, checkpoint52, checkpoint50d_post, checkpoint50b_pre, checkpoint51f_post, checkpoint48b_post, checkpoint51d_post, checkpoint48c_pre, checkpoint51t_post, checkpoint51n_post, checkpoint51s_post, checkpoint48d_pre, checkpoint51j_post, checkpoint47i_post, checkpoint51n_pre, checkpoint47d_post, checkpoint48d_post, checkpoint48f_post, checkpoint52b_pre, checkpoint51l_pre, checkpoint48h_post, checkpoint51q_post, checkpoint51b_pre, checkpoint47g_post, checkpoint52b_post, checkpoint52c_post, checkpoint51h_pre, checkpoint48a_post, checkpoint50f_post, checkpoint50a_post, checkpoint50f_pre, checkpoint47j_post, branch-exfmods-tag, branchpoint-genmake2, checkpoint51r_post, checkpoint48c_post, checkpoint51i_post, checkpoint51b_post, checkpoint51c_post, checkpoint52d_post, checkpoint50g_post, checkpoint52a_pre, checkpoint50h_post, checkpoint50e_pre, checkpoint50i_post, checkpoint51i_pre, checkpoint47f_post, checkpoint50e_post, branch-netcdf, checkpoint50d_pre, checkpoint51e_post, checkpoint48, checkpoint49, checkpoint51o_post, checkpoint51f_pre, checkpoint48g_post, checkpoint47h_post, checkpoint52a_post, checkpoint51g_post, ecco_c52_e35, checkpoint50b_post, checkpoint51m_post, checkpoint51a_post, checkpoint51p_post, checkpoint51u_post
Branch point for: branch-exfmods-curt, branch-genmake2, branch-nonh, tg2-branch, netcdf-sm0, checkpoint51n_branch
Changes since 1.17: +2 -1 lines 
o fixed bug in the oceanic pressure coordinates code: vertical viscosity
  at the bottom boundary had an erroneous half slip boundary condition
o added bottom drag and no slip boundary condition capabilities to
  oceanic pressure coordinates code
1 mlosch 1.18 C $Header: /u/gcmpack/MITgcm/model/src/impldiff.F,v 1.17 2001/09/26 18:09:15 cnh 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 adcroft 1.1
57 heimbach 1.12 C-- Initialise
58 heimbach 1.16 DO k=1,Nr
59     DO j=jMin,jMax
60     DO i=iMin,iMax
61     gYNm1(i,j,k,bi,bj) = 0. _d 0
62     ENDDO
63     ENDDO
64     ENDDO
65 heimbach 1.12
66     C-- Old aLower
67 heimbach 1.14 DO j=jMin,jMax
68     DO i=iMin,iMax
69 heimbach 1.12 a(i,j,1) = 0. _d 0
70     ENDDO
71     ENDDO
72     DO k=2,Nr
73 heimbach 1.14 DO j=jMin,jMax
74     DO i=iMin,iMax
75 heimbach 1.12 a(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
76     & *KappaRX(i,j, k )*recip_drC( k )
77 mlosch 1.18 IF (recip_hFac(i,j,k-1,bi,bj).EQ.0.) a(i,j,k)=0.
78 heimbach 1.12 ENDDO
79     ENDDO
80     ENDDO
81    
82     C-- Old aUpper
83     DO k=1,Nr-1
84 heimbach 1.14 DO j=jMin,jMax
85     DO i=iMin,iMax
86 heimbach 1.12 c(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
87     & *KappaRX(i,j,k+1)*recip_drC(k+1)
88 adcroft 1.13 IF (recip_hFac(i,j,k+1,bi,bj).EQ.0.) c(i,j,k)=0.
89 heimbach 1.12 ENDDO
90     ENDDO
91     ENDDO
92 heimbach 1.14 DO j=jMin,jMax
93     DO i=iMin,iMax
94 heimbach 1.12 c(i,j,Nr) = 0. _d 0
95     ENDDO
96     ENDDO
97    
98     C-- Old aCenter
99     DO k=1,Nr
100 heimbach 1.14 DO j=jMin,jMax
101     DO i=iMin,iMax
102 heimbach 1.12 b(i,j,k) = 1. _d 0 - c(i,j,k) - a(i,j,k)
103     ENDDO
104     ENDDO
105     ENDDO
106    
107     C-- Old and new gam, bet are the same
108     DO k=1,Nr
109 heimbach 1.14 DO j=jMin,jMax
110     DO i=iMin,iMax
111 heimbach 1.12 bet(i,j,k) = 0. _d 0
112     gam(i,j,k) = 0. _d 0
113     ENDDO
114     ENDDO
115     ENDDO
116    
117 heimbach 1.10 C-- Only need do anything if Nr>1
118     IF (Nr.GT.1) THEN
119    
120 heimbach 1.12 k = 1
121 cnh 1.5 C-- Beginning of forward sweep (top level)
122 adcroft 1.1 DO j=jMin,jMax
123     DO i=iMin,iMax
124 heimbach 1.12 IF (b(i,j,1).NE.0.) bet(i,j,1) = 1. _d 0 / b(i,j,1)
125 adcroft 1.1 ENDDO
126     ENDDO
127 heimbach 1.10
128 adcroft 1.1 ENDIF
129 heimbach 1.9
130 cnh 1.5 C-- Middle of forward sweep
131 cnh 1.6 IF (Nr.GT.2) THEN
132 heimbach 1.10
133 heimbach 1.12 CADJ loop = sequential
134     DO k=2,Nr
135 heimbach 1.9
136 adcroft 1.1 DO j=jMin,jMax
137     DO i=iMin,iMax
138 heimbach 1.12 gam(i,j,k) = c(i,j,k-1)*bet(i,j,k-1)
139     IF ( ( b(i,j,k) - a(i,j,k)*gam(i,j,k) ) .NE. 0.)
140     & bet(i,j,k) = 1. _d 0 / ( b(i,j,k) - a(i,j,k)*gam(i,j,k) )
141 adcroft 1.1 ENDDO
142     ENDDO
143 heimbach 1.9
144 adcroft 1.1 ENDDO
145 heimbach 1.10
146 adcroft 1.1 ENDIF
147 heimbach 1.10
148 heimbach 1.11
149 heimbach 1.12 DO j=jMin,jMax
150     DO i=iMin,iMax
151     gYNm1(i,j,1,bi,bj) = gXNm1(i,j,1,bi,bj)*bet(i,j,1)
152 heimbach 1.10 ENDDO
153 heimbach 1.12 ENDDO
154     DO k=2,Nr
155 heimbach 1.10 DO j=jMin,jMax
156     DO i=iMin,iMax
157 heimbach 1.12 gYnm1(i,j,k,bi,bj) = bet(i,j,k)*
158     & (gXnm1(i,j,k,bi,bj) - a(i,j,k)*gYnm1(i,j,k-1,bi,bj))
159 heimbach 1.9 ENDDO
160     ENDDO
161 heimbach 1.12 ENDDO
162 heimbach 1.9
163    
164 heimbach 1.12 C-- Backward sweep
165     CADJ loop = sequential
166     DO k=Nr-1,1,-1
167     DO j=jMin,jMax
168     DO i=iMin,iMax
169     gYnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
170     & -gam(i,j,k+1)*gYnm1(i,j,k+1,bi,bj)
171     ENDDO
172 adcroft 1.1 ENDDO
173     ENDDO
174 heimbach 1.9
175 heimbach 1.12 DO k=1,Nr
176 adcroft 1.1 DO j=jMin,jMax
177     DO i=iMin,iMax
178 heimbach 1.12 gXnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
179 adcroft 1.1 ENDDO
180     ENDDO
181     ENDDO
182    
183     RETURN
184     END
  ViewVC Help
Powered by ViewVC 1.1.22