forked from altMITgcm/MITgcm
-
Notifications
You must be signed in to change notification settings - Fork 0
/
aim_aim2sioce.F
241 lines (207 loc) · 8.19 KB
/
aim_aim2sioce.F
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
C $Header: /u/gcmpack/MITgcm/pkg/aim_v23/aim_aim2sioce.F,v 1.9 2013/05/02 20:10:14 jmc Exp $
C $Name: $
#include "AIM_OPTIONS.h"
#ifdef ALLOW_THSICE
#include "THSICE_OPTIONS.h"
#endif
CBOP
C !ROUTINE: AIM_AIM2SIOCE
C !INTERFACE:
SUBROUTINE AIM_AIM2SIOCE(
I land_frc, siceFrac,
O prcAtm, snowPrc,
I bi, bj, myTime, myIter, myThid)
C !DESCRIPTION: \bv
C *==========================================================*
C | S/R AIM_AIM2SIOCE
C | o Interface between AIM and thSIce pkg or (coupled) ocean
C *==========================================================*
C | o compute surface fluxes over ocean (ice-free + ice covered)
C | for diagnostics, thsice package and (slab, coupled) ocean
C *==========================================================*
C \ev
C !USES:
IMPLICIT NONE
C == Global variables ===
C-- size for MITgcm & Physics package :
#include "AIM_SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "FFIELDS.h"
C-- Physics package
#include "AIM_PARAMS.h"
#include "com_physcon.h"
#include "com_physvar.h"
#ifdef ALLOW_THSICE
#include "THSICE_SIZE.h"
#include "THSICE_PARAMS.h"
#include "THSICE_VARS.h"
#endif
C updated fields (in commom blocks):
C if using thSIce:
C Qsw(inp) :: SW radiation through the sea-ice down to the ocean (+=up)
C Qsw(out) :: SW radiation down to the ocean (ice-free + ice-covered)(+=up)
C Qnet(out) :: Net heat flux out of the ocean (ice-free ocean only)(+=up)
C and the Ice-Covered contribution will be added in S/R THSICE_STEP_FWD
C EmPmR(out) :: Net fresh water flux out off the ocean (ice-free ocean only)
C and the Ice-Covered contribution will be added in S/R THSICE_STEP_FWD
C sHeating(in/out) :: air - seaice surface heat flux left to melt the ice
C icFrwAtm :: Evaporation over sea-ice [kg/m2/s] (>0 if evaporate)
C icFlxSW :: net SW heat flux through the ice to the ocean [W/m2] (+=dw)
C if not using thSIce:
C Qsw(out) :: SW radiation down to the ocean (ice-free + ice-covered)(+=up)
C Qnet(out) :: Net heat flux out of the ocean (ice-free + ice-covered)(+=up)
C EmPmR(out) :: Net fresh water flux out off the ocean (ice-free + ice-covered)
C !INPUT/OUTPUT PARAMETERS:
C == Routine arguments ==
C land_frc :: land fraction [0-1]
C siceFrac :: sea-ice fraction (relative to full grid-cell) [0-1]
C prcAtm :: total precip from the atmosphere [kg/m2/s]
C snowPrc :: snow precip over sea-ice [kg/m2/s]
C bi,bj :: Tile indices
C myTime :: Current time of simulation ( s )
C myIter :: Current iteration number in simulation
C myThid :: My Thread Id number
_RS land_frc(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RL siceFrac(sNx,sNy)
_RL prcAtm (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
_RL snowPrc(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
INTEGER bi, bj, myIter, myThid
_RL myTime
CEOP
#ifdef ALLOW_AIM
C == Local variables ==
C i,j,I2 :: loop counters
C convPrcEvp :: units conversion factor for Precip & Evap:
C :: from AIM units (g/m2/s) to model EmPmR units ( kg/m2/s )
_RL convPrcEvp
_RL icFrac, opFrac
INTEGER i,j,I2
C-- Initialisation :
C-- Atmospheric Physics Fluxes
C from g/m2/s to kg/m2/s :
convPrcEvp = 1. _d -3
DO j=1,sNy
DO i=1,sNx
IF ( land_frc(i,j,bi,bj).GE.1. _d 0 ) THEN
C- Full Land grid-cell: set all fluxes to zero (this has no effect on the
C model integration and just put this to get meaningfull diagnostics)
prcAtm(i,j) = 0. _d 0
Qnet(i,j,bi,bj) = 0. _d 0
EmPmR(i,j,bi,bj)= 0. _d 0
Qsw(i,j,bi,bj) = 0. _d 0
ELSE
I2 = i+(j-1)*sNx
C- Total Precip (no distinction between ice-covered / ice-free fraction):
prcAtm(i,j) = ( PRECNV(I2,myThid)
& + PRECLS(I2,myThid) )
C- Net surface heat flux over ice-free ocean (+=down)
C note: with aim_splitSIOsFx=F, ice-free & ice covered contribution are
C already merged together and Qnet is the mean heat flux over the grid box.
Qnet(i,j,bi,bj) =
& SSR(I2,2,myThid)
& - SLR(I2,2,myThid)
& - SHF(I2,2,myThid)
& - EVAP(I2,2,myThid)*ALHC
C- E-P over ice-free ocean [m/s]: (same as above is aim_splitSIOsFx=F)
EmPmR(i,j,bi,bj) = ( EVAP(I2,2,myThid)
& - prcAtm(i,j) ) * convPrcEvp
C- Net short wave (ice-free ocean) into the ocean (+=down)
Qsw(i,j,bi,bj) = SSR(I2,2,myThid)
ENDIF
ENDDO
ENDDO
#ifdef ALLOW_THSICE
IF ( useThSIce ) THEN
DO j=1,sNy
DO i=1,sNx
I2 = i+(j-1)*sNx
C- Mixed-Layer Ocean: (for thsice slab_ocean and coupler)
C NOTE: masking is now applied much earlier, during initialisation
c IF (land_frc(i,j,bi,bj).EQ.1. _d 0) hOceMxL(i,j,bi,bj) = 0.
C- Evaporation over sea-ice: (for thsice)
icFrwAtm(i,j,bi,bj) = EVAP(I2,3,myThid)*convPrcEvp
C- short-wave downward heat flux (ice-free ocean + ice-covered):
C note: at this point we already called THSICE_IMPL_TEMP to solve for
C seaice temp and SW flux through the ice. SW is not modified after, and
C can therefore combine the open-ocean & ice-covered ocean SW fluxes.
icFrac = iceMask(i,j,bi,bj)
opFrac = 1. _d 0 - icFrac
Qsw(i,j,bi,bj) = icFrac*icFlxSW(i,j,bi,bj)
& + opFrac*Qsw(i,j,bi,bj)
ENDDO
ENDDO
IF ( aim_energPrecip ) THEN
C-- Add energy flux related to Precip. (snow, T_rain) over sea-ice
DO j=1,sNy
DO i=1,sNx
IF ( iceMask(i,j,bi,bj).GT.0. _d 0 ) THEN
I2 = i+(j-1)*sNx
IF ( EnPrec(I2,myThid).GE.0. _d 0 ) THEN
C- positive => add to surface heating
sHeating(i,j,bi,bj) = sHeating(i,j,bi,bj)
& + EnPrec(I2,myThid)*prcAtm(i,j)
snowPrc(i,j) = 0. _d 0
ELSE
C- negative => make snow
snowPrc(i,j) = prcAtm(i,j)*convPrcEvp
ENDIF
ELSE
snowPrc(i,j) = 0. _d 0
ENDIF
ENDDO
ENDDO
ENDIF
ELSEIF ( aim_splitSIOsFx ) THEN
#else /* ALLOW_THSICE */
IF ( aim_splitSIOsFx ) THEN
#endif /* ALLOW_THSICE */
C- aim_splitSIOsFx=T: fluxes over sea-ice (3) & ice-free ocean (2) were
C computed separately and here we merge the 2 fractions
DO j=1,sNy
DO i=1,sNx
I2 = i+(j-1)*sNx
IF ( siceFrac(i,j) .GT. 0. ) THEN
icFrac = siceFrac(i,j)/(1. _d 0 - land_frc(i,j,bi,bj))
opFrac = 1. _d 0 - icFrac
C- Net surface heat flux over sea-ice + ice-free ocean (+=down)
Qnet(i,j,bi,bj) = Qnet(i,j,bi,bj)*opFrac
& + ( SSR(I2,3,myThid)
& - SLR(I2,3,myThid)
& - SHF(I2,3,myThid)
& - EVAP(I2,3,myThid)*ALHC
& )*icFrac
C- E-P over sea-ice + ice-free ocean [m/s]:
EmPmR(i,j,bi,bj) = EmPmR(i,j,bi,bj)*opFrac
& + ( EVAP(I2,3,myThid)
& - prcAtm(i,j) ) * convPrcEvp * icFrac
C- Net short wave (ice-free ocean) into the ocean (+=down)
Qsw(i,j,bi,bj) = opFrac*Qsw(i,j,bi,bj)
& + icFrac*SSR(I2,3,myThid)
ENDIF
ENDDO
ENDDO
C-- end of If useThSIce / elseif aim_splitSIOsFx blocks
ENDIF
IF ( aim_energPrecip ) THEN
C-- Ice free fraction: Add energy flux related to Precip. (snow, T_rain):
DO j=1,sNy
DO i=1,sNx
I2 = i+(j-1)*sNx
Qnet(i,j,bi,bj) = Qnet(i,j,bi,bj)
& + EnPrec(I2,myThid)*prcAtm(i,j)
ENDDO
ENDDO
ENDIF
DO j=1,sNy
DO i=1,sNx
C- Total Precip : convert units
prcAtm(i,j) = prcAtm(i,j) * convPrcEvp
C- Oceanic convention: Heat flux are > 0 upward ; reverse sign.
Qsw(i,j,bi,bj) = -Qsw(i,j,bi,bj)
Qnet(i,j,bi,bj)= -Qnet(i,j,bi,bj)
ENDDO
ENDDO
#endif /* ALLOW_AIM */
RETURN
END