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

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

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

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