/[MITgcm]/MITgcm/pkg/mom_vecinv/mom_vi_hdissip.F
ViewVC logotype

Contents of /MITgcm/pkg/mom_vecinv/mom_vi_hdissip.F

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


Revision 1.30 - (show annotations) (download)
Wed Jun 7 01:55:15 2006 UTC (17 years, 10 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint59, checkpoint58y_post, checkpoint58t_post, checkpoint58m_post, checkpoint58w_post, checkpoint58o_post, checkpoint58p_post, checkpoint58q_post, mitgcm_mapl_00, checkpoint58r_post, checkpoint58n_post, checkpoint59a, checkpoint59c, checkpoint59b, checkpoint58k_post, checkpoint58v_post, checkpoint58l_post, checkpoint58x_post, checkpoint58h_post, checkpoint58j_post, checkpoint58i_post, checkpoint58u_post, checkpoint58s_post
Changes since 1.29: +13 -13 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 C $Header: /u/gcmpack/MITgcm/pkg/mom_vecinv/mom_vi_hdissip.F,v 1.29 2006/04/05 21:14:45 mlosch Exp $
2 C $Name: $
3
4 #include "MOM_VECINV_OPTIONS.h"
5
6 SUBROUTINE MOM_VI_HDISSIP(
7 I bi,bj,k,
8 I hDiv,vort3,tension,strain,KE,
9 I hFacZ,dStar,zStar,
10 I viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
11 I harmonic,biharmonic,useVariableViscosity,
12 O uDissip,vDissip,
13 I myThid)
14
15 cph(
16 cph The following line was commented in the argument list
17 cph TAMC cannot digest commented lines within continuing lines
18 c I viscAh_Z,viscAh_D,viscA4_Z,viscA4_D,
19 cph)
20
21 IMPLICIT NONE
22 C
23 C Calculate horizontal dissipation terms
24 C [del^2 - del^4] (u,v)
25 C
26
27 C == Global variables ==
28 #include "SIZE.h"
29 #include "GRID.h"
30 #include "EEPARAMS.h"
31 #include "PARAMS.h"
32
33 C == Routine arguments ==
34 INTEGER bi,bj,k
35 _RL hDiv(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
36 _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
37 _RL tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
38 _RL strain(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
39 _RL KE(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
40 _RS hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
41 _RL dStar(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
42 _RL zStar(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
43 _RL uDissip(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
44 _RL vDissip(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
45 _RL viscAh_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
46 _RL viscAh_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
47 _RL viscA4_Z(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
48 _RL viscA4_D(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
49 LOGICAL harmonic, biharmonic, useVariableViscosity
50 INTEGER myThid
51
52 C == Local variables ==
53 INTEGER I,J
54 _RL Zip,Zij,Zpj,Dim,Dij,Dmj,uD2,vD2,uD4,vD4
55
56 C - Laplacian terms
57 IF (harmonic) THEN
58 DO j=2-Oly,sNy+Oly-1
59 DO i=2-Olx,sNx+Olx-1
60
61 Dim=hDiv( i ,j-1)
62 Dij=hDiv( i , j )
63 Dmj=hDiv(i-1, j )
64 Zip=hFacZ( i ,j+1)*vort3( i ,j+1)
65 Zij=hFacZ( i , j )*vort3( i , j )
66 Zpj=hFacZ(i+1, j )*vort3(i+1, j )
67
68 C This bit scales the harmonic dissipation operator to be proportional
69 C to the grid-cell area over the time-step. viscAh is then non-dimensional
70 C and should be less than 1/8, for example viscAh=0.01
71 IF (useVariableViscosity) THEN
72 Dij=Dij*viscAh_D(i,j)
73 Dim=Dim*viscAh_D(i,j-1)
74 Dmj=Dmj*viscAh_D(i-1,j)
75 Zij=Zij*viscAh_Z(i,j)
76 Zip=Zip*viscAh_Z(i,j+1)
77 Zpj=Zpj*viscAh_Z(i+1,j)
78 uD2 = (
79 & cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)
80 & -_recip_hFacW(i,j,k,bi,bj)*( Zip-Zij )*recip_DYG(i,j,bi,bj) )
81 #ifdef ISOTROPIC_COS_SCALING
82 & *cosFacU(j,bi,bj)
83 #endif /* ISOTROPIC_COS_SCALING */
84 vD2 = (
85 & _recip_hFacS(i,j,k,bi,bj)*( Zpj-Zij )*recip_DXG(i,j,bi,bj)
86 & *cosFacV(j,bi,bj)
87 & +( Dij-Dim )*recip_DYC(i,j,bi,bj) )
88 #ifdef ISOTROPIC_COS_SCALING
89 & *cosFacV(j,bi,bj)
90 #endif /* ISOTROPIC_COS_SCALING */
91 ELSE
92 uD2 = viscAhD*
93 & cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)
94 & - viscAhZ*_recip_hFacW(i,j,k,bi,bj)*
95 & ( Zip-Zij )*recip_DYG(i,j,bi,bj)
96 #ifdef ISOTROPIC_COS_SCALING
97 & *cosFacU(j,bi,bj)
98 #endif /* ISOTROPIC_COS_SCALING */
99 vD2 = viscAhZ*_recip_hFacS(i,j,k,bi,bj)*
100 & cosFacV(j,bi,bj)*( Zpj-Zij )*recip_DXG(i,j,bi,bj)
101 & + viscAhD* ( Dij-Dim )*recip_DYC(i,j,bi,bj)
102 #ifdef ISOTROPIC_COS_SCALING
103 & *cosFacV(j,bi,bj)
104 #endif /* ISOTROPIC_COS_SCALING */
105 ENDIF
106
107 uDissip(i,j) = uD2
108 vDissip(i,j) = vD2
109
110 ENDDO
111 ENDDO
112 ELSE
113 DO j=2-Oly,sNy+Oly-1
114 DO i=2-Olx,sNx+Olx-1
115 uDissip(i,j) = 0.
116 vDissip(i,j) = 0.
117 ENDDO
118 ENDDO
119 ENDIF
120
121 C - Bi-harmonic terms
122 IF (biharmonic) THEN
123 DO j=2-Oly,sNy+Oly-1
124 DO i=2-Olx,sNx+Olx-1
125
126 #ifdef MOM_VI_ORIGINAL_VISCA4
127 Dim=dyF( i ,j-1,bi,bj)*dStar( i ,j-1)
128 Dij=dyF( i , j ,bi,bj)*dStar( i , j )
129 Dmj=dyF(i-1, j ,bi,bj)*dStar(i-1, j )
130
131 Zip=dxV( i ,j+1,bi,bj)*hFacZ( i ,j+1)*zStar( i ,j+1)
132 Zij=dxV( i , j ,bi,bj)*hFacZ( i , j )*zStar( i , j )
133 Zpj=dxV(i+1, j ,bi,bj)*hFacZ(i+1, j )*zStar(i+1, j )
134 #else
135 Dim=dStar( i ,j-1)
136 Dij=dStar( i , j )
137 Dmj=dStar(i-1, j )
138
139 Zip=hFacZ( i ,j+1)*zStar( i ,j+1)
140 Zij=hFacZ( i , j )*zStar( i , j )
141 Zpj=hFacZ(i+1, j )*zStar(i+1, j )
142 #endif
143
144 C This bit scales the harmonic dissipation operator to be proportional
145 C to the grid-cell area over the time-step. viscAh is then non-dimensional
146 C and should be less than 1/8, for example viscAh=0.01
147 IF (useVariableViscosity) THEN
148 Dij=Dij*viscA4_D(i,j)
149 Dim=Dim*viscA4_D(i,j-1)
150 Dmj=Dmj*viscA4_D(i-1,j)
151 Zij=Zij*viscA4_Z(i,j)
152 Zip=Zip*viscA4_Z(i,j+1)
153 Zpj=Zpj*viscA4_Z(i+1,j)
154
155 #ifdef MOM_VI_ORIGINAL_VISCA4
156 uD4 = recip_rAw(i,j,bi,bj)*(
157 & ( (Dij-Dmj)*cosFacU(j,bi,bj) )
158 & -_recip_hFacW(i,j,k,bi,bj)*( Zip-Zij )
159 #ifdef ISOTROPIC_COS_SCALING
160 & *cosFacU(j,bi,bj)
161 #endif /* ISOTROPIC_COS_SCALING */
162 & )
163 vD4 = recip_rAs(i,j,bi,bj)*(
164 & _recip_hFacS(i,j,k,bi,bj)*( (Zpj-Zij)*cosFacV(j,bi,bj) )
165 & + ( Dij-Dim )
166 #ifdef ISOTROPIC_COS_SCALING
167 & *cosFacV(j,bi,bj)
168 #endif /* ISOTROPIC_COS_SCALING */
169 & )
170 ELSE
171 uD4 = recip_rAw(i,j,bi,bj)*(
172 & viscA4*( (Dij-Dmj)*cosFacU(j,bi,bj) )
173 & -_recip_hFacW(i,j,k,bi,bj)*viscA4*( Zip-Zij )
174 #ifdef ISOTROPIC_COS_SCALING
175 & *cosFacU(j,bi,bj)
176 #endif /* ISOTROPIC_COS_SCALING */
177 & )
178 vD4 = recip_rAs(i,j,bi,bj)*(
179 & _recip_hFacS(i,j,k,bi,bj)*viscA4*( (Zpj-Zij)*cosFacV(j,bi,bj) )
180 & + viscA4*( Dij-Dim )
181 #ifdef ISOTROPIC_COS_SCALING
182 & *cosFacV(j,bi,bj)
183 #endif /* ISOTROPIC_COS_SCALING */
184 & )
185 #else /* MOM_VI_ORIGINAL_VISCA4 */
186 uD4 = (
187 & cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)
188 & -_recip_hFacW(i,j,k,bi,bj)*( Zip-Zij )*recip_DYG(i,j,bi,bj) )
189 #ifdef ISOTROPIC_COS_SCALING
190 & *cosFacU(j,bi,bj)
191 #endif /* ISOTROPIC_COS_SCALING */
192 vD4 = (
193 & _recip_hFacS(i,j,k,bi,bj)*( Zpj-Zij )*recip_DXG(i,j,bi,bj)
194 & *cosFacV(j,bi,bj)
195 & +( Dij-Dim )*recip_DYC(i,j,bi,bj) )
196 #ifdef ISOTROPIC_COS_SCALING
197 & *cosFacV(j,bi,bj)
198 #endif /* ISOTROPIC_COS_SCALING */
199 ELSE
200 uD4 = viscA4D*
201 & cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)
202 & - viscA4Z*_recip_hFacW(i,j,k,bi,bj)*
203 & ( Zip-Zij )*recip_DYG(i,j,bi,bj)
204 #ifdef ISOTROPIC_COS_SCALING
205 & *cosFacU(j,bi,bj)
206 #endif /* ISOTROPIC_COS_SCALING */
207 vD4 = viscA4Z*_recip_hFacS(i,j,k,bi,bj)*
208 & cosFacV(j,bi,bj)*( Zpj-Zij )*recip_DXG(i,j,bi,bj)
209 & + viscA4D* ( Dij-Dim )*recip_DYC(i,j,bi,bj)
210 #ifdef ISOTROPIC_COS_SCALING
211 & *cosFacV(j,bi,bj)
212 #endif /* ISOTROPIC_COS_SCALING */
213 #endif /* MOM_VI_ORIGINAL_VISCA4 */
214 ENDIF
215
216 uDissip(i,j) = uDissip(i,j) - uD4
217 vDissip(i,j) = vDissip(i,j) - vD4
218
219 ENDDO
220 ENDDO
221 ENDIF
222
223 RETURN
224 END

  ViewVC Help
Powered by ViewVC 1.1.22