|
||||
File indexing completed on 2018-03-02 18:38:12 UTC
view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC88ac963ab8 Jean*0001 #include "BULK_FORCE_OPTIONS.h" 0002 0003 CBOP 0004 C !ROUTINE: BULKF_FORMULA_AIM 0005 C !INTERFACE: 0006 SUBROUTINE BULKF_FORMULA_AIM( 0007 I Tsurf, SLRD, 0008 I T1, T0, Q0, Vsurf, 0009 O SHF, EVAP, SLRU, 0010 O dEvp, sFlx, 0011 I iceornot, myThid ) 0012 0013 C !DESCRIPTION: \bv 0014 C *==========================================================* 0015 C | S/R BULKF_FORMULA_AIM 0016 C | o compute surface flux over ocean and sea-ice, 0017 C | using AIM surface flux formulation 0018 C *==========================================================* 0019 C *==========================================================* 0020 C \ev 0021 0022 C !USES: 0023 IMPLICIT NONE 0024 0025 C Resolution parameters 0026 0027 #include "EEPARAMS.h" 0028 #include "SIZE.h" 0029 #include "PARAMS.h" 0030 #include "BULKF_PARAMS.h" 0031 0032 C !INPUT/OUTPUT PARAMETERS: 0033 C == Routine Arguments == 0034 C-- Input: 0035 C FMASK :: fractional land-sea mask (2-dim) 0036 C Tsurf :: surface temperature (2-dim) 0037 C SSR :: sfc sw radiation (net flux) (2-dim) 0038 C SLRD :: sfc lw radiation (downward flux)(2-dim) 0039 C T1 :: near-surface air temperature (from Pot.temp) 0040 C T0 :: near-surface air temperature (2-dim) 0041 C Q0 :: near-surface sp. humidity [g/kg](2-dim) 0042 C Vsurf :: surface wind speed [m/s] (2-dim,input) 0043 C-- Output: 0044 C SHF :: sensible heat flux (2-dim) 0045 C EVAP :: evaporation [g/(m^2 s)] (2-dim) 0046 C SLRU :: sfc lw radiation (upward flux) (2-dim) 0047 C Shf0 :: sensible heat flux over freezing surf. 0048 C dShf :: sensible heat flux derivative relative to surf. temp 0049 C Evp0 :: evaporation computed over freezing surface (Ts=0.oC) 0050 C dEvp :: evaporation derivative relative to surf. temp 0051 C Slr0 :: upward long wave radiation over freezing surf. 0052 C dSlr :: upward long wave rad. derivative relative to surf. temp 0053 C sFlx :: net heat flux (+=down) except SW, function of surf. temp Ts: 0054 C 0: Flux(Ts=0.oC) ; 1: Flux(Ts^n) ; 2: d.Flux/d.Ts(Ts^n) 0055 C TSFC :: surface temperature (clim.) (2-dim) 0056 C TSKIN :: skin surface temperature (2-dim) 0057 C-- Input: 0058 C iceornot :: 0=open water, 1=ice cover 0059 C myThid :: Thread number for this instance of the routine 0060 C-- 0061 INTEGER NGP 0062 PARAMETER ( NGP = 1 ) 0063 c _RL PSA(NGP), FMASK(NGP), EMISloc 0064 _RL Tsurf(NGP) 0065 c _RL SSR(NGP) 0066 _RL SLRD(NGP) 0067 _RL T1(NGP), T0(NGP), Q0(NGP), Vsurf(NGP) 0068 0069 _RL SHF(NGP), EVAP(NGP), SLRU(NGP) 0070 _RL dEvp(NGP), sFlx(NGP,0:2) 0071 c _RL Shf0(NGP), dShf(NGP), Evp0(NGP), dEvp(NGP) 0072 c _RL Slr0(NGP), dSlr(NGP), sFlx(NGP,0:2) 0073 c _RL TSFC(NGP), TSKIN(NGP) 0074 0075 INTEGER iceornot 0076 INTEGER myThid 0077 CEOP 0078 0079 #ifdef ALLOW_FORMULA_AIM 0080 0081 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| 0082 C FWIND0 = ratio of near-sfc wind to lowest-level wind 0083 C CHS = heat exchange coefficient over sea 0084 C VGUST = wind speed for sub-grid-scale gusts 0085 C DTHETA = Potential temp. gradient for stability correction 0086 C dTstab = potential temp. increment for stability function derivative 0087 C FSTAB = Amplitude of stability correction (fraction) 0088 C P0 = reference pressure [Pa=N/m2] 0089 C GG = gravity accel. [m/s2] 0090 C RD = gas constant for dry air [J/kg/K] 0091 C CP = specific heat at constant pressure [J/kg/K] 0092 C ALHC = latent heat of condensation [J/g] 0093 C ALHF = latent heat of freezing [J/g] 0094 C SBC = Stefan-Boltzmann constant 0095 C EMISloc :: longwave surface emissivity 0096 c _RL FWIND0, CHS, VGUST, DTHETA, dTstab, FSTAB c8e7a08f49 Jean*0097 _RL P0, ALHC, ALHF, RD, CP, SBC, EMISloc 88ac963ab8 Jean*0098 EQUIVALENCE ( ocean_emissivity , EMISloc ) 0099 EQUIVALENCE ( Lvap , ALHC ) 0100 EQUIVALENCE ( Lfresh , ALHF ) 0101 EQUIVALENCE ( Rgas , RD ) 0102 EQUIVALENCE ( cpair , CP ) 0103 EQUIVALENCE ( stefan , SBC ) 0104 0105 C-- Local variables: 0106 C PSA :: norm. surface pressure [p/p0] (2-dim) 0107 C DENVV :: surface flux (sens,lat.) coeff. (=Rho*|V|) [kg/m2/s] 0108 C CDENVV :: surf. heat flux (sens.,lat.) coeff including stability effect 0109 C ALHevp :: Latent Heat of evaporation 0110 _RL PSA(NGP) c8e7a08f49 Jean*0111 _RL DENVV(NGP), PRD 88ac963ab8 Jean*0112 _RL Shf0(NGP), dShf(NGP), Evp0(NGP) 0113 _RL Slr0(NGP), dSlr(NGP) 0114 _RL TSFC(NGP), TSKIN(NGP) 0115 _RL CDENVV(NGP), RDTH, FSSICE 0116 _RL ALHevp, Fstb0, dTstb, dFstb 0117 _RL QSAT0(NGP,2) 0118 _RL QDUMMY(1), RDUMMY(1), TS2 0119 INTEGER J 0120 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| 0121 0122 PSA(1) = 1. _d 0 c8e7a08f49 Jean*0123 P0 = 1. _d +5 88ac963ab8 Jean*0124 C--- 0125 0126 ALHevp = ALHC 0127 C Evap of snow/ice: account for Latent Heat of freezing : c8e7a08f49 Jean*0128 c IF ( aim_energPrecip .OR. useThSIce ) ALHevp = ALHC + ALHF 0129 IF ( iceornot.GE.1 ) ALHevp = ALHC + ALHF 88ac963ab8 Jean*0130 0131 C 1.4 Density * wind speed (including gustiness factor) 0132 0133 PRD = P0/RD c8e7a08f49 Jean*0134 c VG2 = VGUST*VGUST 0135 c factWind2 = FWIND0*FWIND0 88ac963ab8 Jean*0136 0137 DO J=1,NGP 0138 c SPEED0(J)=SQRT(factWind2*Vsurf2(J)+VG2) 0139 c DENVV(J)=(PRD*PSA(J)/T0(J))*SPEED0(J) c8e7a08f49 Jean*0140 C-- assuming input file "WspeedFile" contains the time-average "SPEED0" 0141 C from AIM output (aimPhytave: fields # 15 ; aimDiag: WINDS ) : 88ac963ab8 Jean*0142 DENVV(J)=(PRD*PSA(J)/T0(J))*Vsurf(J) 0143 ENDDO 0144 0145 C 1.5 Define effective skin temperature to compensate for 0146 C non-linearity of heat/moisture fluxes during the daily cycle 0147 0148 DO J=1,NGP 0149 TSKIN(J) = Tsurf(J) + celsius2K 0150 TSFC(J)=273.16 _d 0 0151 ENDDO 0152 0153 C-- 2. Computation of fluxes over land and sea 0154 0155 C 2.1 Stability correction 0156 0157 RDTH = FSTAB/DTHETA 0158 0159 DO J=1,NGP 0160 FSSICE=1.+MIN(DTHETA,MAX(-DTHETA,TSKIN(J)-T1(J)))*RDTH 0161 CDENVV(J)=CHS*DENVV(J)*FSSICE 0162 ENDDO 0163 0164 IF ( dTstab.GT.0. _d 0 ) THEN 0165 C- account for stability function derivative relative to Tsurf: 0166 C note: to avoid discontinuity in the derivative (because of min,max), compute 0167 C the derivative using the discrete form: F(Ts+dTstab)-F(Ts-dTstab)/2.dTstab 0168 DO J=1,NGP 0169 Fstb0 = 1.+MIN(DTHETA,MAX(-DTHETA,TSFC(J) -T1(J)))*RDTH 0170 Shf0(J) = CHS*DENVV(J)*Fstb0 0171 dTstb = ( DTHETA+dTstab-ABS(TSKIN(J)-T1(J)) )/dTstab 0172 dFstb = RDTH*MIN(1. _d 0, MAX(0. _d 0, dTstb*0.5 _d 0)) 0173 dShf(J) = CHS*DENVV(J)*dFstb 0174 ENDDO 0175 ENDIF 0176 0177 C 2.2 Evaporation 0178 0179 CALL BULKF_SH2RH_AIM( 2, NGP, TSKIN, PSA, 1. _d 0, 0180 & QDUMMY, dEvp, QSAT0(1,1), myThid ) 0181 CALL BULKF_SH2RH_AIM( 0, NGP, TSFC, PSA, 1. _d 0, 0182 & QDUMMY, RDUMMY,QSAT0(1,2), myThid ) 0183 0184 IF ( dTstab.GT.0. _d 0 ) THEN 0185 C- account for stability function derivative relative to Tsurf: 0186 DO J=1,NGP 0187 EVAP(J) = CDENVV(J)*(QSAT0(J,1)-Q0(J)) 0188 Evp0(J) = Shf0(J)*(QSAT0(J,2)-Q0(J)) 0189 dEvp(J) = CDENVV(J)*dEvp(J) 0190 & + dShf(J)*(QSAT0(J,1)-Q0(J)) 0191 ENDDO 0192 ELSE 0193 DO J=1,NGP 0194 EVAP(J) = CDENVV(J)*(QSAT0(J,1)-Q0(J)) 0195 Evp0(J) = CDENVV(J)*(QSAT0(J,2)-Q0(J)) 0196 dEvp(J) = CDENVV(J)*dEvp(J) 0197 ENDDO 0198 ENDIF 0199 0200 C 2.3 Sensible heat flux 0201 0202 IF ( dTstab.GT.0. _d 0 ) THEN 0203 C- account for stability function derivative relative to Tsurf: 0204 DO J=1,NGP 0205 SHF(J) = CDENVV(J)*CP*(TSKIN(J)-T0(J)) 0206 Shf0(J) = Shf0(J)*CP*(TSFC(J) -T0(J)) 0207 dShf(J) = CDENVV(J)*CP 0208 & + dShf(J)*CP*(TSKIN(J)-T0(J)) 0209 dShf(J) = MAX( dShf(J), 0. _d 0 ) 0210 C-- do not allow negative derivative vs Ts of Sensible+Latent H.flux: 0211 C a) quiet unrealistic ; 0212 C b) garantee positive deriv. of total H.flux (needed for implicit solver) 0213 dEvp(J) = MAX( dEvp(J), -dShf(J)/ALHevp ) 0214 ENDDO 0215 ELSE 0216 DO J=1,NGP 0217 SHF(J) = CDENVV(J)*CP*(TSKIN(J)-T0(J)) 0218 Shf0(J) = CDENVV(J)*CP*(TSFC(J) -T0(J)) 0219 dShf(J) = CDENVV(J)*CP 0220 ENDDO 0221 ENDIF 0222 0223 C 2.4 Emission of lw radiation from the surface 0224 0225 DO J=1,NGP 0226 TS2 = TSFC(J)*TSFC(J) 0227 Slr0(J) = SBC*TS2*TS2 0228 TS2 = TSKIN(J)*TSKIN(J) 0229 SLRU(J) = SBC*TS2*TS2 0230 dSlr(J) = 4. _d 0 *SBC*TS2*TSKIN(J) 0231 ENDDO 0232 0233 C-- Compute net surface heat flux and its derivative ./. surf. temp. 0234 DO J=1,NGP 0235 sFlx(J,0)= ( SLRD(J) - EMISloc*Slr0(J) ) 0236 & - ( Shf0(J) + ALHevp*Evp0(J) ) 0237 sFlx(J,1)= ( SLRD(J) - EMISloc*SLRU(J) ) 0238 & - ( SHF(J) + ALHevp*EVAP(J) ) 0239 sFlx(J,2)= -EMISloc*dSlr(J) 0240 & - ( dShf(J) + ALHevp*dEvp(J) ) 0241 ENDDO 0242 0243 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| 0244 #endif /* ALLOW_FORMULA_AIM */ 0245 0246 RETURN 0247 END
[ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
This page was automatically generated from https://github.com/MITgcm/MITgcm by the 2.2.1-MITgcm-0.1 LXR engine. The LXR team |