/[MITgcm]/MITgcm/pkg/thsice/thsice_main.F
ViewVC logotype

Contents of /MITgcm/pkg/thsice/thsice_main.F

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


Revision 1.7 - (show annotations) (download)
Sun Apr 9 17:35:30 2006 UTC (18 years, 2 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint58e_post, checkpoint58d_post
Changes since 1.6: +51 -1 lines
Starting thsice adjoint

1 C $Header: /u/gcmpack/MITgcm/pkg/thsice/thsice_main.F,v 1.6 2006/01/22 15:58:39 jmc Exp $
2 C $Name: $
3
4 #include "THSICE_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: THSICE_MAIN
8 C !INTERFACE:
9 SUBROUTINE THSICE_MAIN(
10 I myTime, myIter, myThid )
11 C !DESCRIPTION: \bv
12 C *==========================================================*
13 C | S/R THSICE_MAIN
14 C | o Therm_SeaIce main routine.
15 C | step forward Thermodynamic_SeaIce variables and modify
16 C | ocean surface forcing accordingly.
17 C *==========================================================*
18
19 C !USES:
20 IMPLICIT NONE
21
22 C === Global variables ===
23 #include "SIZE.h"
24 #include "EEPARAMS.h"
25 #include "PARAMS.h"
26 #include "GRID.h"
27 #include "SURFACE.h"
28 #include "DYNVARS.h"
29 #include "FFIELDS.h"
30 #include "THSICE_PARAMS.h"
31 #include "THSICE_VARS.h"
32 #ifdef ALLOW_AUTODIFF_TAMC
33 # include "tamc.h"
34 # include "tamc_keys.h"
35 #endif
36
37 C !INPUT/OUTPUT PARAMETERS:
38 C === Routine arguments ===
39 C myIter :: iteration counter for this thread
40 C myTime :: time counter for this thread
41 C myThid :: thread number for this instance of the routine.
42 _RL myTime
43 INTEGER myIter
44 INTEGER myThid
45 CEOP
46
47 #ifdef ALLOW_THSICE
48 C !LOCAL VARIABLES:
49 C === Local variables ===
50 INTEGER i,j
51 INTEGER bi,bj
52 INTEGER iMin, iMax
53 INTEGER jMin, jMax
54 _RL prcAtm(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
55 _RL evpAtm(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
56 _RL flxSW (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
57
58 _RL tauFac
59
60 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
61
62 IF ( stressReduction.GT. 0. _d 0 ) THEN
63 C- needs new Ice Fraction in halo region to apply wind-stress reduction
64 iMin = 1-Olx
65 iMax = sNx+Olx-1
66 jMin = 1-Oly
67 jMax = sNy+Oly-1
68 #ifdef ATMOSPHERIC_LOADING
69 ELSEIF ( useRealFreshWaterFlux ) THEN
70 C- needs sea-ice loading in part of the halo regions for grad.Phi0surf
71 C to be valid at the boundaries ( d/dx 1:sNx+1 ; d/dy 1:sNy+1 )
72 iMin = 0
73 iMax = sNx+1
74 jMin = 0
75 jMax = sNy+1
76 #endif /* ATMOSPHERIC_LOADING */
77 ELSE
78 iMin = 1
79 iMax = sNx
80 jMin = 1
81 jMax = sNy
82 ENDIF
83
84 DO bj=myByLo(myThid),myByHi(myThid)
85 DO bi=myBxLo(myThid),myBxHi(myThid)
86
87 #ifdef ALLOW_AUTODIFF_TAMC
88 act1 = bi - myBxLo(myThid)
89 max1 = myBxHi(myThid) - myBxLo(myThid) + 1
90 act2 = bj - myByLo(myThid)
91 max2 = myByHi(myThid) - myByLo(myThid) + 1
92 act3 = myThid - 1
93 max3 = nTx*nTy
94 act4 = ikey_dynamics - 1
95 iicekey = (act1 + 1) + act2*max1
96 & + act3*max1*max2
97 & + act4*max1*max2*max3
98 #endif /* ALLOW_AUTODIFF_TAMC */
99
100 C-- Mixed layer thickness: take the 1rst layer
101 #ifdef NONLIN_FRSURF
102 IF ( staggerTimeStep .AND. nonlinFreeSurf.GT.0 ) THEN
103 IF ( select_rStar.GT.0 ) THEN
104 DO j = jMin, jMax
105 DO i = iMin, iMax
106 hOceMxL(i,j,bi,bj) = drF(1)*h0FacC(i,j,1,bi,bj)
107 & *rStarFacC(i,j,bi,bj)
108 ENDDO
109 ENDDO
110 ELSE
111 DO j = jMin, jMax
112 DO i = iMin, iMax
113 IF ( ksurfC(i,j,bi,bj).EQ.1 ) THEN
114 hOceMxL(i,j,bi,bj) = drF(1)*hFac_surfC(i,j,bi,bj)
115 ELSE
116 hOceMxL(i,j,bi,bj) = drF(1)*hfacC(i,j,1,bi,bj)
117 ENDIF
118 ENDDO
119 ENDDO
120 ENDIF
121 ELSE
122 #else /* ndef NONLIN_FRSURF */
123 IF (.TRUE.) THEN
124 #endif /* NONLIN_FRSURF */
125 DO j = jMin, jMax
126 DO i = iMin, iMax
127 hOceMxL(i,j,bi,bj) = drF(1)*hfacC(i,j,1,bi,bj)
128 ENDDO
129 ENDDO
130 ENDIF
131
132 #ifdef ALLOW_AUTODIFF_TAMC
133 CADJ STORE uvel (:,:,1,bi,bj) = comlev1_bibj, key=iicekey, byte=isbyte
134 CADJ STORE vvel (:,:,1,bi,bj) = comlev1_bibj, key=iicekey, byte=isbyte
135 #endif
136
137 DO j = jMin, jMax
138 DO i = iMin, iMax
139 tOceMxL(i,j,bi,bj) = theta(i,j,1,bi,bj)
140 sOceMxL(i,j,bi,bj) = salt (i,j,1,bi,bj)
141 v2ocMxL(i,j,bi,bj) =
142 & ( uvel(i,j,1,bi,bj)*uvel(i,j,1,bi,bj)
143 & + uvel(i+1,j,1,bi,bj)*uvel(i+1,j,1,bi,bj)
144 & + vvel(i,j+1,1,bi,bj)*vvel(i,j+1,1,bi,bj)
145 & + vvel(i,j,1,bi,bj)*vvel(i,j,1,bi,bj)
146 & )*0.5 _d 0
147 prcAtm(i,j) = 0.
148 evpAtm(i,j) = 0.
149 flxSW (i,j) = 0.
150 snowPrc(i,j,bi,bj) = 0. _d 0
151 siceAlb(i,j,bi,bj) = 0. _d 0
152 ENDDO
153 ENDDO
154
155 #ifdef ALLOW_AUTODIFF_TAMC
156 CADJ STORE iceMask = comlev1, key = iicekey
157 CADJ STORE iceHeight = comlev1, key = iicekey
158 CADJ STORE snowHeight = comlev1, key = iicekey
159 CADJ STORE Tsrf = comlev1, key = iicekey
160 CADJ STORE Qice1 = comlev1, key = iicekey
161 CADJ STORE Qice2 = comlev1, key = iicekey
162 CADJ STORE snowAge = comlev1, key = iicekey
163
164 CADJ STORE sHeating = comlev1, key = iicekey
165 CADJ STORE flxCndBt = comlev1, key = iicekey
166 CADJ STORE snowPrc = comlev1, key = iicekey
167
168 CADJ STORE hOceMxL = comlev1, key = iicekey
169 CADJ STORE tOceMxL = comlev1, key = iicekey
170 CADJ STORE sOceMxL = comlev1, key = iicekey
171 CADJ STORE v2ocMxL = comlev1, key = iicekey
172
173 CADJ STORE empmr = comlev1, key = iicekey
174 CADJ STORE qnet = comlev1, key = iicekey
175 #endif
176
177 #ifdef ALLOW_BULK_FORCE
178 IF ( useBulkforce ) THEN
179 CALL THSICE_GET_PRECIP(
180 I iceMask,
181 O prcAtm, snowPrc(1-Olx,1-Oly,bi,bj), flxSW,
182 I iMin,iMax,jMin,jMax, bi,bj, myThid )
183 ENDIF
184 #endif
185
186
187 CALL THSICE_STEP_FWD(
188 I bi, bj, iMin, iMax, jMin, jMax,
189 I prcAtm,
190 U evpAtm, flxSW,
191 I myTime, myIter, myThid )
192
193 CALL THSICE_AVE(
194 I evpAtm, flxSW,
195 I bi,bj, myTime, myIter, myThid )
196
197 c ENDDO
198 c ENDDO
199
200 c IF ( .FALSE. ) THEN
201 #ifdef ALLOW_AUTODIFF_TAMC
202 CADJ STORE fu(:,:,bi,bj) = comlev1_bibj, key=iicekey, byte=isbyte
203 CADJ STORE fv(:,:,bi,bj) = comlev1_bibj, key=iicekey, byte=isbyte
204 #endif
205 IF ( stressReduction.GT. 0. _d 0 ) THEN
206 DO j = jMin, jMax
207 DO i = iMin+1,iMax
208 tauFac = stressReduction
209 & *(iceMask(i-1,j,bi,bj)+iceMask(i,j,bi,bj))*0.5 _d 0
210 fu(i,j,bi,bj) = (1. _d 0 - tauFac)*fu(i,j,bi,bj)
211 ENDDO
212 ENDDO
213 DO j = jMin+1, jMax
214 DO i = iMin, iMax
215 tauFac = stressReduction
216 & *(iceMask(i,j-1,bi,bj)+iceMask(i,j,bi,bj))*0.5 _d 0
217 fv(i,j,bi,bj) = (1. _d 0 - tauFac)*fv(i,j,bi,bj)
218 ENDDO
219 ENDDO
220 ENDIF
221
222 C-- end bi,bj loop
223 ENDDO
224 ENDDO
225
226 #ifdef ATMOSPHERIC_LOADING
227 c IF (useRealFreshWaterFlux) _EXCH_XY_RS(sIceLoad, myThid)
228 #endif
229
230 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
231 #endif /*ALLOW_THSICE*/
232
233 RETURN
234 END

  ViewVC Help
Powered by ViewVC 1.1.22