/[MITgcm]/MITgcm/pkg/gmredi/gmredi_slope_limit.F
ViewVC logotype

Annotation of /MITgcm/pkg/gmredi/gmredi_slope_limit.F

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


Revision 1.1 - (hide annotations) (download)
Wed Jun 21 19:42:54 2000 UTC (24 years ago) by adcroft
Branch: MAIN
CVS Tags: checkpoint29, branch-atmos-merge-freeze, branch-atmos-merge-start, branch-atmos-merge-shapiro, checkpoint33, checkpoint32, checkpoint31, checkpoint30, 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
Packaged GM/Redi routines.

1 adcroft 1.1 C $Header: $
2    
3     #include "GMREDI_OPTIONS.h"
4    
5     CStartOfInterface
6     SUBROUTINE GMREDI_SLOPE_LIMIT(
7     I dSigmadRReal,
8     I depthZ,
9     U SlopeX, SlopeY,
10     O dRdSigmaLtd,
11     I bi,bj, myThid )
12     C /==========================================================\
13     C | SUBROUTINE GMREDI_SLOPE_LIMIT |
14     C | o Calculate slopes for use in GM/Redi tensor |
15     C |==========================================================|
16     C | On entry: |
17     C | dSigmadRReal conatins the d/dz Sigma |
18     C | SlopeX/Y contains X/Y gradients of sigma |
19     C | depthZ conatins the height (m) of level |
20     C | On exit: |
21     C | SlopeX/Y contains X/Y slopes |
22     C | dRdSigmaLtd conatins the effective dSig/dz |
23     C \==========================================================/
24     IMPLICIT NONE
25    
26     C == Global variables ==
27     #include "SIZE.h"
28     #include "EEPARAMS.h"
29     #include "GMREDI.h"
30     #include "PARAMS.h"
31    
32     C == Routine arguments ==
33     C
34     _RL SlopeX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
35     _RL SlopeY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
36     _RL dSigmaDrReal(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
37     _RL dRdSigmaLtd(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
38     _RL depthZ
39     INTEGER bi,bj,myThid
40     CEndOfInterface
41    
42     #ifdef ALLOW_GMREDI
43    
44     C == Local variables ==
45     _RL Small_Number
46     PARAMETER(Small_Number=1.E-12)
47     _RL gradSmod,stratLimit,f1,Smod,f2,Rnondim,Cspd,Lrho
48     _RL dSigmaDrLtd
49     _RL fpi
50     PARAMETER(fpi=3.141592653589793047592)
51     INTEGER i,j
52    
53    
54     if (GM_taper_scheme.EQ.'orig') then
55    
56     C- Original implementation in mitgcmuv
57     C (this turns out to be the same as Cox slope clipping)
58     DO j=1-Oly+1,sNy+Oly-1
59     DO i=1-Olx+1,sNx+Olx-1
60    
61     gradSmod=SlopeX(i,j)*SlopeX(i,j)
62     & +SlopeY(i,j)*SlopeY(i,j)
63     if (gradSmod .NE. 0.) gradSmod=sqrt(gradSmod)
64    
65     stratLimit=-Small_Number-gradSmod*GM_rMaxSlope
66     if (dSigmaDrReal(i,j).LT.stratLimit) then
67     dSigmaDrLtd=dSigmaDrReal(i,j)
68     else
69     dSigmaDrLtd=stratLimit
70     endif
71     dRdSigmaLtd(i,j)=1./dSigmaDrLtd
72     SlopeX(i,j)=-SlopeX(i,j)*dRdSigmaLtd(i,j)
73     SlopeY(i,j)=-SlopeY(i,j)*dRdSigmaLtd(i,j)
74    
75     ENDDO
76     ENDDO
77    
78     elseif (GM_taper_scheme.EQ.'clipping') then
79    
80     C- Cox 1987 "Slope clipping"
81     DO j=1-Oly+1,sNy+Oly-1
82     DO i=1-Olx+1,sNx+Olx-1
83    
84     gradSmod=SlopeX(i,j)*SlopeX(i,j)
85     & +SlopeY(i,j)*SlopeY(i,j)
86     if (gradSmod .NE. 0.) gradSmod=sqrt(gradSmod)
87    
88     dSigmaDrLtd=-(Small_Number+gradSmod*GM_rMaxSlope)
89     if (dSigmaDrReal(i,j).LT.dSigmaDrLtd)
90     & dSigmaDrLtd=dSigmaDrReal(i,j)
91     dRdSigmaLtd(i,j)=1./dSigmaDrLtd
92     SlopeX(i,j)=-SlopeX(i,j)*dRdSigmaLtd(i,j)
93     SlopeY(i,j)=-SlopeY(i,j)*dRdSigmaLtd(i,j)
94    
95     ENDDO
96     ENDDO
97    
98     elseif (GM_taper_scheme.EQ.'gkw91') then
99    
100     C- Gerdes, Koberle and Willebrand, Clim. Dyn. 1991
101     DO j=1-Oly+1,sNy+Oly-1
102     DO i=1-Olx+1,sNx+Olx-1
103    
104     gradSmod=SlopeX(i,j)*SlopeX(i,j)
105     & +SlopeY(i,j)*SlopeY(i,j)
106     if (gradSmod .NE. 0.) gradSmod=sqrt(gradSmod)
107    
108     dSigmaDrLtd=dSigmaDrReal(i,j)
109     if (dSigmaDrLtd.NE.0.) then
110     dRdSigmaLtd(i,j)=1./dSigmaDrLtd
111     else
112     dRdSigmaLtd(i,j)=0.
113     endif
114     if (gradSmod.LE.GM_maxSlope*abs(dSigmaDrReal(i,j))) then
115     C If the slope is < Smax then calculate slopes properly
116     SlopeX(i,j)=-SlopeX(i,j)*dRdSigmaLtd(i,j)
117     SlopeY(i,j)=-SlopeY(i,j)*dRdSigmaLtd(i,j)
118     else
119     C If the slope is > Smax then taper slopes
120     SlopeX(i,j)=-SlopeX(i,j)*dSigmaDrLtd
121     & *((GM_maxSlope/gradSmod)**2)
122     SlopeY(i,j)=-SlopeY(i,j)*dSigmaDrLtd
123     & *((GM_maxSlope/gradSmod)**2)
124     endif
125    
126     ENDDO
127     ENDDO
128    
129     elseif (GM_taper_scheme.EQ.'dm95') then
130    
131     C- Danabasoglu and McWilliams, J. Clim. 1995
132     DO j=1-Oly+1,sNy+Oly-1
133     DO i=1-Olx+1,sNx+Olx-1
134    
135     dSigmaDrLtd=dSigmaDrReal(i,j)
136     if (dSigmaDrLtd.NE.0.) then
137     dRdSigmaLtd(i,j)=1./dSigmaDrLtd
138     else
139     dRdSigmaLtd(i,j)=0.
140     endif
141     SlopeX(i,j)=-SlopeX(i,j)*dRdSigmaLtd(i,j)
142     SlopeY(i,j)=-SlopeY(i,j)*dRdSigmaLtd(i,j)
143     Smod=SlopeX(i,j)*SlopeX(i,j)+SlopeY(i,j)*SlopeY(i,j)
144     if (Smod.NE.0.) Smod=sqrt(Smod)
145     f1=0.5*(1.+tanh( (GM_Scrit-Smod)/GM_Sd ))
146     SlopeX(i,j)=SlopeX(i,j)*f1
147     SlopeY(i,j)=SlopeY(i,j)*f1
148    
149     ENDDO
150     ENDDO
151    
152     elseif (GM_taper_scheme.EQ.'ldd97') then
153    
154     C- Large, Danabasoglu and Doney, JPO 1997
155     DO j=1-Oly+1,sNy+Oly-1
156     DO i=1-Olx+1,sNx+Olx-1
157    
158     dSigmaDrLtd=dSigmaDrReal(i,j)
159     if (dSigmaDrLtd.NE.0.) then
160     dRdSigmaLtd(i,j)=1./dSigmaDrLtd
161     else
162     dRdSigmaLtd(i,j)=0.
163     endif
164     SlopeX(i,j)=-SlopeX(i,j)*dRdSigmaLtd(i,j)
165     SlopeY(i,j)=-SlopeY(i,j)*dRdSigmaLtd(i,j)
166     Smod=SlopeX(i,j)*SlopeX(i,j)+SlopeY(i,j)*SlopeY(i,j)
167     if (Smod.NE.0.) Smod=sqrt(Smod)
168     f1=0.5*(1.+tanh( (GM_Scrit-Smod)/GM_Sd ))
169     Cspd=2.
170     Lrho=100.e3
171     if (FCori(i,j,bi,bj).NE.0.) Lrho=Cspd/abs(Fcori(i,j,bi,bj))
172     Lrho=min(Lrho,100.e3)
173     Lrho=max(Lrho,15.e3)
174     if (Smod.NE.0.) then
175     Rnondim=depthZ/(Lrho*Smod)
176     else
177     Rnondim=0.
178     endif
179     f2=0.5*(1.+sin( fpi*(Rnondim-0.5)))
180     SlopeX(i,j)=SlopeX(i,j)*f1
181     SlopeY(i,j)=SlopeY(i,j)*f1
182    
183     ENDDO
184     ENDDO
185    
186     elseif (GM_taper_scheme.EQ.' ') then
187    
188     C- No tapering/clipping selected so ...
189     DO j=1-Oly+1,sNy+Oly-1
190     DO i=1-Olx+1,sNx+Olx-1
191    
192     dSigmaDrLtd=dSigmaDrReal(i,j)
193     if (dSigmaDrLtd.NE.0.) then
194     dRdSigmaLtd(i,j)=1./dSigmaDrLtd
195     else
196     dRdSigmaLtd(i,j)=0.
197     endif
198     SlopeX(i,j)=-SlopeX(i,j)*dRdSigmaLtd(i,j)
199     SlopeY(i,j)=-SlopeY(i,j)*dRdSigmaLtd(i,j)
200    
201     ENDDO
202     ENDDO
203    
204     endif
205    
206    
207     #endif /* ALLOW_GMREDI */
208    
209     RETURN
210     END

  ViewVC Help
Powered by ViewVC 1.1.22