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_ADVECT_THICKNESS_X ( myThid , |
9 |
O hflux_x , |
10 |
O h , |
11 |
I time_step ) |
12 |
|
13 |
IMPLICIT NONE |
14 |
|
15 |
C O hflux_x ! flux per unit width across face |
16 |
C O h |
17 |
C I time_step |
18 |
|
19 |
C === Global variables === |
20 |
#include "SIZE.h" |
21 |
#include "GRID.h" |
22 |
#include "EEPARAMS.h" |
23 |
#include "PARAMS.h" |
24 |
#include "STREAMICE.h" |
25 |
|
26 |
INTEGER myThid |
27 |
_RL hflux_x (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy) |
28 |
_RL h (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy) |
29 |
_RL time_step |
30 |
|
31 |
#ifdef ALLOW_STREAMICE |
32 |
|
33 |
C LOCAL VARIABLES |
34 |
|
35 |
INTEGER i, j, bi, bj, Gi, Gj, k |
36 |
_RL uface, phi |
37 |
_RL stencil (-1:1) |
38 |
LOGICAL H0_valid(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy) |
39 |
! there are valid cells to calculate a |
40 |
! slope-limited 2nd order flux |
41 |
_RL SLOPE_LIMITER |
42 |
! _RL total_vol_out |
43 |
external SLOPE_LIMITER |
44 |
|
45 |
! total_vol_out = 0.0 |
46 |
|
47 |
DO bj=myByLo(myThid),myByHi(myThid) |
48 |
DO bi=myBxLo(myThid),myBxHi(myThid) |
49 |
DO j=1-oly,sNy+oly |
50 |
DO i=1-olx,sNx+olx |
51 |
H0_valid(i,j,bi,bj)=.false. |
52 |
ENDDO |
53 |
ENDDO |
54 |
ENDDO |
55 |
ENDDO |
56 |
|
57 |
DO bj=myByLo(myThid),myByHi(myThid) |
58 |
DO bi=myBxLo(myThid),myBxHi(myThid) |
59 |
DO j=1-2,sNy+2 |
60 |
Gj = (myYGlobalLo-1)+(bj-1)*sNy+j |
61 |
IF ((Gj .ge. 1) .and. (Gj .le. Ny)) THEN |
62 |
DO i=1-1,sNx+2 |
63 |
C THESE ARRAY BOUNDS INSURE THAT AFTER THIS STEP, |
64 |
C VALUES WILL BE RELIABLE 2 GRID CELLS OUT IN THE |
65 |
C X DIRECTION AND 3 CELLS OUT IN THE Y DIR |
66 |
IF ((STREAMICE_hmask(i,j,bi,bj).eq.1.0) .or. |
67 |
& ((STREAMICE_hmask(i-1,j,bi,bj).eq.1.0) .and. |
68 |
& (STREAMICE_hmask(i,j,bi,bj).ne.1.0))) THEN |
69 |
|
70 |
Gi = (myXGlobalLo-1)+(bi-1)*sNx+i |
71 |
|
72 |
IF (STREAMICE_ufacemask(i,j,bi,bj).eq.4.0) THEN |
73 |
hflux_x (i,j,bi,bj) = u_flux_bdry_SI (i,j,bi,bj) |
74 |
ELSE |
75 |
|
76 |
uface = .5 * |
77 |
& (U_streamice(i,j,bi,bj)+U_streamice(i,j+1,bi,bj)) |
78 |
|
79 |
|
80 |
IF (uface .gt. 0. _d 0) THEN |
81 |
DO k=-1,1 |
82 |
stencil (k) = h(i+k-1,j,bi,bj) |
83 |
ENDDO |
84 |
IF ((STREAMICE_hmask(i,j,bi,bj).eq.1.0) .and. |
85 |
& (STREAMICE_hmask(i-2,j,bi,bj).eq.1.0)) |
86 |
& H0_valid(i,j,bi,bj)=.true. |
87 |
|
88 |
IF ((STREAMICE_hmask(i-1,j,bi,bj).eq.3.0)) |
89 |
& THEN ! we are at western bdry and there is a thick. bdry cond |
90 |
hflux_x (i,j,bi,bj) = h(i-1,j,bi,bj) * uface |
91 |
! PRINT *, "BOUNDARY FLUX RIGHT", hflux_x (i,j,bi,bj), |
92 |
! & h(i-1,j,bi,bj),uface,i,j |
93 |
ELSEIF (H0_valid(i,j,bi,bj)) THEN |
94 |
phi = SLOPE_LIMITER ( |
95 |
& stencil(0)-stencil(-1), |
96 |
& stencil(1)-stencil(0)) |
97 |
hflux_x (i,j,bi,bj) = uface * |
98 |
& (stencil(0) - phi * .5 * (stencil(0)-stencil(1))) |
99 |
ELSE ! one of the two cells needed for a HO scheme is missing, use FO scheme |
100 |
hflux_x (i,j,bi,bj) = uface * stencil(0) |
101 |
ENDIF |
102 |
|
103 |
ELSE ! uface <= 0 |
104 |
|
105 |
DO k=-1,1 |
106 |
stencil (k) = h(i-k,j,bi,bj) |
107 |
ENDDO |
108 |
IF ((STREAMICE_hmask(i-1,j,bi,bj).eq.1.0) .and. |
109 |
& (STREAMICE_hmask(i+1,j,bi,bj).eq.1.0)) |
110 |
& H0_valid(i,j,bi,bj)=.true. |
111 |
|
112 |
IF ((STREAMICE_hmask(i,j,bi,bj).eq.3.0)) |
113 |
& THEN ! we are at western bdry and there is a thick. bdry cond |
114 |
hflux_x (i,j,bi,bj) = h(i,j,bi,bj) * uface |
115 |
! PRINT *, "BOUNDARY FLUX LEFT", hflux_x(i,j,bi,bj), |
116 |
! & h(i,j,bi,bj),uface |
117 |
ELSEIF (H0_valid(i,j,bi,bj)) THEN |
118 |
phi = SLOPE_LIMITER ( |
119 |
& stencil(0)-stencil(-1), |
120 |
& stencil(1)-stencil(0)) |
121 |
hflux_x (i,j,bi,bj) = uface * |
122 |
& (stencil(0) - phi * .5 * (stencil(0)-stencil(1))) |
123 |
ELSE ! one of the two cells needed for a HO scheme is missing, use FO scheme |
124 |
hflux_x (i,j,bi,bj) = uface * stencil(0) |
125 |
ENDIF |
126 |
|
127 |
ENDIF |
128 |
|
129 |
ENDIF |
130 |
|
131 |
! if (streamice_ufacemask(i,j,bi,bj).eq.2.0) THEN |
132 |
! total_vol_out = total_vol_out + |
133 |
! & hflux_x (i,j,bi,bj) * time_step |
134 |
! ENDIF |
135 |
|
136 |
ENDIF |
137 |
ENDDO |
138 |
ENDIF |
139 |
ENDDO |
140 |
ENDDO |
141 |
ENDDO |
142 |
|
143 |
C X-FLUXES AT CELL BOUNDARIES CALCULATED; NOW TAKE FLUX DIVERGENCE TO INCREMENT THICKNESS |
144 |
|
145 |
DO bj=myByLo(myThid),myByHi(myThid) |
146 |
DO bi=myBxLo(myThid),myBxHi(myThid) |
147 |
DO j=1-3,sNy+3 |
148 |
Gj = (myYGlobalLo-1)+(bj-1)*sNy+j |
149 |
IF ((Gj .ge. 1) .and. (Gj .le. Ny)) THEN |
150 |
DO i=1-2,sNx+2 |
151 |
IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN |
152 |
h(i,j,bi,bj) = h(i,j,bi,bj) - time_step * |
153 |
& (hflux_x(i+1,j,bi,bj)*dyG(i+1,j,bi,bj) - |
154 |
& hflux_x(i,j,bi,bj)*dyG(i,j,bi,bj)) * |
155 |
& recip_rA (i,j,bi,bj) |
156 |
ENDIF |
157 |
ENDDO |
158 |
ENDIF |
159 |
ENDDO |
160 |
ENDDO |
161 |
ENDDO |
162 |
|
163 |
|
164 |
|
165 |
#endif |
166 |
RETURN |
167 |
END SUBROUTINE STREAMICE_ADVECT_THICKNESS_X |
168 |
|