/[MITgcm]/MITgcm/pkg/aim_v23/phy_suflux_sice.F
ViewVC logotype

Annotation of /MITgcm/pkg/aim_v23/phy_suflux_sice.F

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


Revision 1.5 - (hide annotations) (download)
Thu Jul 22 22:58:38 2004 UTC (19 years, 10 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57t_post, checkpoint57o_post, checkpoint57v_post, checkpoint57m_post, checkpoint55c_post, checkpoint54e_post, checkpoint57s_post, checkpoint57k_post, checkpoint55d_pre, checkpoint57d_post, checkpoint57g_post, checkpoint57b_post, checkpoint57c_pre, checkpoint55j_post, checkpoint56b_post, checkpoint57i_post, checkpoint57y_post, checkpoint57e_post, checkpoint55h_post, checkpoint57g_pre, checkpoint55b_post, checkpoint54d_post, checkpoint56c_post, checkpoint57y_pre, checkpoint55, checkpoint57f_pre, checkpoint57a_post, checkpoint54f_post, checkpoint55g_post, checkpoint55f_post, checkpoint57r_post, checkpoint58, checkpoint57a_pre, checkpoint55i_post, checkpoint57, checkpoint56, eckpoint57e_pre, checkpoint57h_done, checkpoint57x_post, checkpoint57n_post, checkpoint57w_post, checkpoint57p_post, checkpint57u_post, checkpoint57f_post, checkpoint58a_post, checkpoint57q_post, checkpoint57z_post, checkpoint57c_post, checkpoint55e_post, checkpoint55a_post, checkpoint57j_post, checkpoint57h_pre, checkpoint57l_post, checkpoint57h_post, checkpoint56a_post, checkpoint55d_post
Changes since 1.4: +35 -17 lines
fix a bug in SH & Lat flux derivative vs Ts ;
ensure positive deriv. vs Ts of SH+Lat heat-flux ;

1 jmc 1.5 C $Header: /u/gcmpack/MITgcm/pkg/aim_v23/phy_suflux_sice.F,v 1.4 2004/06/24 23:43:11 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "AIM_OPTIONS.h"
5    
6     CBOP
7     C !ROUTINE: SUFLUX_SICE
8     C !INTERFACE:
9     SUBROUTINE SUFLUX_SICE(
10     I PSA, FMASK, EMISloc,
11     I Tsurf, dTskin, SSR, SLRD,
12 jmc 1.4 I T1, T0, Q0, DENVV,
13 jmc 1.1 O SHF, EVAP, SLRU,
14 jmc 1.4 O Shf0, dShf, Evp0, dEvp, Slr0, dSlr, sFlx,
15 jmc 1.1 O TSFC, TSKIN,
16     I bi,bj,myThid)
17    
18     C !DESCRIPTION: \bv
19     C *==========================================================*
20     C | S/R SUFLUX_SICE
21     C | o compute surface flux over sea-ice
22     C *==========================================================*
23     C | o contains part of original S/R SUFLUX (Speedy code)
24     C *==========================================================*
25     C \ev
26    
27     C !USES:
28     IMPLICIT NONE
29    
30     C Resolution parameters
31    
32     C-- size for MITgcm & Physics package :
33     #include "AIM_SIZE.h"
34     #include "EEPARAMS.h"
35    
36     C-- Physics package
37     #include "AIM_PARAMS.h"
38    
39     C Physical constants + functions of sigma and latitude
40     #include "com_physcon.h"
41    
42     C Surface flux constants
43     #include "com_sflcon.h"
44    
45     C !INPUT/OUTPUT PARAMETERS:
46     C == Routine Arguments ==
47     C-- Input:
48     C PSA :: norm. surface pressure [p/p0] (2-dim)
49     C FMASK :: fractional land-sea mask (2-dim)
50     C EMISloc:: longwave surface emissivity
51     C Tsurf :: surface temperature (2-dim)
52     C dTskin :: temp. correction for daily-cycle heating [K]
53     C SSR :: sfc sw radiation (net flux) (2-dim)
54     C SLRD :: sfc lw radiation (downward flux)(2-dim)
55 jmc 1.4 C T1 :: near-surface air temperature (from Pot.temp)
56 jmc 1.1 C T0 :: near-surface air temperature (2-dim)
57     C Q0 :: near-surface sp. humidity [g/kg](2-dim)
58 jmc 1.4 C DENVV :: surface flux (sens,lat.) coeff. (=Rho*|V|) [kg/m2/s]
59 jmc 1.1 C-- Output:
60     C SHF :: sensible heat flux (2-dim)
61     C EVAP :: evaporation [g/(m^2 s)] (2-dim)
62     C SLRU :: sfc lw radiation (upward flux) (2-dim)
63 jmc 1.4 C Shf0 :: sensible heat flux over freezing surf.
64     C dShf :: sensible heat flux derivative relative to surf. temp
65 jmc 1.1 C Evp0 :: evaporation computed over freezing surface (Ts=0.oC)
66     C dEvp :: evaporation derivative relative to surf. temp
67     C Slr0 :: upward long wave radiation over freezing surf.
68     C dSlr :: upward long wave rad. derivative relative to surf. temp
69 jmc 1.2 C sFlx :: net heat flux (+=down) except SW, function of surf. temp Ts:
70 jmc 1.1 C 0: Flux(Ts=0.oC) ; 1: Flux(Ts^n) ; 2: d.Flux/d.Ts(Ts^n)
71     C TSFC :: surface temperature (clim.) (2-dim)
72     C TSKIN :: skin surface temperature (2-dim)
73     C-- Input:
74     C bi,bj :: tile index
75     C myThid :: Thread number for this instance of the routine
76     C--
77     _RL PSA(NGP), FMASK(NGP), EMISloc
78     _RL Tsurf(NGP), dTskin(NGP)
79     _RL SSR(NGP), SLRD(NGP)
80 jmc 1.4 _RL T1(NGP), T0(NGP), Q0(NGP), DENVV(NGP)
81 jmc 1.1
82     _RL SHF(NGP), EVAP(NGP), SLRU(NGP)
83 jmc 1.4 _RL Shf0(NGP), dShf(NGP), Evp0(NGP), dEvp(NGP)
84     _RL Slr0(NGP), dSlr(NGP), sFlx(NGP,0:2)
85 jmc 1.1 _RL TSFC(NGP), TSKIN(NGP)
86    
87     INTEGER bi,bj,myThid
88     CEOP
89    
90     #ifdef ALLOW_AIM
91    
92     C-- Local variables:
93 jmc 1.4 C CDENVV :: surf. heat flux (sens.,lat.) coeff including stability effect
94 jmc 1.5 C ALHevp :: Latent Heat of evaporation
95 jmc 1.4 _RL CDENVV(NGP), RDTH, FSSICE
96 jmc 1.5 _RL ALHevp, Fstb0, dTstb, dFstb
97 jmc 1.1 _RL QSAT0(NGP,2)
98     _RL QDUMMY(1), RDUMMY(1), TS2
99     INTEGER J
100    
101     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
102    
103 jmc 1.5 ALHevp = ALHC
104     C Evap of snow/ice: account for Latent Heat of freezing :
105     IF ( aim_energPrecip ) ALHevp = ALHC + ALHF
106    
107 jmc 1.1 C 1.5 Define effective skin temperature to compensate for
108     C non-linearity of heat/moisture fluxes during the daily cycle
109    
110     DO J=1,NGP
111     c TSKIN(J) = Tsurf(J) + dTskin(J)
112     c TSFC(J)=273.16 _d 0 + dTskin(J)
113     TSKIN(J) = Tsurf(J)
114     TSFC(J)=273.16 _d 0
115     ENDDO
116    
117     C-- 2. Computation of fluxes over land and sea
118    
119 jmc 1.4 C 2.1 Stability correction
120 jmc 1.1
121 jmc 1.4 RDTH = FSTAB/DTHETA
122 jmc 1.1
123     DO J=1,NGP
124 jmc 1.4 FSSICE=1.+MIN(DTHETA,MAX(-DTHETA,TSKIN(J)-T1(J)))*RDTH
125     CDENVV(J)=CHS*DENVV(J)*FSSICE
126 jmc 1.1 ENDDO
127    
128 jmc 1.4 IF ( dTstab.GT.0. _d 0 ) THEN
129     C- account for stability function derivative relative to Tsurf:
130     C note: to avoid discontinuity in the derivative (because of min,max), compute
131     C the derivative using the discrete form: F(Ts+dTstab)-F(Ts-dTstab)/2.dTstab
132     DO J=1,NGP
133     Fstb0 = 1.+MIN(DTHETA,MAX(-DTHETA,TSFC(J) -T1(J)))*RDTH
134 jmc 1.5 Shf0(J) = CHS*DENVV(J)*Fstb0
135 jmc 1.4 dTstb = ( DTHETA+dTstab-ABS(TSKIN(J)-T1(J)) )/dTstab
136     dFstb = RDTH*MIN(1. _d 0, MAX(0. _d 0, dTstb*0.5 _d 0))
137 jmc 1.5 dShf(J) = CHS*DENVV(J)*dFstb
138 jmc 1.4 ENDDO
139 jmc 1.5 C- deBug part:
140     c J = 6 + (17-1)*sNx
141     c IF ( bi.EQ.3 .AND. J.LE.NGP )
142     c & WRITE(6,1020)'SUFLUX_SICE: Stab=',Shf0(J),CDENVV(J),dShf(J)
143 jmc 1.4 ENDIF
144    
145     C 2.2 Evaporation
146 jmc 1.1
147     CALL SHTORH (2, NGP, TSKIN, PSA, 1. _d 0, QDUMMY, dEvp,
148     & QSAT0(1,1), myThid)
149     CALL SHTORH (0, NGP, TSFC, PSA, 1. _d 0, QDUMMY, RDUMMY,
150     & QSAT0(1,2), myThid)
151    
152 jmc 1.4 IF ( dTstab.GT.0. _d 0 ) THEN
153     C- account for stability function derivative relative to Tsurf:
154     DO J=1,NGP
155     EVAP(J) = CDENVV(J)*(QSAT0(J,1)-Q0(J))
156     Evp0(J) = Shf0(J)*(QSAT0(J,2)-Q0(J))
157     dEvp(J) = CDENVV(J)*dEvp(J)
158     & + dShf(J)*(QSAT0(J,1)-Q0(J))
159     ENDDO
160     ELSE
161     DO J=1,NGP
162 jmc 1.1 EVAP(J) = CDENVV(J)*(QSAT0(J,1)-Q0(J))
163     Evp0(J) = CDENVV(J)*(QSAT0(J,2)-Q0(J))
164     dEvp(J) = CDENVV(J)*dEvp(J)
165 jmc 1.4 ENDDO
166     ENDIF
167    
168     C 2.3 Sensible heat flux
169    
170     IF ( dTstab.GT.0. _d 0 ) THEN
171     C- account for stability function derivative relative to Tsurf:
172     DO J=1,NGP
173     SHF(J) = CDENVV(J)*CP*(TSKIN(J)-T0(J))
174     Shf0(J) = Shf0(J)*CP*(TSFC(J) -T0(J))
175     dShf(J) = CDENVV(J)*CP
176     & + dShf(J)*CP*(TSKIN(J)-T0(J))
177 jmc 1.5 dShf(J) = MAX( dShf(J), 0. _d 0 )
178     C-- do not allow negative derivative vs Ts of Sensible+Latent H.flux:
179     C a) quiet unrealistic ;
180     C b) garantee positive deriv. of total H.flux (needed for implicit solver)
181     dEvp(J) = MAX( dEvp(J), -dShf(J)/ALHevp )
182 jmc 1.4 ENDDO
183     ELSE
184     DO J=1,NGP
185     SHF(J) = CDENVV(J)*CP*(TSKIN(J)-T0(J))
186     Shf0(J) = CDENVV(J)*CP*(TSFC(J) -T0(J))
187     dShf(J) = CDENVV(J)*CP
188     ENDDO
189     ENDIF
190 jmc 1.1
191     C 2.4 Emission of lw radiation from the surface
192    
193     DO J=1,NGP
194     TS2 = TSFC(J)*TSFC(J)
195     Slr0(J) = SBC*TS2*TS2
196     TS2 = TSKIN(J)*TSKIN(J)
197     SLRU(J) = SBC*TS2*TS2
198     dSlr(J) = 4. _d 0 *SBC*TS2*TSKIN(J)
199     ENDDO
200    
201     C-- Compute net surface heat flux and its derivative ./. surf. temp.
202     DO J=1,NGP
203 jmc 1.4 sFlx(J,0)= ( SLRD(J) - EMISloc*Slr0(J) )
204 jmc 1.5 & - ( Shf0(J) + ALHevp*Evp0(J) )
205 jmc 1.4 sFlx(J,1)= ( SLRD(J) - EMISloc*SLRU(J) )
206 jmc 1.5 & - ( SHF(J) + ALHevp*EVAP(J) )
207 jmc 1.4 sFlx(J,2)= -EMISloc*dSlr(J)
208 jmc 1.5 & - ( dShf(J) + ALHevp*dEvp(J) )
209 jmc 1.1 ENDDO
210 jmc 1.5
211     C- deBug part: -----------------
212     c1010 FORMAT(A,I3,2F10.3,F10.4)
213     c1020 FORMAT(A,1P4E11.3)
214     c J = 6 + (17-1)*sNx
215     c IF ( bi.EQ.3 .AND. J.LE.NGP ) THEN
216     c WRITE(6,1010) 'SUFLUX_SICE: 1,sFlx=', 1,
217     c & sFlx(J,0),sFlx(J,1),sFlx(J,2)
218     c WRITE(6,1010) 'SUFLUX_SICE: 0,Evap=', 0,Evp0(J),EVAP(J),dEvp(J)
219     c WRITE(6,1010) 'SUFLUX_SICE: -,LWup=',-1,Slr0(J),SLRU(J),dSlr(J)
220     c WRITE(6,1010) 'SUFLUX_SICE: -, SHF=',-1,Shf0(J),SHF(J), dShf(J)
221     c WRITE(6,1010) 'SUFLUX_SICE: -, LAT=',-1,
222     c & ALHevp*Evp0(J),ALHevp*EVAP(J),ALHevp*dEvp(J)
223     c ENDIF
224 jmc 1.1
225     C-- 3. Adjustment of skin temperature and fluxes over land
226     C-- based on energy balance (to be implemented)
227     C <= done separately for each surface type (land,ocean,sea-ice)
228    
229     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
230     #endif /* ALLOW_AIM */
231    
232     RETURN
233     END

  ViewVC Help
Powered by ViewVC 1.1.22