/[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.7 - (show annotations) (download)
Tue Feb 15 00:21:12 2005 UTC (19 years, 3 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57m_post, checkpoint57g_pre, checkpoint57s_post, checkpoint57g_post, checkpoint57r_post, checkpoint57i_post, checkpoint57n_post, checkpoint57l_post, checkpoint57f_post, checkpoint57h_pre, checkpoint57h_post, checkpoint57e_post, checkpoint57p_post, checkpoint57q_post, eckpoint57e_pre, checkpoint57h_done, checkpoint57j_post, checkpoint57f_pre, checkpoint57o_post, checkpoint57k_post
Changes since 1.6: +3 -2 lines 
return filter tendency as output (in temporary array)
1 C $Header: /u/gcmpack/MITgcm/pkg/shap_filt/shap_filt_tracer_s2.F,v 1.6 2004/12/04 00:19: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
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 = nShapTr - nShapTrPhys
61
62 DO bj=myByLo(myThid),myByHi(myThid)
63 DO bi=myBxLo(myThid),myBxHi(myThid)
64 DO K=1,kSize
65 DO J=1-Oly,sNy+Oly
66 DO I=1-Olx,sNx+Olx
67 tmpFld(i,j,k,bi,bj)=field(i,j,k,bi,bj)
68 ENDDO
69 ENDDO
70 ENDDO
71 ENDDO
72 ENDDO
73
74
75 C ( d_xx +d_yy )^n tmpFld
76
77 C-- Computational Filter
78 DO N=1,nShapComput
79
80 IF (kSize.EQ.Nr) THEN
81 _EXCH_XYZ_R8( tmpFld, myThid )
82 ELSE
83 _EXCH_XY_R8( tmpFld, myThid )
84 ENDIF
85
86 DO bj=myByLo(myThid),myByHi(myThid)
87 DO bi=myBxLo(myThid),myBxHi(myThid)
88 DO K=1,kSize
89
90 DO J=1,sNy
91 DO I=1,sNx
92 tmpGrd(i,j) =
93 & ( tmpFld(i+1,j,k,bi,bj)-tmpFld( i ,j,k,bi,bj) )
94 & *_maskW(i+1,j,k,bi,bj)
95 & -( tmpFld( i ,j,k,bi,bj)-tmpFld(i-1,j,k,bi,bj) )
96 & *_maskW( i ,j,k,bi,bj)
97 & +( tmpFld(i,j+1,k,bi,bj)-tmpFld(i, j ,k,bi,bj) )
98 & *_maskS(i,j+1,k,bi,bj)
99 & -( tmpFld(i, j ,k,bi,bj)-tmpFld(i,j-1,k,bi,bj) )
100 & *_maskS(i, j ,k,bi,bj)
101 ENDDO
102 ENDDO
103
104 DO J=1,sNy
105 DO I=1,sNx
106 tmpFld(i,j,k,bi,bj) = -0.125*tmpGrd(i,j)
107 ENDDO
108 ENDDO
109
110 ENDDO
111 ENDDO
112 ENDDO
113 C end loop N=1,nShapComput
114 ENDDO
115
116 DO N=1,nShapTrPhys
117 C-- Physical space Filter
118
119 IF (kSize.EQ.Nr) THEN
120 _EXCH_XYZ_R8( tmpFld, myThid )
121 ELSE
122 _EXCH_XY_R8( tmpFld, myThid )
123 ENDIF
124
125 DO bj=myByLo(myThid),myByHi(myThid)
126 DO bi=myBxLo(myThid),myBxHi(myThid)
127 DO K=1,kSize
128
129 DO J=1,sNy
130 DO I=1,sNx
131 tmpGrd(i,j) =
132 & ( tmpFld(i+1,j,k,bi,bj)-tmpFld( i ,j,k,bi,bj) )
133 & *_hFacW(i+1,j,k,bi,bj)
134 & *DYG(i+1,j,bi,bj)
135 & *recip_DXC(i+1,j,bi,bj)
136 & -( tmpFld( i ,j,k,bi,bj)-tmpFld(i-1,j,k,bi,bj) )
137 & *_hFacW( i ,j,k,bi,bj)
138 & *DYG( i ,j,bi,bj)
139 & *recip_DXC( i ,j,bi,bj)
140 & +( tmpFld(i,j+1,k,bi,bj)-tmpFld(i, j ,k,bi,bj) )
141 & *_hFacS(i,j+1,k,bi,bj)
142 & *DXG(i,j+1,bi,bj)
143 & *recip_DYC(i,j+1,bi,bj)
144 & -( tmpFld(i, j ,k,bi,bj)-tmpFld(i,j-1,k,bi,bj) )
145 & *_hFacS(i, j ,k,bi,bj)
146 & *DXG(i, j ,bi,bj)
147 & *recip_DYC(i, j ,bi,bj)
148 ENDDO
149 ENDDO
150
151 IF (Shap_TrLength.EQ.0.) THEN
152 DO J=1,sNy
153 DO I=1,sNx
154 tmpFld(i,j,k,bi,bj) = -0.125*tmpGrd(i,j)
155 & *recip_hFacC(i,j,k,bi,bj)
156 ENDDO
157 ENDDO
158 ELSE
159 DO J=1,sNy
160 DO I=1,sNx
161 tmpFld(i,j,k,bi,bj) = -0.125*tmpGrd(i,j)
162 & *recip_hFacC(i,j,k,bi,bj)*recip_rA(i,j,bi,bj)
163 & *Shap_TrLength*Shap_TrLength
164 ENDDO
165 ENDDO
166 ENDIF
167
168 ENDDO
169 ENDDO
170 ENDDO
171 C end loop N=1,nShapTrPhys
172 ENDDO
173
174 C F <- [1 - (d_xx+d_yy)^n *deltaT/tau].F
175 DO bj=myByLo(myThid),myByHi(myThid)
176 DO bi=myBxLo(myThid),myBxHi(myThid)
177 DO K=1,kSize
178 DO J=1,sNy
179 DO I=1,sNx
180 field(i,j,k,bi,bj)=field(i,j,k,bi,bj)
181 & -tmpFld(i,j,k,bi,bj)*dTtracerLev(1)/Shap_Trtau
182 tmpFld(i,j,k,bi,bj)= -tmpFld(i,j,k,bi,bj)/Shap_Trtau
183 ENDDO
184 ENDDO
185 ENDDO
186 ENDDO
187 ENDDO
188
189 IF (kSize.EQ.Nr) THEN
190 _EXCH_XYZ_R8( field, myThid )
191 ELSEIF (kSize.EQ.1) THEN
192 _EXCH_XY_R8( field, myThid )
193 ELSE
194 STOP 'S/R SHAP_FILT_TRACER_S4: kSize is wrong'
195 ENDIF
196
197 ENDIF
198 #endif /* ALLOW_SHAP_FILT */
199
200 RETURN
201 END
  ViewVC Help
Powered by ViewVC 1.1.22