/[MITgcm]/MITgcm_contrib/dgoldberg/streamice/streamice_adv_flux_fl_y.F
ViewVC logotype

Annotation of /MITgcm_contrib/dgoldberg/streamice/streamice_adv_flux_fl_y.F

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


Revision 1.1 - (hide annotations) (download)
Sat Jun 8 22:15:33 2013 UTC (12 years, 1 month ago) by dgoldberg
Branch: MAIN
new advected scalar; new advection scheme for thickness update; corresponding TAF directives

1 dgoldberg 1.1 C $Header: /u/gcmpack/MITgcm_contrib/dgoldberg/streamice/streamice_advect_thickness_x.F,v 1.4 2012/09/18 17:06:48 dgoldberg Exp $
2     C $Name: $
3    
4     #include "STREAMICE_OPTIONS.h"
5    
6     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7    
8     SUBROUTINE STREAMICE_ADV_FLUX_FL_Y ( myThid ,
9     I VADV ,
10     I TRAC ,
11     I BC_FACEMASK,
12     I BC_YVALUES,
13     O YFLUX,
14     I time_step )
15    
16     IMPLICIT NONE
17    
18     C O hflux_x ! flux per unit width across face
19     C O h
20     C I time_step
21    
22     C === Global variables ===
23     #include "SIZE.h"
24     #include "GRID.h"
25     #include "EEPARAMS.h"
26     #include "PARAMS.h"
27     #include "STREAMICE.h"
28     !#include "GAD_FLUX_LIMITER.h"
29    
30     #ifdef ALLOW_AUTODIFF_TAMC
31     # include "tamc.h"
32     #endif
33    
34    
35     INTEGER myThid
36     _RL VADV (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
37     _RL TRAC (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
38     _RL BC_FACEMASK (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
39     _RL BC_YVALUES (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
40     _RL YFLUX (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
41     _RL time_step
42    
43     #ifdef ALLOW_STREAMICE
44    
45     C LOCAL VARIABLES
46    
47     INTEGER i, j, bi, bj, Gi, Gj, k
48     _RL vface, phi, cfl, Cr, rdenom, d0, d1, psi
49     _RL stencil (-1:1)
50     LOGICAL H0_valid(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
51     ! there are valid cells to calculate a
52     ! slope-limited 2nd order flux
53     _RL SLOPE_LIMITER
54     ! _RL total_vol_out
55     external SLOPE_LIMITER
56    
57     ! total_vol_out = 0.0
58    
59    
60     DO bj=myByLo(myThid),myByHi(myThid)
61     DO bi=myBxLo(myThid),myBxHi(myThid)
62     DO j=1-oly,sNy+oly
63     DO i=1-olx,sNx+olx
64     H0_valid(i,j,bi,bj)=.false.
65     ENDDO
66     ENDDO
67     ENDDO
68     ENDDO
69    
70     DO bj=myByLo(myThid),myByHi(myThid)
71     DO bi=myBxLo(myThid),myBxHi(myThid)
72     DO j=1-1,sNy+2
73     Gj = (myYGlobalLo-1)+(bj-1)*sNy+j
74     DO i=1-1,sNx+1
75     C THESE ARRAY BOUNDS INSURE THAT AFTER THIS STEP,
76     C VALUES WILL BE RELIABLE 1 GRID CELLS OUT IN THE
77     C X DIRECTION AND 1 CELLS OUT IN THE Y DIR
78     Gi = (myXGlobalLo-1)+(bi-1)*sNx+i
79     IF (((Gi .ge. 1) .and. (Gi .le. Nx)).or.
80     & STREAMICE_EW_PERIODIC ) THEN
81     IF ((STREAMICE_hmask(i,j,bi,bj).eq.1.0) .or.
82     & ((STREAMICE_hmask(i,j-1,bi,bj).eq.1.0) .and.
83     & (STREAMICE_hmask(i,j,bi,bj).ne.1.0))) THEN
84    
85    
86     vface = VADV(i,j,bi,bj)
87     cfl = ABS(vface) * time_step * recip_dyC(i,j,bi,bj)
88    
89     ! IF (BC_FACEMASK(i,j,bi,bj).eq.3.0) THEN
90     ! YFLUX (i,j,bi,bj) = BC_YVALUES(i,j,bi,bj) * vface
91    
92     IF (BC_FACEMASK(i,j,bi,bj).eq.3.0 .and.
93     & vface.gt.0 .and.
94     & STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN
95     YFLUX (i,j,bi,bj) = BC_YVALUES(i,j,bi,bj) * vface
96     ELSEIF
97     & (BC_FACEMASK(i,j,bi,bj).eq.3.0 .and.
98     & vface.le.0 .and.
99     & STREAMICE_hmask(i,j-1,bi,bj).eq.1.0) THEN
100     YFLUX (i,j,bi,bj) = BC_YVALUES(i,j,bi,bj) * vface
101    
102    
103     ELSE
104    
105     IF (vface .gt. 0. _d 0) THEN
106     DO k=-1,1
107     stencil (k) = TRAC(i,j+k-1,bi,bj)
108     ENDDO
109     IF ((STREAMICE_hmask(i,j,bi,bj).eq.1.0) .and.
110     & (STREAMICE_hmask(i,j-2,bi,bj).eq.1.0))
111     & H0_valid(i,j,bi,bj)=.true.
112    
113     IF (((Gj.eq.1).and.(STREAMICE_hmask(i,j-1,bi,bj).eq.3.0))
114     & .and.(.not.STREAMICE_NS_PERIODIC))
115     & THEN ! we are at western bdry and there is a thick. bdry cond
116    
117     YFLUX (i,j,bi,bj) = TRAC(i,j-1,bi,bj) * vface
118    
119     ELSEIF (H0_valid(i,j,bi,bj)) THEN
120    
121    
122     rdenom = (stencil(1)-stencil(0))
123     IF (rdenom .ne. 0.) THEN
124     Cr = (stencil(0)-stencil(-1))/rdenom
125     ELSE
126     Cr = 1.E20 * (stencil(0)-stencil(-1))
127     ENDIF
128    
129     IF (STREAMICE_ADV_SCHEME.ne.'DST3') THEN
130     ! phi = SLOPE_LIMITER(stencil(0)-stencil(-1),
131     ! & stencil(1)-stencil(0))
132     phi = SLOPE_LIMITER (Cr)
133     ELSE
134     d0 = (2.-cfl)*(1.-cfl)/6.0
135     d1 = (1.-cfl**2)/6.0
136     psi = d0+d1*Cr
137     phi = MAX(0. _d 0,MIN(MIN(1. _d 0,psi),
138     & Cr*(1. _d 0 -CFL)/(CFL+1. _d -20) ))
139     ENDIF
140    
141     IF (STREAMICE_ADV_SCHEME.ne.'DST3') THEN
142     YFLUX (i,j,bi,bj) = vface *
143     & (stencil(0) + phi * .5 * (1.0-cfl) *
144     & (stencil(1)-stencil(0)))
145     ELSE
146     YFLUX (i,j,bi,bj) = vface *
147     & (stencil(0) + phi *
148     & (stencil(1)-stencil(0)))
149     ENDIF
150    
151     ELSE ! one of the two cells needed for a HO scheme is missing, use FO scheme
152    
153     YFLUX (i,j,bi,bj) = vface * stencil(0)
154    
155     ENDIF
156    
157     ELSE ! vface <= 0
158    
159     DO k=-1,1
160     stencil (k) = TRAC(i,j-k,bi,bj)
161     ENDDO
162     IF ((STREAMICE_hmask(i-1,j,bi,bj).eq.1.0) .and.
163     & (STREAMICE_hmask(i+1,j,bi,bj).eq.1.0))
164     & H0_valid(i,j,bi,bj)=.true.
165    
166     IF (((Gj.eq.Ny).and.(STREAMICE_hmask(i,j+1,bi,bj).eq.3.0))
167     & .and.(.not.STREAMICE_NS_PERIODIC))
168     & THEN ! we are at western bdry and there is a thick. bdry cond
169    
170     YFLUX (i,j,bi,bj) = TRAC(i,j+1,bi,bj) * vface
171    
172     ELSEIF (H0_valid(i,j,bi,bj)) THEN
173    
174     rdenom = (stencil(1)-stencil(0))
175     IF (rdenom .ne. 0.) THEN
176     Cr = (stencil(0)-stencil(-1))/rdenom
177     ELSE
178     Cr = 1.E20 * (stencil(0)-stencil(-1))
179     ENDIF
180    
181     IF (STREAMICE_ADV_SCHEME.ne.'DST3') THEN
182     ! phi = SLOPE_LIMITER(stencil(0)-stencil(-1),
183     ! & stencil(1)-stencil(0))
184     phi = SLOPE_LIMITER (Cr)
185     ELSE
186     d0 = (2.-cfl)*(1.-cfl)/6.0
187     d1 = (1.-cfl**2)/6.0
188     psi = d0+d1*Cr
189     phi = MAX(0. _d 0,MIN(MIN(1. _d 0,psi),
190     & Cr*(1. _d 0 -CFL)/(CFL+1. _d -20) ))
191     ENDIF
192    
193     IF (STREAMICE_ADV_SCHEME.ne.'DST3') THEN
194     YFLUX (i,j,bi,bj) = vface *
195     & (stencil(0) + phi * .5 * (1.0-cfl) *
196     & (stencil(1)-stencil(0)))
197     ELSE
198     YFLUX (i,j,bi,bj) = vface *
199     & (stencil(0) + phi *
200     & (stencil(1)-stencil(0)))
201     ENDIF
202    
203     ELSE ! one of the two cells needed for a HO scheme is missing, use FO scheme
204    
205     yflux (i,j,bi,bj) = vface * stencil(0)
206    
207     ENDIF
208    
209     ENDIF
210    
211     ENDIF
212    
213     ENDIF
214     ENDIF
215     ENDDO
216     ENDDO
217     ENDDO
218     ENDDO
219    
220     !C X-FLUXES AT CELL BOUNDARIES CALCULATED; NOW TAKE FLUX DIVERGENCE TO INCREMENT THICKNESS
221     !
222     ! DO bj=myByLo(myThid),myByHi(myThid)
223     ! DO bi=myBxLo(myThid),myBxHi(myThid)
224     ! DO j=1-3,sNy+3
225     ! Gj = (myYGlobalLo-1)+(bj-1)*sNy+j
226     ! IF ((Gj .ge. 1) .and. (Gj .le. Ny)) THEN
227     ! DO i=1-2,sNx+2
228     ! IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN
229     ! h(i,j,bi,bj) = h(i,j,bi,bj) - time_step *
230     ! & (hflux_x(i+1,j,bi,bj)*dyG(i+1,j,bi,bj) -
231     ! & hflux_x(i,j,bi,bj)*dyG(i,j,bi,bj)) *
232     ! & recip_rA (i,j,bi,bj)
233     ! ENDIF
234     ! ENDDO
235     ! ENDIF
236     ! ENDDO
237     ! ENDDO
238     ! ENDDO
239    
240    
241    
242     #endif
243     RETURN
244     END SUBROUTINE STREAMICE_ADV_FLUX_FL_Y
245    

  ViewVC Help
Powered by ViewVC 1.1.22