Changeset 650
- Timestamp:
- May 4, 2012, 6:35:58 PM (13 years ago)
- Location:
- trunk/LMDZ.GENERIC
- Files:
-
- 15 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.GENERIC/README
r649 r650 661 661 - rcp, cpp and mugaz can now be computed using gases.def in newstart 662 662 - Correction of a bug arising in gcm.F when the solar days are long (thanks Melanie V.) 663 663 - Corrected the temperature used to differentiate sublimation and evaporation in watersat_grad 664 - Minor name changes in watercommon 665 - Better physical parametrization of the effective radius of liquid and icy water cloud particles in callcorrk 666 (for radfixed=true) 667 - Added consistency check in inifis 668 - Moved 1d water initialization from physiqu to rcm1d -
trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90
r622 r650 172 172 173 173 ! are particle radii fixed? 174 logical radfixed174 logical, save :: radfixed 175 175 real maxrad, minrad 176 176 177 ! Local water cloud optical properties 178 real, parameter :: rad_froid=35.0e-6 179 real, parameter :: rad_chaud=10.0e-6 180 real, parameter :: coef_chaud=0.13 181 real, parameter :: coef_froid=0.09 182 real zfice 183 177 184 real CBRT 178 185 external CBRT … … 183 190 REAL*8 muvarrad(L_LEVELS) 184 191 185 radfixed=.false.186 187 192 !======================================================================= 188 193 ! Initialization on first call … … 195 200 qsiaer(:,:,:)=0.0 196 201 giaer(:,:,:) =0.0 202 radfixed=.false. 197 203 198 204 199 205 if(firstcall) then 206 if(kastprof)then 207 radfixed=.true. 208 endif 200 209 201 210 call system('rm -f surf_vals_long.out') … … 305 314 enddo 306 315 307 if(kastprof)then 308 radfixed=.true. 309 endif 316 310 317 311 318 if(radfixed)then … … 320 327 do ig=1,ngrid 321 328 !reffrad(ig,l,2) = 2.e-5 ! H2O ice 322 reffrad(ig,l,2) = 5.e-6 ! H2O ice 329 reffrad(ig,l,2) = rad_chaud ! H2O ice 330 331 zfice = 1.0 - (pt(ig,l)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds) 332 zfice = MIN(MAX(zfice,0.0),1.0) 333 reffrad(ig,l,2)= rad_chaud * (1.-zfice) + rad_froid * zfice 334 nueffrad(ig,l,2) = coef_chaud * (1.-zfice) + coef_froid * zfice 323 335 enddo 324 336 enddo -
trunk/LMDZ.GENERIC/libf/phystd/callkeys.h
r596 r650 11 11 & , Nmix_co2, Nmix_h2o & 12 12 & , dusttau & 13 & , nonideal &14 13 & , meanOLR & 15 14 & , specOLR & -
trunk/LMDZ.GENERIC/libf/phystd/evap.F
r253 r650 61 61 ! ignoring qevap term creates huge difference when qevap large!!! 62 62 63 zdelta = MAX(0.,SIGN(1.,T o-tevap(ig,l))) ! what is this?63 zdelta = MAX(0.,SIGN(1.,T_h2O_ice_liq-tevap(ig,l))) ! what is this? 64 64 ! for division between water / liquid 65 65 dqevap(ig,l) = MAX(0.0,qevap(ig,l,igcm_h2o_ice))/dtime -
trunk/LMDZ.GENERIC/libf/phystd/hydrol.F90
r368 r650 3 3 albedo0,albedo,mu0,pdtsurf,pdtsurf_hyd,hice) 4 4 5 use watercommon_h, only: T o, RLFTT, rhowater, mx_eau_sol5 use watercommon_h, only: T_h2O_ice_liq, RLFTT, rhowater, mx_eau_sol 6 6 7 7 implicit none … … 189 189 hice(ig) = 0. 190 190 191 if(twater(ig) .lt. T o)then192 E=min((T o+Tsaldiff-twater(ig))*pcapcal(ig),RLFTT*rhowater*maxicethick)191 if(twater(ig) .lt. T_h2O_ice_liq)then 192 E=min((T_h2O_ice_liq+Tsaldiff-twater(ig))*pcapcal(ig),RLFTT*rhowater*maxicethick) 193 193 hice(ig) = E/(RLFTT*rhowater) 194 194 hice(ig) = max(hice(ig),0.0) … … 221 221 222 222 ! melt the snow 223 if(ztsurf(ig).gt.T o)then223 if(ztsurf(ig).gt.T_h2O_ice_liq)then 224 224 if(zqsurf(ig,iice).gt.1.0e-8)then 225 225 226 a = (ztsurf(ig)-T o)*pcapcal(ig)/RLFTT226 a = (ztsurf(ig)-T_h2O_ice_liq)*pcapcal(ig)/RLFTT 227 227 b = zqsurf(ig,iice) 228 228 fsnoi = min(a,b) … … 240 240 if(zqsurf(ig,iliq).gt.1.0e-8)then 241 241 242 a = -(ztsurf(ig)-T o)*pcapcal(ig)/RLFTT242 a = -(ztsurf(ig)-T_h2O_ice_liq)*pcapcal(ig)/RLFTT 243 243 b = zqsurf(ig,iliq) 244 244 -
trunk/LMDZ.GENERIC/libf/phystd/inifis.F
r622 r650 219 219 call getin("calladj",calladj) 220 220 write(*,*) " calladj = ",calladj 221 222 223 ! Test of incompatibility224 if (enertest.and.nonideal) then225 print*,'Energy conservation calculations currently226 & assume ideal gas!'227 call abort228 endif229 221 230 222 write(*,*) "call CO2 condensation ?" … … 232 224 call getin("co2cond",co2cond) 233 225 write(*,*) " co2cond = ",co2cond 234 226 ! Test of incompatibility 227 if (co2cond.and.(.not.tracer)) then 228 print*,'We need a CO2 ice tracer to condense CO2' 229 call abort 230 endif 231 235 232 write(*,*) "CO2 supersaturation level ?" 236 233 co2supsat=1.0 ! default value -
trunk/LMDZ.GENERIC/libf/phystd/moistadj.F90
r253 r650 1 1 subroutine moistadj(t, pq, pplev, pplay, dtmana, dqmana, ptimestep, rneb) 2 2 3 use watercommon_h, only: T o, RLVTT, RCPD3 use watercommon_h, only: T_h2O_ice_liq, RLVTT, RCPD 4 4 5 5 implicit none … … 124 124 DO i = 1, ngridmx 125 125 v_cptt(i,k) = RCPD * local_t(i,k) 126 ENDDO127 ENDDO128 129 DO k = 1, nlayermx130 DO i = 1, ngridmx131 v_cptt(i,k) = RCPD * local_t(i,k)132 126 v_t = local_t(i,k) 133 127 v_p = pplay(i,k) … … 142 136 ! DO i = 1, ngridmx 143 137 ! v_t = local_t(i,k) 144 ! IF (v_t.LT.T o) THEN138 ! IF (v_t.LT.T_h2O_ice_liq) THEN 145 139 ! print*,'RHs=',q(i,k) / v_qs(i,k) 146 140 ! ELSE -
trunk/LMDZ.GENERIC/libf/phystd/physiq.F90
r622 r650 8 8 9 9 use radinc_h, only : naerkind,L_NSPECTI,L_NSPECTV 10 use watercommon_h, only : mx_eau_sol, To, RLVTT, mH2O10 use watercommon_h, only : RLVTT 11 11 use gases_h 12 12 use radcommon_h, only: sigma … … 370 370 real vdifcncons(ngridmx) 371 371 real cadjncons(ngridmx) 372 real qzero1D373 save qzero1D374 372 375 373 ! double precision qsurf_hist(ngridmx,nqmx) … … 544 542 ! initialise variables for water cycle 545 543 546 qzero1D=0.0547 548 544 if(enertest)then 549 545 watertest = .true. 550 endif551 552 if(ngrid.eq.1)then553 qzero1D = 0.0554 qsurf(1,igcm_h2o_vap) = qzero1D555 do l=1, nlayer556 pq(1,l,igcm_h2o_vap)=0.0557 pq(1,l,igcm_h2o_ice)=0.0558 enddo559 546 endif 560 547 -
trunk/LMDZ.GENERIC/libf/phystd/rain.F90
r622 r650 2 2 3 3 4 use watercommon_h, only: T o, RLVTT, RCPD, RCPV, RV, RVTMP24 use watercommon_h, only: T_h2O_ice_liq,T_h2O_ice_clouds, RLVTT, RCPD, RCPV, RV, RVTMP2 5 5 6 6 implicit none … … 82 82 83 83 REAL zoliq(ngridmx) 84 REAL ztglace85 84 REAL zdz(ngridmx),zrho(ngridmx),ztot(ngridmx), zrhol(ngridmx) 86 85 REAL zchau(ngridmx),zfroi(ngridmx),zfrac(ngridmx),zneb(ngridmx) … … 140 139 ENDDO 141 140 ENDDO 142 143 ! Determine the cold clouds by their temperature144 ztglace = To - 15.0145 141 146 142 ! Initialise the outputs … … 219 215 DO i = 1, ngridmx 220 216 ttemp = zt(i,k) 221 IF (ttemp .ge. T o) THEN217 IF (ttemp .ge. T_h2O_ice_liq) THEN 222 218 lconvert=rainthreshold 223 219 ELSEIF (ttemp .gt. t_crit) THEN … … 245 241 zrho(i) = pplay(i,k) / ( zt(i,k) * R ) 246 242 zdz(i) = (pplev(i,k)-pplev(i,k+1)) / (zrho(i)*g) 247 zfrac(i) = (zt(i,k)- ztglace) / (To-ztglace)243 zfrac(i) = (zt(i,k)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds) 248 244 zfrac(i) = MAX(zfrac(i), 0.0) 249 245 zfrac(i) = MIN(zfrac(i), 1.0) … … 291 287 call abort 292 288 endif 293 IF (t(i,1) .LT. T o) THEN289 IF (t(i,1) .LT. T_h2O_ice_liq) THEN 294 290 dqssnow(i) = zrfl(i) 295 291 dqsrain(i) = 0.0 -
trunk/LMDZ.GENERIC/libf/phystd/rcm1d.F
r601 r650 73 73 ! REAL co2ice ! co2ice layer (kg.m-2) !not used anymore 74 74 integer :: i_co2_ice=0 ! tracer index of co2 ice 75 integer :: i_h2o_ice=0 ! tracer index of co2 ice 75 integer :: i_h2o_ice=0 ! tracer index of h2o ice 76 integer :: i_h2o_vap=0 ! tracer index of h2o vapor 76 77 REAL emis ! surface layer 77 78 REAL q2(nlayermx+1) ! Turbulent Kinetic Energy … … 111 112 REAL cloudfrac(ngridmx,nlayermx) 112 113 REAL hice(ngridmx),totcloudfrac(ngridmx) 114 REAL qzero1D !initial water amount on the ground 113 115 114 116 c======================================================================= … … 228 230 i_co2_ice=0 229 231 i_h2o_ice=0 232 i_h2o_vap=0 230 233 do iq=1,nq 231 234 if (tnom(iq)=="co2_ice") then … … 233 236 elseif (tnom(iq)=="h2o_ice") then 234 237 i_h2o_ice=iq 238 elseif (tnom(iq)=="h2o_vap") then 239 i_h2o_vap=iq 235 240 endif 236 241 enddo … … 275 280 area(1)=1.E+0 276 281 aire(1)=area(1) !JL+EM to have access to the area in the diagfi.nc files. area in comgeomfi.h and aire in comgeom.h 277 !!! GEOPOTENTIAL. useless in 1D because control b usurface pressure282 !!! GEOPOTENTIAL. useless in 1D because control by surface pressure 278 283 phisfi(1)=0.E+0 279 284 !!! LATITUDE. can be set. … … 493 498 qsurf(iq) = 0. 494 499 ENDDO 500 501 if (water) then 502 qzero1D = 0.0 503 qsurf(i_h2o_vap) = qzero1D 504 endif 495 505 496 506 c Initialisation pour prendre en compte les vents en 1-D -
trunk/LMDZ.GENERIC/libf/phystd/turbdiff.F90
r600 r650 7 7 pdqdif,pdqsdif,lastcall) 8 8 9 use watercommon_h, only : RLVTT, T o, RCPD, mx_eau_sol9 use watercommon_h, only : RLVTT, T_h2O_ice_liq, RCPD, mx_eau_sol 10 10 use radcommon_h, only : sigma 11 11 … … 609 609 else !dqsdif_total(ig).ge.0 610 610 !If water vapour is condensing, we must decide whether it forms ice or liquid. 611 if(ztsrf(ig).gt.T o)then611 if(ztsrf(ig).gt.T_h2O_ice_liq)then 612 612 pdqsdif(ig,iice_surf)=0.0 613 613 pdqsdif(ig,iliq_surf)=dqsdif_total(ig)/ptimestep -
trunk/LMDZ.GENERIC/libf/phystd/vdifc.F
r600 r650 7 7 & pdqdif,pdqsdif,lastcall) 8 8 9 use watercommon_h, only : RLVTT, T o, RCPD, mx_eau_sol9 use watercommon_h, only : RLVTT, T_h2O_ice_liq, RCPD, mx_eau_sol 10 10 use radcommon_h, only : sigma 11 11 … … 644 644 ! If water vapour is condensing, we must decide whether it forms ice or liquid. 645 645 if(dqsdif_total(ig).gt.0)then ! a bug was here! 646 if(ztsrf2(ig).gt.T o)then646 if(ztsrf2(ig).gt.T_h2O_ice_liq)then 647 647 pdqsdif(ig,iice)=0.0 648 648 pdqsdif(ig,ivap)=dqsdif_total(ig)/ptimestep … … 704 704 Wdiff=Wtot/ptimestep+pdqsdif(ig,ivap)+pdqsdif(ig,iice) 705 705 706 if(ztsrf2(ig).gt.T o)then706 if(ztsrf2(ig).gt.T_h2O_ice_liq)then 707 707 pdqsdif(ig,ivap)=pdqsdif(ig,ivap)-Wdiff 708 708 else -
trunk/LMDZ.GENERIC/libf/phystd/watercommon_h.F90
r253 r650 4 4 5 5 real, parameter :: T_coup = 234.0 6 real, parameter :: To = 273.16 6 real, parameter :: T_h2O_ice_liq = 273.16 7 real, parameter :: T_h2O_ice_clouds = T_h2O_ice_liq-15. 7 8 real, parameter :: mH2O = 18.01528 8 9 9 10 ! benjamin additions 10 real, parameter :: RLVTT = 2.257E+6 ! Latent heat of sublimation (J kg-1)11 real, parameter :: RLVTT = 2.257E+6 ! Latent heat of vaporization (J kg-1) 11 12 real, parameter :: RLSTT = 2.257E+6 ! 2.591E+6 in reality ! Latent heat of sublimation (J kg-1) 12 13 !real, parameter :: RLVTT = 0.0 13 14 !real, parameter :: RLSTT = 0.0 14 15 15 real, parameter :: RLFTT = 3.334E+5 ! Latent heat of fusion (J kg-1) 16 real, parameter :: RLFTT = 3.334E+5 ! Latent heat of fusion (J kg-1) ! entails an energy sink but better description of albedo 16 17 !real, parameter :: RLFTT = 0.0 17 18 real, parameter :: rhowater = 1.0E+3 ! mass of water (kg/m^3) 19 real, parameter :: capcal_h2o_liq = 4181.3 ! specific heat capacity of liquid water J/kg/K 18 20 real, parameter :: mx_eau_sol = 150 ! mass of water (kg/m^2) 19 21 … … 21 23 ! lw=(2834.3-0.28*(zt(ig,l)-To)-0.004*(zt(ig,l)-To)**2)*1.e+3 22 24 23 real epsi, RCPD, RCPV, RV, RVTMP2 24 save epsi, RCPD, RCPV, RV, RVTMP2 25 real, save :: epsi, RCPD, RCPV, RV, RVTMP2 25 26 26 27 end module watercommon_h -
trunk/LMDZ.GENERIC/libf/phystd/watersat.F90
r253 r650 1 1 subroutine watersat(T,p,qsat) 2 2 3 use watercommon_h, only: T o, epsi3 use watercommon_h, only: T_h2O_ice_liq, epsi 4 4 implicit none 5 5 … … 28 28 ! x by epsi converts to mass mixing ratio 29 29 30 if (T.lt.T o) then ! solid / vapour30 if (T.lt.T_h2O_ice_liq) then ! solid / vapour 31 31 qsat = 100.0 * epsi * 10**(2.07023 - 0.00320991 & 32 32 * T - 2484.896 / T + 3.56654 * alog10(T)) -
trunk/LMDZ.GENERIC/libf/phystd/watersat_grad.F90
r135 r650 1 1 subroutine watersat_grad(T,qsat,dqsat) 2 2 3 use watercommon_h, only: T_ coup, RLVTT, RCPD3 use watercommon_h, only: T_h2O_ice_liq, RLVTT, RCPD,T_coup 4 4 implicit none 5 5 … … 7 7 ! Purpose 8 8 ! ------- 9 ! Compute the water mass mixing ratio at saturation (kg/kg)10 ! for a given pressure (Pa) andtemperature (K)9 ! Compute the L/cp*d (q_sat)/d T 10 ! for a given temperature (K) 11 11 ! 12 12 ! Authors … … 22 22 real dqsat 23 23 24 if (T.lt.T_coup) then ! solid / vapour 24 ! if (T.lt.T_coup) then ! solid / vapour !why use T_coup?????????? JL12 25 if (T.lt.T_h2O_ice_liq) then ! solid / vapour 25 26 dqsat = RLVTT/RCPD*qsat*(3.56654/T & 26 27 +2484.896*LOG(10.)/T**2 &
Note: See TracChangeset
for help on using the changeset viewer.