/[MITgcm]/MITgcm_contrib/lab_sea_test/budget.F
ViewVC logotype

Contents of /MITgcm_contrib/lab_sea_test/budget.F

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


Revision 1.1 - (show annotations) (download)
Mon Jul 12 01:00:20 2004 UTC (19 years, 9 months ago) by dimitri
Branch: MAIN
CVS Tags: HEAD
added my_min_max for pkg/seaice routines

1 C $Header: /u/gcmpack/MITgcm/pkg/seaice/budget.F,v 1.11 2003/10/09 04:19:20 edhill Exp $
2 C $Name: $
3
4 #include "SEAICE_OPTIONS.h"
5
6 CStartOfInterface
7 SUBROUTINE BUDGET(UG, TICE, HICE1, FICE1, KOPEN, bi, bj)
8 C /==========================================================\
9 C | SUBROUTINE budget |
10 C | o Calculate ice growth rate |
11 C | see Hibler, MWR, 108, 1943-1973, 1980 |
12 C |==========================================================|
13 C \==========================================================/
14 IMPLICIT NONE
15
16 C === Global variables ===
17 #include "SIZE.h"
18 #include "EEPARAMS.h"
19 #include "FFIELDS.h"
20 #include "SEAICE_PARAMS.h"
21 #include "SEAICE_FFIELDS.h"
22
23 C Subset of variables from SEAICE.h
24 _RL HEFF (1-OLx:sNx+OLx,1-OLy:sNy+OLy,3,nSx,nSy)
25 _RL HSNOW (1-OLx:sNx+OLx,1-OLy:sNy+OLy, nSx,nSy)
26 _RL QNETO (1-OLx:sNx+OLx,1-OLy:sNy+OLy, nSx,nSy)
27 _RL QNETI (1-OLx:sNx+OLx,1-OLy:sNy+OLy, nSx,nSy)
28 _RL QSWO (1-OLx:sNx+OLx,1-OLy:sNy+OLy, nSx,nSy)
29 _RL QSWI (1-OLx:sNx+OLx,1-OLy:sNy+OLy, nSx,nSy)
30 COMMON/TRANS/HEFF,HSNOW
31 COMMON/QFLUX/QNETO,QNETI,QSWO,QSWI
32
33 C === Routine arguments ===
34 _RL UG (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
35 _RL TICE (1-OLx:sNx+OLx, 1-OLy:sNy+OLy, nSx,nSy)
36 _RL HICE1 (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
37 _RL FICE1 (1-OLx:sNx+OLx, 1-OLy:sNy+OLy, nSx,nSy)
38 INTEGER KOPEN
39 INTEGER bi, bj
40 CEndOfInterface
41
42 #ifdef ALLOW_SEAICE
43
44 C === Local variables ===
45 C i,j,k,bi,bj - Loop counters
46
47 INTEGER i, j
48 INTEGER ITER
49 _RL QS1, C1, C2, C3, C4, C5, TB, D1, D1W, D1I, D3
50 _RL TMELT, TMELTP, XKI, XKS, HCUT, ASNOW, XIO
51
52 _RL HICE (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
53 _RL ALB (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
54 _RL A1 (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
55 _RL A2 (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
56 _RL A3 (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
57 _RL B (1-OLx:sNx+OLx, 1-OLy:sNy+OLy)
58
59 _RL mymin_R8, mymax_R8
60 external mymin_R8, mymax_R8
61
62 C IF KOPEN LT 0, THEN DO OPEN WATER BUDGET
63 C NOW DEFINE ASSORTED CONSTANTS
64 C SATURATION VAPOR PRESSURE CONSTANT
65 QS1=0.622 _d +00/1013.0 _d +00
66 C MAYKUTS CONSTANTS FOR SAT. VAP. PRESSURE TEMP. POLYNOMIAL
67 C1=2.7798202 _d -06
68 C2=-2.6913393 _d -03
69 C3=0.97920849 _d +00
70 C4=-158.63779 _d +00
71 C5=9653.1925 _d +00
72 C FREEZING TEMPERATURE OF SEAWATER
73 TB=271.2 _d +00
74 C SENSIBLE HEAT CONSTANT
75 D1=SEAICE_sensHeat
76 C WATER LATENT HEAT CONSTANT
77 D1W=SEAICE_latentWater
78 C ICE LATENT HEAT CONSTANT
79 D1I=SEAICE_latentIce
80 C STEFAN BOLTZMAN CONSTANT TIMES 0.97 EMISSIVITY
81 D3=SEAICE_emissivity
82 C MELTING TEMPERATURE OF ICE
83 TMELT=273.16 _d +00
84 TMELTP=273.159 _d +00
85 C ICE CONDUCTIVITY
86 XKI=SEAICE_iceConduct
87 C SNOW CONDUCTIVITY
88 XKS=SEAICE_snowConduct
89 C CUTOFF SNOW THICKNESS
90 HCUT=SEAICE_snowThick
91 C PENETRATION SHORTWAVE RADIATION FACTOR
92 XIO=SEAICE_shortwave
93
94 DO J=1,sNy
95 DO I=1,sNx
96 TICE(I,J,bi,bj)=MYMIN_R8(273.16 _d 0+MAX_TICE,TICE(I,J,bi,bj))
97 ATEMP(I,J,bi,bj)=MYMAX_R8(273.16 _d 0+MIN_ATEMP,ATEMP(I,J,bi,bj))
98 LWDOWN(I,J,bi,bj)=MYMAX_R8(MIN_LWDOWN,LWDOWN(I,J,bi,bj))
99 ENDDO
100 ENDDO
101
102 C NOW DECIDE IF OPEN WATER OR ICE
103 IF(KOPEN.LE.0) THEN
104
105 C NOW DETERMINE OPEN WATER HEAT BUD. ASSUMING TICE=WATER TEMP.
106 C WATER ALBEDO IS ASSUMED TO BE THE CONSTANT SEAICE_waterAlbedo
107 DO J=1,sNy
108 DO I=1,sNx
109 #ifdef SEAICE_EXTERNAL_FLUXES
110 FICE1(I,J,bi,bj)=QNET(I,J,bi,bj)+Qsw(I,J,bi,bj)
111 QSWO(I,J,bi,bj)=Qsw(I,J,bi,bj)
112 #else /* SEAICE_EXTERNAL_FLUXES undefined */
113 ALB(I,J)=SEAICE_waterAlbedo
114 A1(I,J)=(ONE-ALB(I,J))*SWDOWN(I,J,bi,bj)
115 & +LWDOWN(I,J,bi,bj)*0.97 _d 0
116 & +D1*UG(I,J)*ATEMP(I,J,bi,bj)+D1W*UG(I,J)*AQH(I,J,bi,bj)
117 B(I,J)=QS1*6.11 _d +00*EXP(17.2694 _d +00
118 & *(TICE(I,J,bi,bj)-TMELT)
119 & /(TICE(I,J,bi,bj)-TMELT+237.3 _d +00))
120 A2(I,J)=-D1*UG(I,J)*TICE(I,J,bi,bj)-D1W*UG(I,J)*B(I,J)
121 & -D3*(TICE(I,J,bi,bj)**4)
122 FICE1(I,J,bi,bj)=-A1(I,J)-A2(I,J)
123 QSWO(I,J,bi,bj)=-(ONE-ALB(I,J))*SWDOWN(I,J,bi,bj)
124 #endif /* SEAICE_EXTERNAL_FLUXES */
125 QNETO(I,J,bi,bj)=FICE1(I,J,bi,bj)-QSWO(I,J,bi,bj)
126 ENDDO
127 ENDDO
128
129 ELSE
130
131 C COME HERE IF ICE COVER
132 C FIRST PUT MINIMUM ON ICE THICKNESS
133 DO J=1,sNy
134 DO I=1,sNx
135 HICE(I,J)=MYMAX_R8(HICE1(I,J),0.05 _d +00)
136 HICE(I,J)=MYMIN_R8(HICE(I,J),9.0 _d +00)
137 ENDDO
138 ENDDO
139 C NOW DECIDE ON ALBEDO
140 DO J=1,sNy
141 DO I=1,sNx
142 ALB(I,J)=SEAICE_dryIceAlb
143 IF(TICE(I,J,bi,bj).GT.TMELTP) ALB(I,J)=SEAICE_wetIceAlb
144 ASNOW=SEAICE_drySnowAlb
145 IF(TICE(I,J,bi,bj).GT.TMELTP) ASNOW=SEAICE_wetSnowAlb
146 IF(HSNOW(I,J,bi,bj).GT.HCUT) THEN
147 ALB(I,J)=ASNOW
148 ELSE
149 ALB(I,J)=ALB(I,J)+(HSNOW(I,J,bi,bj)/HCUT)*(ASNOW-ALB(I,J))
150 IF(ALB(I,J).GT.ASNOW) ALB(I,J)=ASNOW
151 END IF
152 ENDDO
153 ENDDO
154 C NOW DETERMINE FIXED FORCING TERM IN HEAT BUDGET
155 DO J=1,sNy
156 DO I=1,sNx
157 IF(HSNOW(I,J,bi,bj).GT.0.0) THEN
158 C NO SW PENETRATION WITH SNOW
159 A1(I,J)=(ONE-ALB(I,J))*SWDOWN(I,J,bi,bj)
160 & +LWDOWN(I,J,bi,bj)*0.97 _d 0
161 & +D1*UG(I,J)*ATEMP(I,J,bi,bj)+D1I*UG(I,J)*AQH(I,J,bi,bj)
162 ELSE
163 C SW PENETRATION UNDER ICE
164 A1(I,J)=(ONE-ALB(I,J))*SWDOWN(I,J,bi,bj)
165 & *(ONE-XIO*EXP(-1.5 _d 0*HICE(I,J)))
166 & +LWDOWN(I,J,bi,bj)*0.97 _d 0
167 & +D1*UG(I,J)*ATEMP(I,J,bi,bj)+D1I*UG(I,J)*AQH(I,J,bi,bj)
168 ENDIF
169 ENDDO
170 ENDDO
171 C NOW COMPUTE OTHER TERMS IN HEAT BUDGET
172 C COME HERE AT START OF ITERATION
173
174 crg check wether a2 is needed in the list of variables
175 cdm Ralf, the line below causes following error message
176 cdm INTERNAL ERROR: cannot find var clone to ada2
177 cdm c$taf loop = iteration TICE,A2
178 cdm iterative solver for ice growth rate
179 cdm inputs: TICE ice temperature
180 cdm UG forcing
181 cdm HSNOW snow thickness
182 cdm HICE ice thickness
183 cdm outputs: A2 is needed for FICE1, which is ice growth rate
184 cdm TICE
185 DO ITER=1,IMAX_TICE
186
187 DO J=1,sNy
188 DO I=1,sNx
189 B(I,J)=QS1*(C1*TICE(I,J,bi,bj)**4+C2*TICE(I,J,bi,bj)**3
190 & +C3*TICE(I,J,bi,bj)**2+C4*TICE(I,J,bi,bj)+C5)
191 A2(I,J)=-D1*UG(I,J)*TICE(I,J,bi,bj)-D1I*UG(I,J)*B(I,J)
192 & -D3*(TICE(I,J,bi,bj)**4)
193 B(I,J)=XKS/(HSNOW(I,J,bi,bj)/HICE(I,J)+XKS/XKI)/HICE(I,J)
194 A3(I,J)=4.0 _d +00*D3*(TICE(I,J,bi,bj)**3)+B(I,J)+D1*UG(I,J)
195 B(I,J)=B(I,J)*(TB-TICE(I,J,bi,bj))
196 cdm
197 cdm if(TICE(I,J,bi,bj).le.206.)
198 cdm & print '(A,3i4,f12.2)','### ITER,I,J,TICE',
199 cdm & ITER,I,J,TICE(I,J,bi,bj)
200 cdm
201 ENDDO
202 ENDDO
203 C NOW DECIDE IF IT IS TIME TO ESTIMATE GROWTH RATES
204 C NOW DETERMINE NEW ICE TEMPERATURE
205 DO J=1,sNy
206 DO I=1,sNx
207 TICE(I,J,bi,bj)=TICE(I,J,bi,bj)
208 & +(A1(I,J)+A2(I,J)+B(I,J))/A3(I,J)
209 TICE(I,J,bi,bj)=MYMAX_R8(273.16 _d 0+MIN_TICE,TICE(I,J,bi,bj))
210 ENDDO
211 ENDDO
212 C NOW SET ICE TEMP TO MIN OF TMELT/ITERATION RESULT
213 DO J=1,sNy
214 DO I=1,sNx
215 TICE(I,J,bi,bj)=MYMIN_R8(TICE(I,J,bi,bj),TMELT)
216 ENDDO
217 ENDDO
218
219 C END OF ITERATION
220 ENDDO
221
222 DO J=1,sNy
223 DO I=1,sNx
224 FICE1(I,J,bi,bj)=-A1(I,J)-A2(I,J)
225 IF(HSNOW(I,J,bi,bj).GT.0.0) THEN
226 C NO SW PENETRATION WITH SNOW
227 QSWI(I,J,bi,bj)=ZERO
228 ELSE
229 C SW PENETRATION UNDER ICE
230 QSWI(I,J,bi,bj)=-(ONE-ALB(I,J))*SWDOWN(I,J,bi,bj)
231 & *XIO*EXP(-1.5 _d 0*HICE(I,J))
232 ENDIF
233 ENDDO
234 ENDDO
235
236 END IF
237
238 #endif /* ALLOW_SEAICE */
239
240 RETURN
241 END

  ViewVC Help
Powered by ViewVC 1.1.22