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

Contents of /MITgcm/pkg/mom_common/mom_hdissip.F

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


Revision 1.3 - (show annotations) (download)
Thu Apr 14 14:22:51 2005 UTC (19 years, 1 month ago) by baylor
Branch: MAIN
CVS Tags: checkpoint57o_post, checkpoint57m_post, checkpoint57k_post, checkpoint57g_post, checkpoint57i_post, checkpoint57g_pre, checkpoint57r_post, checkpoint57h_done, checkpoint57n_post, checkpoint57p_post, checkpoint57q_post, checkpoint57j_post, checkpoint57h_pre, checkpoint57l_post, checkpoint57h_post
Changes since 1.2: +15 -7 lines
Changed Smagorinsky length scale to harmonic mean of cell area, e.g.,
lenA2=(2*rA(i,j,bi,bj)/(dxF(I,J,bi,bj)+dyF(I,J,bi,bj)))**2

Should be more robust for anisotropic grids.

1 C $Header: /u/gcmpack/MITgcm/pkg/mom_common/mom_hdissip.F,v 1.2 2005/03/10 03:45:11 baylor Exp $
2 C $Name: $
3
4 #include "MOM_COMMON_OPTIONS.h"
5
6 SUBROUTINE MOM_HDISSIP(
7 I bi,bj,k,
8 I tension,strain,hFacZ,viscAt,viscAs,
9 O uDissip,vDissip,
10 I myThid)
11 IMPLICIT NONE
12 C
13 C Calculate horizontal dissipation terms in terms of tension and strain
14 C
15 C Du = d/dx At Tension + d/dy As Strain
16 C Dv = d/dx As Strain - d/dy At Tension
17
18 C == Global variables ==
19 #include "SIZE.h"
20 #include "EEPARAMS.h"
21 #include "GRID.h"
22 #include "PARAMS.h"
23
24 C == Routine arguments ==
25 INTEGER bi,bj,k
26 _RL tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
27 _RL strain(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
28 _RS hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
29 _RL viscAt
30 _RL viscAs
31 _RL uDissip(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
32 _RL vDissip(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
33 INTEGER myThid
34
35 C == Local variables ==
36 INTEGER I,J
37 _RL viscA_t(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
38 _RL viscA_s(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
39 _RL ASmag, smagfac
40 _RL vg2Min, vg2Max, AlinMax, AlinMin
41 _RL lenA2, lenAz2
42
43 IF (deltaTmom.NE.0.) THEN
44 vg2Min=viscAhGridMin/deltaTmom
45 vg2Max=viscAhGridMax/deltaTmom
46 ELSE
47 vg2Min=0.
48 vg2Max=0.
49 ENDIF
50
51 C - Calculate Smagorinsky Coefficients
52 smagfac=(viscC2smag/pi)**2
53 DO j=2-Oly,sNy+Oly-1
54 DO i=2-Olx,sNx+Olx-1
55 IF (viscC2smag.NE.0.) THEN
56 C Geometric Mean is used as length scale
57 lenA2=(2*rA(i,j,bi,bj)/
58 & (dxF(I,J,bi,bj)+dyF(I,J,bi,bj)))**2
59
60 Asmag=smagfac*lenA2
61 & *sqrt(tension(i,j)**2
62 & +0.25*(strain(i+1, j )**2+strain( i ,j+1)**2
63 & +strain(i-1, j )**2+strain( i ,j-1)**2))
64 viscA_t(i,j)=min(viscAhMax,max(viscAt,Asmag))
65
66 IF (vg2Max.GT.0.) THEN
67 AlinMax=vg2Max*lenA2
68 viscA_t(i,j)=min(AlinMax,viscA_t(i,j))
69 ENDIF
70 AlinMin=vg2Min*lenA2
71 viscA_t(i,j)=max(AlinMin,viscA_t(i,j))
72
73 C Geometric Mean is used as length scale
74 lenAz2=(2*rAz(i,j,bi,bj)/
75 & (dxV(I,J,bi,bj)+dyU(I,J,bi,bj)))**2
76 Asmag=smagfac*lenAz2
77 & *sqrt(strain(i,j)**2
78 & +0.25*(tension( i , j )**2+tension( i ,j+1)**2
79 & +tension(i+1, j )**2+tension(i+1,j+1)**2))
80 viscA_s(i,j)=min(viscAhMax,max(viscAs,Asmag))
81
82 IF (vg2Max.GT.0.) THEN
83 AlinMax=vg2Max*lenAz2
84 viscA_s(i,j)=min(AlinMax,viscA_s(i,j))
85 ENDIF
86 AlinMin=vg2Min*lenAz2
87 viscA_s(i,j)=max(AlinMin,viscA_s(i,j))
88
89 ELSE
90 viscA_t(i,j)=viscAt
91 viscA_s(i,j)=viscAs
92 ENDIF
93 ENDDO
94 ENDDO
95
96 C - Laplacian and bi-harmonic terms
97 DO j=2-Oly,sNy+Oly-1
98 DO i=2-Olx,sNx+Olx-1
99
100 uDissip(i,j) =
101 & recip_dyg(i,j,bi,bj)*recip_dyg(i,j,bi,bj)
102 & *recip_dxc(i,j,bi,bj)
103 & *(
104 & dyf( i , j ,bi,bj)*dyf( i , j ,bi,bj)
105 & *viscA_t( i , j )*tension( i , j )
106 & -dyf(i-1, j ,bi,bj)*dyf(i-1, j ,bi,bj)
107 & *viscA_t(i-1, j )*tension(i-1, j )
108 & )
109 & +recip_dxc(i,j,bi,bj)*recip_dxc(i,j,bi,bj)
110 & *recip_dyg(i,j,bi,bj)
111 & *(
112 & dxv( i ,j+1,bi,bj)*dxv( i ,j+1,bi,bj)
113 & *viscA_s(i,j+1)*strain( i ,j+1)
114 & -dxv( i , j ,bi,bj)*dxv( i , j ,bi,bj)
115 & *viscA_s(i, j )*strain( i , j )
116 & )
117
118 vDissip(i,j) =
119 & recip_dyc(i,j,bi,bj)*recip_dyc(i,j,bi,bj)
120 & *recip_dxg(i,j,bi,bj)
121 & *(
122 & dyu(i+1, j ,bi,bj)*dyu(i+1, j ,bi,bj)
123 & *viscA_s(i+1,j)*strain(i+1,j)
124 & -dyu( i , j ,bi,bj)*dyu( i , j ,bi,bj)
125 & *viscA_s( i ,j)*strain( i ,j)
126 & )
127 & -recip_dxg(i,j,bi,bj)*recip_dxg(i,j,bi,bj)
128 & *recip_dyc(i,j,bi,bj)
129 & *(
130 & dxf( i , j ,bi,bj)*dxf( i , j ,bi,bj)
131 & *viscA_t(i, j )*tension(i, j )
132 & -dxf( i ,j-1,bi,bj)*dxf( i ,j-1,bi,bj)
133 & *viscA_t(i,j-1)*tension(i,j-1)
134 & )
135
136 ENDDO
137 ENDDO
138
139 RETURN
140 END

  ViewVC Help
Powered by ViewVC 1.1.22