/[MITgcm]/MITgcm/pkg/seaice/seaice_advdiff.F
ViewVC logotype

Contents of /MITgcm/pkg/seaice/seaice_advdiff.F

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


Revision 1.1 - (show annotations) (download)
Thu Feb 16 10:41:48 2006 UTC (18 years, 3 months ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint58b_post
 add a few new advection schemes to seaice:
 ENUM_UPWIND_1RST, ENUM_DST2, ENUM_FLUX_LIMIT, ENUM_DST3,
 ENUM_DST3_FLUX_LIMIT
 Default is still the old one

1 C $Header: $
2 C $Name: $
3
4 #include "SEAICE_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: SEAICE_ADVDIFF
8
9 C !INTERFACE: ==========================================================
10 SUBROUTINE SEAICE_ADVDIFF(
11 I myTime, myIter, myThid )
12
13 C !DESCRIPTION: \bv
14 C /===========================================================\
15 C | SUBROUTINE SEAICE_ADVDIFF |
16 C | o driver for different advection routines |
17 C | calls an adaption of gad_advection to call different |
18 C | advection routines of pkg/generic_advdiff |
19 C \===========================================================/
20 IMPLICIT NONE
21 c \ev
22
23 C !USES: ===============================================================
24 C === Global variables ===
25 #include "SIZE.h"
26 #include "EEPARAMS.h"
27 #include "PARAMS.h"
28 #include "GRID.h"
29 #include "GAD.h"
30 #include "SEAICE_PARAMS.h"
31 #include "SEAICE.h"
32
33 #ifdef ALLOW_AUTODIFF_TAMC
34 # include "tamc.h"
35 #endif
36
37 C !INPUT PARAMETERS: ===================================================
38 C === Routine arguments ===
39 C UICE/VICE - ice velocity
40 C HEFF - scalar field to be advected
41 C HEFFM - mask for scalar field
42 C myTime - current time
43 C myIter - iteration number
44 C myThid - Thread no. that called this routine.
45 CML _RL UICE (1-OLx:sNx+OLx,1-OLy:sNy+OLy,3,nSx,nSy)
46 CML _RL VICE (1-OLx:sNx+OLx,1-OLy:sNy+OLy,3,nSx,nSy)
47 CML _RL HEFF (1-OLx:sNx+OLx,1-OLy:sNy+OLy,3,nSx,nSy)
48 CML _RL HEFFM (1-OLx:sNx+OLx,1-OLy:sNy+OLy, nSx,nSy)
49 _RL myTime
50 INTEGER myIter
51 INTEGER myThid
52 CEndOfInterface
53
54 #ifdef ALLOW_SEAICE
55 C !LOCAL VARIABLES: ====================================================
56 C === Local variables ===
57 C i,j,bi,bj - Loop counters
58 C uc/vc - current ice velocity on C-grid
59 C fld - copy of scalar field
60 C gfld - tendency of scalar field
61 INTEGER i, j, bi, bj, k3
62 LOGICAL SEAICEmultiDimAdvection
63
64 _RL uc (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
65 _RL vc (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
66 _RL fld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
67 _RL gfld (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
68 CEOP
69
70 SEAICEmultidimadvection = .TRUE.
71 IF ( SEAICEadvScheme.EQ.ENUM_CENTERED_2ND
72 & .OR.SEAICEadvScheme.EQ.ENUM_UPWIND_3RD
73 & .OR.SEAICEadvScheme.EQ.ENUM_CENTERED_4TH ) THEN
74 SEAICEmultiDimAdvection = .FALSE.
75 ENDIF
76
77
78 IF ( SEAICEmultiDimAdvection ) THEN
79 C This has to be done to comply with the time stepping in advect.F:
80 C Making sure that the following routines see the different
81 C time levels correctly
82 C At the end of the routine ADVECT,
83 C timelevel 1 is updated with advection contribution
84 C and diffusion contribution
85 C (which was computed in DIFFUS on timelevel 3)
86 C timelevel 2 is the previous timelevel 1
87 C timelevel 3 is the total diffusion tendency * deltaT
88 C (empty if no diffusion)
89
90 DO bj=myByLo(myThid),myByHi(myThid)
91 DO bi=myBxLo(myThid),myBxHi(myThid)
92 DO j=1-OLy,sNy+OLy
93 DO i=1-OLx,sNx+OLx
94 HEFF(I,J,3,bi,bj) = 0. _d 0 !HEFF(I,J,2,bi,bj)
95 HEFF(I,J,2,bi,bj) = HEFF(I,J,1,bi,bj)
96 AREA(I,J,3,bi,bj) = 0. _d 0 !AREA(I,J,2,bi,bj)
97 AREA(I,J,2,bi,bj) = AREA(I,J,1,bi,bj)
98 ENDDO
99 ENDDO
100 ENDDO
101 ENDDO
102
103 DO bj=myByLo(myThid),myByHi(myThid)
104 DO bi=myBxLo(myThid),myBxHi(myThid)
105 C average seaice velocity to c-grid
106 DO j=1-Oly,sNy+Oly-1
107 DO i=1-Olx,sNx+Olx-1
108 uc(I,J,bi,bj)=.5 _d 0*(UICE(I,J,1,bi,bj)+UICE(I,J+1,1,bi,bj))
109 vc(I,J,bi,bj)=.5 _d 0*(VICE(I,J,1,bi,bj)+VICE(I+1,J,1,bi,bj))
110 ENDDO
111 ENDDO
112 ENDDO
113 ENDDO
114 C Do we need this? I am afraid so.
115 CALL EXCH_UV_XY_RL(uc,vc,.TRUE.,myThid)
116
117 C Thickness (Volume)
118 DO bj=myByLo(myThid),myByHi(myThid)
119 DO bi=myBxLo(myThid),myBxHi(myThid)
120 DO j=1-Oly,sNy+Oly
121 DO i=1-Olx,sNx+Olx
122 fld (I,J,bi,bj) = HEFF(I,J,1,bi,bj)
123 gfld(I,J,bi,bj) = 0. _d 0
124 ENDDO
125 ENDDO
126 ENDDO
127 ENDDO
128
129 DO bj=myByLo(myThid),myByHi(myThid)
130 DO bi=myBxLo(myThid),myBxHi(myThid)
131 C advection
132 CALL SEAICE_ADVECTION(
133 I SEAICEadvScheme, GAD_HEFF,
134 I uc, vc, fld,
135 O gfld,
136 I bi, bj, myTime, myIter, myThid)
137 C now do the "explicit" time step
138 DO j=1,sNy
139 DO i=1,sNx
140 HEFF(I,J,1,bi,bj) = HEFFM(I,J,bi,bj) * (
141 & HEFF(I,J,1,bi,bj) + SEAICE_deltaTtherm * gFld(I,J,bi,bj)
142 & )
143 ENDDO
144 ENDDO
145 ENDDO
146 ENDDO
147 IF ( DIFF1 .GT. 0. _d 0 ) THEN
148 C Do we need this?
149 CALL SEAICE_EXCH( HEFF, myThid )
150 C diffusion
151 CALL SEAICE_DIFFUSION(
152 U HEFF,
153 I HEFFM, SEAICE_deltaTtherm, myTime, myIter, myThid )
154 ENDIF
155
156 C Fractional Area
157 DO bj=myByLo(myThid),myByHi(myThid)
158 DO bi=myBxLo(myThid),myBxHi(myThid)
159 DO j=1-Oly,sNy+Oly
160 DO i=1-Olx,sNx+Olx
161 fld (I,J,bi,bj) = AREA(I,J,1,bi,bj)
162 gfld(I,J,bi,bj) = 0. _d 0
163 ENDDO
164 ENDDO
165 ENDDO
166 ENDDO
167
168 DO bj=myByLo(myThid),myByHi(myThid)
169 DO bi=myBxLo(myThid),myBxHi(myThid)
170 C advection
171 CALL SEAICE_ADVECTION(
172 I SEAICEadvScheme, GAD_AREA,
173 I uc, vc, fld,
174 O gfld,
175 I bi, bj, myTime, myIter, myThid)
176 C now do the "explicit" time step
177 DO j=1,sNy
178 DO i=1,sNx
179 AREA(I,J,1,bi,bj) = HEFFM(I,J,bi,bj) * (
180 & AREA(I,J,1,bi,bj) + SEAICE_deltaTtherm * gFld(I,J,bi,bj)
181 & )
182 ENDDO
183 ENDDO
184 ENDDO
185 ENDDO
186
187 IF ( DIFF1 .GT. 0. _d 0 ) THEN
188 C Do we need this? Probably not, but it is done in ADVECT
189 CALL SEAICE_EXCH( AREA, myThid )
190 C diffusion
191 CALL SEAICE_DIFFUSION(
192 U AREA,
193 I HEFFM, SEAICE_deltaTtherm, myTime, myIter, myThid )
194 ENDIF
195
196 ELSE
197 C if not multiDimAdvection
198 CALL ADVECT( UICE, VICE, HEFF, HEFFM, myThid )
199 CALL ADVECT( UICE, VICE, AREA, HEFFM, myThid )
200 ENDIF
201
202 #endif /* ALLOW_SEAICE */
203
204 RETURN
205 END

  ViewVC Help
Powered by ViewVC 1.1.22