/[MITgcm]/MITgcm/pkg/shap_filt/shap_filt_tracer_s2.F
ViewVC logotype

Contents of /MITgcm/pkg/shap_filt/shap_filt_tracer_s2.F

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


Revision 1.8 - (show annotations) (download)
Tue Sep 27 22:11:06 2005 UTC (18 years, 8 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57t_post
Changes since 1.7: +13 -12 lines
use shap-funct S2 with nShap_Phys=nShap instead of shap-funct S2G
 and get the same results.

1 C $Header: /u/gcmpack/MITgcm/pkg/shap_filt/shap_filt_tracer_s2.F,v 1.7 2005/02/15 00:21:12 jmc Exp $
2 C $Name: $
3
4 #include "SHAP_FILT_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: SHAP_FILT_TRACER_S2
8 C !INTERFACE:
9 SUBROUTINE SHAP_FILT_TRACER_S2(
10 U field, tmpFld,
11 I nShapTr, kSize, myTime, myThid )
12 C !DESCRIPTION: \bv
13 C *==========================================================*
14 C | S/R SHAP_FILT_TRACER_S2
15 C | o Applies Shapiro filter to 2D field (cell center).
16 C | o use filtering function "S2" = [1 - (d_xx+d_yy)^n]
17 C | o Options for computational filter (no grid spacing)
18 C | or physical space filter (with grid spacing) or both.
19 C *==========================================================*
20 C \ev
21
22 C !USES:
23 IMPLICIT NONE
24
25 C == Global variables ===
26 #include "SIZE.h"
27 #include "EEPARAMS.h"
28 #include "PARAMS.h"
29 #include "GRID.h"
30 #include "SHAP_FILT.h"
31
32 C !INPUT/OUTPUT PARAMETERS:
33 C == Routine arguments
34 C field :: cell-centered 2D field on which filter applies
35 C tmpFld :: working temporary array
36 C nShapTr :: (total) power of the filter for this tracer
37 C kSize :: length of 3rd Dim : either =1 (2D field) or =Nr (3D field)
38 C myTime :: Current time in simulation
39 C myThid :: Thread number for this instance of SHAP_FILT_TRACER_S2
40 INTEGER nShapTr, kSize
41 _RL field(1-OLx:sNx+OLx,1-OLy:sNy+OLy,kSize,nSx,nSy)
42 _RL tmpFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,kSize,nSx,nSy)
43 _RL myTime
44 INTEGER myThid
45
46 #ifdef ALLOW_SHAP_FILT
47
48 C !LOCAL VARIABLES:
49 C == Local variables ==
50 INTEGER nShapComput, nShapPhysic
51 INTEGER bi,bj,K,I,J,N
52 _RL tmpGrd(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
53 CEOP
54
55 IF (nShapTr.GT.0) THEN
56 C-------
57 C Apply computational filter ^(nShap-nShapPhys) without grid factor
58 C then apply Physical filter ^nShapPhys with grid factors
59 C-------
60 nShapComput = MAX( 0, nShapTr - nShapTrPhys )
61 nShapPhysic = nShapTr - nShapComput
62
63 DO bj=myByLo(myThid),myByHi(myThid)
64 DO bi=myBxLo(myThid),myBxHi(myThid)
65 DO K=1,kSize
66 DO J=1-Oly,sNy+Oly
67 DO I=1-Olx,sNx+Olx
68 tmpFld(i,j,k,bi,bj)=field(i,j,k,bi,bj)
69 ENDDO
70 ENDDO
71 ENDDO
72 ENDDO
73 ENDDO
74
75
76 C ( d_xx +d_yy )^n tmpFld
77
78 C-- Computational Filter
79 DO N=1,nShapComput
80
81 IF (kSize.EQ.Nr) THEN
82 _EXCH_XYZ_R8( tmpFld, myThid )
83 ELSE
84 _EXCH_XY_R8( tmpFld, myThid )
85 ENDIF
86
87 DO bj=myByLo(myThid),myByHi(myThid)
88 DO bi=myBxLo(myThid),myBxHi(myThid)
89 DO K=1,kSize
90
91 DO J=1,sNy
92 DO I=1,sNx
93 tmpGrd(i,j) =
94 & ( tmpFld(i+1,j,k,bi,bj)-tmpFld( i ,j,k,bi,bj) )
95 & *_maskW(i+1,j,k,bi,bj)
96 & -( tmpFld( i ,j,k,bi,bj)-tmpFld(i-1,j,k,bi,bj) )
97 & *_maskW( i ,j,k,bi,bj)
98 & +( tmpFld(i,j+1,k,bi,bj)-tmpFld(i, j ,k,bi,bj) )
99 & *_maskS(i,j+1,k,bi,bj)
100 & -( tmpFld(i, j ,k,bi,bj)-tmpFld(i,j-1,k,bi,bj) )
101 & *_maskS(i, j ,k,bi,bj)
102 ENDDO
103 ENDDO
104
105 DO J=1,sNy
106 DO I=1,sNx
107 tmpFld(i,j,k,bi,bj) = -0.125*tmpGrd(i,j)
108 ENDDO
109 ENDDO
110
111 ENDDO
112 ENDDO
113 ENDDO
114 C end loop N=1,nShapComput
115 ENDDO
116
117 DO N=1,nShapPhysic
118 C-- Physical space Filter
119
120 IF (kSize.EQ.Nr) THEN
121 _EXCH_XYZ_R8( tmpFld, myThid )
122 ELSE
123 _EXCH_XY_R8( tmpFld, myThid )
124 ENDIF
125
126 DO bj=myByLo(myThid),myByHi(myThid)
127 DO bi=myBxLo(myThid),myBxHi(myThid)
128 DO K=1,kSize
129
130 DO J=1,sNy
131 DO I=1,sNx
132 tmpGrd(i,j) =
133 & ( tmpFld(i+1,j,k,bi,bj)-tmpFld( i ,j,k,bi,bj) )
134 & *_hFacW(i+1,j,k,bi,bj)
135 & *DYG(i+1,j,bi,bj)
136 & *recip_DXC(i+1,j,bi,bj)
137 & -( tmpFld( i ,j,k,bi,bj)-tmpFld(i-1,j,k,bi,bj) )
138 & *_hFacW( i ,j,k,bi,bj)
139 & *DYG( i ,j,bi,bj)
140 & *recip_DXC( i ,j,bi,bj)
141 & +( tmpFld(i,j+1,k,bi,bj)-tmpFld(i, j ,k,bi,bj) )
142 & *_hFacS(i,j+1,k,bi,bj)
143 & *DXG(i,j+1,bi,bj)
144 & *recip_DYC(i,j+1,bi,bj)
145 & -( tmpFld(i, j ,k,bi,bj)-tmpFld(i,j-1,k,bi,bj) )
146 & *_hFacS(i, j ,k,bi,bj)
147 & *DXG(i, j ,bi,bj)
148 & *recip_DYC(i, j ,bi,bj)
149 ENDDO
150 ENDDO
151
152 IF (Shap_TrLength.LE.0.) THEN
153 DO J=1,sNy
154 DO I=1,sNx
155 tmpFld(i,j,k,bi,bj) = -0.125*tmpGrd(i,j)
156 & *recip_hFacC(i,j,k,bi,bj)
157 ENDDO
158 ENDDO
159 ELSE
160 DO J=1,sNy
161 DO I=1,sNx
162 tmpFld(i,j,k,bi,bj) = -0.125*tmpGrd(i,j)
163 & *recip_hFacC(i,j,k,bi,bj)*recip_rA(i,j,bi,bj)
164 & *Shap_TrLength*Shap_TrLength
165 ENDDO
166 ENDDO
167 ENDIF
168
169 ENDDO
170 ENDDO
171 ENDDO
172 C end loop N=1,nShapTrPhys
173 ENDDO
174
175 C F <- [1 - (d_xx+d_yy)^n *deltaT/tau].F
176 DO bj=myByLo(myThid),myByHi(myThid)
177 DO bi=myBxLo(myThid),myBxHi(myThid)
178 DO K=1,kSize
179 DO J=1,sNy
180 DO I=1,sNx
181 field(i,j,k,bi,bj)=field(i,j,k,bi,bj)
182 & -tmpFld(i,j,k,bi,bj)*dTtracerLev(1)/Shap_Trtau
183 tmpFld(i,j,k,bi,bj)= -tmpFld(i,j,k,bi,bj)/Shap_Trtau
184 ENDDO
185 ENDDO
186 ENDDO
187 ENDDO
188 ENDDO
189
190 IF (kSize.EQ.Nr) THEN
191 _EXCH_XYZ_R8( field, myThid )
192 ELSEIF (kSize.EQ.1) THEN
193 _EXCH_XY_R8( field, myThid )
194 ELSE
195 STOP 'S/R SHAP_FILT_TRACER_S4: kSize is wrong'
196 ENDIF
197
198 ENDIF
199 #endif /* ALLOW_SHAP_FILT */
200
201 RETURN
202 END

  ViewVC Help
Powered by ViewVC 1.1.22