/[MITgcm]/MITgcm/pkg/mom_common/mom_v_implicit_r.F
ViewVC logotype

Annotation of /MITgcm/pkg/mom_common/mom_v_implicit_r.F

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

Revision 1.4 - (hide annotations) (download)
Wed Jun 7 01:55:15 2006 UTC (18 years ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint58l_post, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, mitgcm_mapl_00, checkpoint58u_post, checkpoint58w_post, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint64, checkpoint65, checkpoint60, checkpoint61, checkpoint62, checkpoint63, checkpoint58r_post, checkpoint58n_post, checkpoint58x_post, checkpoint58t_post, checkpoint58h_post, checkpoint65b, checkpoint65a, checkpoint58q_post, checkpoint59q, checkpoint59p, checkpoint59r, checkpoint58j_post, checkpoint59e, checkpoint59d, checkpoint59g, checkpoint59f, checkpoint59a, checkpoint59c, checkpoint59b, checkpoint59m, checkpoint59l, checkpoint59o, checkpoint59n, checkpoint59i, checkpoint59h, checkpoint59k, checkpoint59j, checkpoint59, checkpoint58i_post, checkpoint58o_post, checkpoint62c, checkpoint62b, checkpoint62a, checkpoint62g, checkpoint62f, checkpoint62e, checkpoint62d, checkpoint62k, checkpoint62j, checkpoint62i, checkpoint62h, checkpoint62o, checkpoint62n, checkpoint62m, checkpoint62l, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62p, checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62z, checkpoint62y, checkpoint62x, checkpoint58y_post, checkpoint58k_post, checkpoint58v_post, checkpoint58s_post, checkpoint61f, checkpoint61g, checkpoint61d, checkpoint61e, checkpoint61b, checkpoint61c, checkpoint58p_post, checkpoint61a, checkpoint61n, checkpoint61o, checkpoint61l, checkpoint61m, checkpoint61j, checkpoint61k, checkpoint61h, checkpoint61i, checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61r, checkpoint61s, checkpoint61p, checkpoint61q, checkpoint61z, checkpoint61x, checkpoint61y, checkpoint58m_post
Changes since 1.3: +11 -11 lines 
Modifications for bottom topography control
o replace hFacC by _hFacC at various places
o replace ALLOW_HFACC_CONTROL by ALLOW_DEPTH_CONTROL
o add non-self-adjoint cg2d_nsa
o update autodiff support routines
o re-initialise hfac after ctrl_depth_ini
o works for 5x5 box, doesnt work for global_ocean.90x40x15
1 heimbach 1.4 C $Header: /u/gcmpack/MITgcm/pkg/mom_common/mom_v_implicit_r.F,v 1.3 2005/09/04 19:22:48 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "MOM_COMMON_OPTIONS.h"
5    
6     CBOP
7     C !ROUTINE: MOM_V_IMPLICIT_R
8     C !INTERFACE:
9     SUBROUTINE MOM_V_IMPLICIT_R(
10     I kappaRV,
11     I bi, bj, myTime, myIter, myThid )
12     C !DESCRIPTION: \bv
13     C *==========================================================*
14     C | S/R MOM_V_IMPLICIT_R
15     C | o Solve implicitly vertical advection & diffusion
16     C | of momentum, meridional component
17     C *==========================================================*
18     C *==========================================================*
19     C \ev
20    
21     C !USES:
22     IMPLICIT NONE
23     C == Global data ==
24     #include "SIZE.h"
25     #include "EEPARAMS.h"
26     #include "PARAMS.h"
27     #include "GRID.h"
28     #include "DYNVARS.h"
29    
30     C !INPUT/OUTPUT PARAMETERS:
31     C == Routine Arguments ==
32     _RL kappaRV(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
33     INTEGER bi, bj
34     _RL myTime
35     INTEGER myIter, myThid
36    
37     C !LOCAL VARIABLES:
38     C == Local variables ==
39     INTEGER iMin,iMax,jMin,jMax
40     INTEGER i,j,k
41     INTEGER diagonalNumber, errCode
42     c _RL a5d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
43     _RL b5d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
44     _RL c5d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
45     _RL d5d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
46     c _RL e5d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
47     _RL rTrans(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
48 jmc 1.2 _RL rCenter, rUpwind, upwindFac
49 jmc 1.3 #ifdef ALLOW_DIAGNOSTICS
50     CHARACTER*8 diagName
51     LOGICAL DIAGNOSTICS_IS_ON
52     EXTERNAL DIAGNOSTICS_IS_ON
53     _RL vf (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
54     #endif
55 jmc 1.1 CEOP
56    
57     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
58     C Solve for V-component :
59     C----------------------------
60    
61     C-- Initialise
62     iMin = 1
63     jMin = 1
64     iMax = sNx
65     jMax = sNy+1
66     DO k=1,Nr
67     DO j=1-Oly,sNy+Oly
68     DO i=1-Olx,sNx+Olx
69     c a5d(i,j,k) = 0. _d 0
70     b5d(i,j,k) = 0. _d 0
71     c5d(i,j,k) = 1. _d 0
72     d5d(i,j,k) = 0. _d 0
73     c e5d(i,j,k) = 0. _d 0
74     ENDDO
75     ENDDO
76     ENDDO
77     diagonalNumber = 1
78    
79     IF ( implicitViscosity .AND. Nr.GT.1 ) THEN
80    
81     C-- set the tri-diagonal matrix to solve the implicit viscosity
82     diagonalNumber = 3
83     C- 1rst lower diagonal :
84     DO k=2,Nr
85     DO j=jMin,jMax
86     DO i=iMin,iMax
87     IF (maskS(i,j,k-1,bi,bj).EQ.1.)
88     & b5d(i,j,k) = -deltaTmom
89 heimbach 1.4 & *_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
90 jmc 1.1 & *kappaRV(i,j, k )*recip_drC( k )
91     ENDDO
92     ENDDO
93     ENDDO
94     C- 1rst upper diagonal :
95     DO k=1,Nr-1
96     DO j=jMin,jMax
97     DO i=iMin,iMax
98     IF (maskS(i,j,k+1,bi,bj).EQ.1.)
99     & d5d(i,j,k) = -deltaTmom
100 heimbach 1.4 & *_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
101 jmc 1.1 & *KappaRV(i,j,k+1)*recip_drC(k+1)
102     ENDDO
103     ENDDO
104     ENDDO
105     C- Main diagonal :
106     DO k=1,Nr
107     DO j=jMin,jMax
108     DO i=iMin,iMax
109     c5d(i,j,k) = 1. _d 0 - b5d(i,j,k) - d5d(i,j,k)
110     ENDDO
111     ENDDO
112     ENDDO
113    
114     C-- end if implicitDiffusion
115     ENDIF
116    
117     IF ( momImplVertAdv .AND. Nr.GT.1 ) THEN
118    
119     diagonalNumber = 3
120     DO k=2,Nr
121    
122     DO j=jMin,jMax
123     DO i=iMin,iMax
124     rTrans(i,j) = 0.5 _d 0 * (
125     & wVel(i, j ,k,bi,bj)*rA(i, j ,bi,bj)
126     & *maskC(i, j ,k-1,bi,bj)
127     & + wVel(i,j-1,k,bi,bj)*rA(i,j-1,bi,bj)
128     & *maskC(i,j-1,k-1,bi,bj)
129     & )
130     ENDDO
131     ENDDO
132    
133     IF ( vectorInvariantMomentum ) THEN
134 jmc 1.2 C- space Centered/Upwind advection scheme, Advective form:
135     IF ( upwindShear ) THEN
136     upwindFac = 1. _d 0
137     ELSE
138     upwindFac = 0. _d 0
139     ENDIF
140 jmc 1.1 DO j=jMin,jMax
141     DO i=iMin,iMax
142 jmc 1.2 rCenter = 0.5 _d 0 *deltaTmom*rTrans(i,j)
143     & *recip_rAs(i,j,bi,bj)*rkSign
144     rUpwind = ABS(rCenter)*upwindFac
145 jmc 1.1 b5d(i,j,k) = b5d(i,j,k)
146 jmc 1.2 & - (rCenter+rUpwind)
147 heimbach 1.4 & *_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
148 jmc 1.1 c5d(i,j,k) = c5d(i,j,k)
149 jmc 1.2 & + (rCenter+rUpwind)
150 heimbach 1.4 & *_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
151 jmc 1.1 c5d(i,j,k-1) = c5d(i,j,k-1)
152 jmc 1.2 & - (rCenter-rUpwind)
153 heimbach 1.4 & *_recip_hFacS(i,j,k-1,bi,bj)*recip_drF(k-1)
154 jmc 1.1 d5d(i,j,k-1) = d5d(i,j,k-1)
155 jmc 1.2 & + (rCenter-rUpwind)
156 heimbach 1.4 & *_recip_hFacS(i,j,k-1,bi,bj)*recip_drF(k-1)
157 jmc 1.1 ENDDO
158     ENDDO
159     ELSE
160     C- space Centered advection scheme, Flux form:
161     DO j=jMin,jMax
162     DO i=iMin,iMax
163 jmc 1.2 rCenter = 0.5 _d 0 *deltaTmom*rTrans(i,j)
164     & *recip_rAs(i,j,bi,bj)*rkSign
165 jmc 1.1 b5d(i,j,k) = b5d(i,j,k)
166 heimbach 1.4 & - rCenter*_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
167 jmc 1.1 c5d(i,j,k) = c5d(i,j,k)
168 heimbach 1.4 & - rCenter*_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)
169 jmc 1.1 c5d(i,j,k-1) = c5d(i,j,k-1)
170 heimbach 1.4 & + rCenter*_recip_hFacS(i,j,k-1,bi,bj)*recip_drF(k-1)
171 jmc 1.1 d5d(i,j,k-1) = d5d(i,j,k-1)
172 heimbach 1.4 & + rCenter*_recip_hFacS(i,j,k-1,bi,bj)*recip_drF(k-1)
173 jmc 1.1 ENDDO
174     ENDDO
175     STOP 'MOM_IMPLICIT_R: Flux Form not yet finished.'
176     ENDIF
177    
178     C-- end k loop
179     ENDDO
180    
181     C-- end if momImplVertAdv
182     ENDIF
183    
184     IF ( diagonalNumber .EQ. 3 ) THEN
185     C-- Solve tri-diagonal system :
186     CALL SOLVE_TRIDIAGONAL( iMin,iMax, jMin,jMax,
187     I b5d, c5d, d5d,
188     U gV,
189     O errCode,
190     I bi, bj, myThid )
191     IF (errCode.GE.1) THEN
192     STOP 'MOM_IMPLICIT_R: error when solving 3-Diag problem.'
193     ENDIF
194     ELSEIF ( diagonalNumber .NE. 1 ) THEN
195     STOP 'MOM_IMPLICIT_R: no solver available.'
196     ENDIF
197    
198 jmc 1.3 #ifdef ALLOW_DIAGNOSTICS
199     C-- Diagnostics of vertical viscous flux:
200     IF ( useDiagnostics .AND. implicitViscosity ) THEN
201     diagName = 'VISrI_Vm'
202     IF ( DIAGNOSTICS_IS_ON(diagName,myThid) ) THEN
203     DO k= 1,Nr
204     IF ( k.EQ.1 ) THEN
205     C- Note: Needs to call DIAGNOSTICS_FILL at level k=1 even if array == 0
206     C otherwise counter is not incremented !!
207     DO j=1-OLy,sNy+OLy
208     DO i=1-OLx,sNx+OLx
209     vf(i,j) = 0. _d 0
210     ENDDO
211     ENDDO
212     ELSE
213     DO j=jMin,jMax
214     DO i=iMin,iMax
215     vf(i,j) =
216     & -KappaRV(i,j,k)*rAs(i,j,bi,bj)*recip_drC(k)
217     & * (gV(i,j,k,bi,bj) - gV(i,j,k-1,bi,bj))*rkSign
218     & *_maskS(i,j,k,bi,bj)
219     & *_maskS(i,j,k-1,bi,bj)
220     ENDDO
221     ENDDO
222     ENDIF
223     CALL DIAGNOSTICS_FILL(vf,diagName, k,1, 2,bi,bj, myThid)
224     ENDDO
225     ENDIF
226     ENDIF
227     #endif /* ALLOW_DIAGNOSTICS */
228    
229 jmc 1.1 RETURN
230     END
  ViewVC Help
Powered by ViewVC 1.1.22