Changeset 4828 for LMDZ6/trunk/libf/phylmd
- Timestamp:
- Feb 17, 2024, 9:25:57 PM (12 months ago)
- Location:
- LMDZ6/trunk/libf/phylmd
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/trunk/libf/phylmd/lmdz_blowing_snow_sublim_sedim.F90
r4821 r4828 110 110 DO i = 1, ngrid 111 111 112 rhoair = pplay(i,k) / zt(i) / RD 113 dz = (paprs(i,k)-paprs(i,k+1)) / (rhoair*RG) 114 ! BS fall velocity assumed to be constant for now 115 zvelo(i) = fallv_bs 116 117 112 118 IF (zt(i) .GT. RTT) THEN 113 119 … … 117 123 118 124 taumeltbs=taumeltbs0*exp(-max(0.,(zt(i)-RTT)/(tbsmelt-RTT))) 119 dqmelt=min(zqbs(i), zqbs(i)/taumeltbs*dtime)125 dqmelt=min(zqbs(i),-1.*zqbs(i)*(exp(-dtime/taumeltbs)-1.)) 120 126 maxdqmelt= max(RCPD*(1.0+RVTMP2*(zq(i)+zqbs(i)))*(zt(i)-RTT)/(RLMLT+RLVTT),0.) 121 dqmelt=min( max(dqmelt,0.),maxdqmelt)127 dqmelt=min(dqmelt,maxdqmelt) 122 128 ! q update, melting + evaporation 123 129 zq(i) = zq(i) + dqmelt … … 130 136 ! negative celcius temperature 131 137 ! Sublimation scheme 132 138 133 139 rhoair=pplay(i,k)/zt(i)/RD 134 140 dqsub = 1./rhoair*coef_eva_bs*(1.0-zq(i)/qsi(i))*((rhoair*zvelo(i)*zqbs(i))**expo_eva_bs)*dtime … … 148 154 149 155 ! Sedimentation scheme 150 151 rhoair = pplay(i,k) / zt(i) / RD 152 dz = (paprs(i,k)-paprs(i,k+1)) / (rhoair*RG) 153 ! BS fall velocity assumed to be constant for now 154 zvelo(i) = fallv_bs 156 !---------------------- 155 157 156 158 sedimn(i) = rhoair*zqbs(i)*zvelo(i) … … 173 175 zt(i) = zt(i) - zqbs(i) * RLSTT/RCPD/(1.0+RVTMP2*(zq(i)+zqbs(i))) 174 176 ELSE 175 zt(i) = zt(i) - zqbs(i) * (RLVTT+RLMLT)/RCPD/(1.0+RVTMP2*(zq(i)+zqbs(i)))177 zt(i) = zt(i) - zqbs(i) * (RLVTT+RLMLT)/RCPD/(1.0+RVTMP2*(zq(i)+zqbs(i))) 176 178 ENDIF 177 179 zqbs(i)=0. -
LMDZ6/trunk/libf/phylmd/surf_landice_mod.F90
r4724 r4828 144 144 REAL :: tau_dens, tau_dens0, tau_densmin, rhomax, rhohard 145 145 REAL, DIMENSION(klon) :: ws1, rhos, ustart, qsalt, hsalt, snowini 146 REAL :: maxerosion ! in kg/s146 REAL :: ustarmin, maxerosion, tau_eq_salt 147 147 148 148 ! End definition … … 364 364 ! Simple blowing snow param 365 365 if (ok_bs) then 366 ustarmin=1.e-3 366 367 ustart0 = 0.211 367 368 rhoice = 920.0 … … 372 373 tau_densmin=21600.0 ! 1/4 days according to in situ obs by C. Amory during blowing snow + 373 374 ! Marshall et al. 1999 (snow densification during rain) 375 tau_eq_salt=10. 374 376 375 377 ! computation of threshold friction velocity … … 408 410 csalt=fluxsalt/(2.8*ustart(j)) 409 411 hsalt(j)=0.08436*ustar(j)**1.27 410 qsalt(j)=1./rhod*csalt*lambdasalt*RG/(max(ustar(j)**2, 1E-6)) &411 * exp(-lambdasalt*RG*hsalt(j)/max(ustar(j)**2, 1E-6))412 qsalt(j)=1./rhod*csalt*lambdasalt*RG/(max(ustar(j)**2,ustarmin)) & 413 * exp(-lambdasalt*RG*hsalt(j)/max(ustar(j)**2,ustarmin)) 412 414 qsalt(j)=max(qsalt(j),0.) 413 415 enddo … … 417 419 ! default formulation from MAR model (Amory et al. 2021, Gallee et al. 2001) 418 420 do j=1, knon 419 esalt=1./(3.25*max(ustar (j),0.001))420 hsalt(j)=0.08436* ustar(j)**1.27421 esalt=1./(3.25*max(ustarmin,ustar(j))) 422 hsalt(j)=0.08436*(max(ustarmin,ustar(j))**1.27) 421 423 qsalt(j)=(max(ustar(j)**2-ustart(j)**2,0.))/(RG*hsalt(j))*esalt 422 424 !ep=qsalt*cdragm(j)*sqrt(u1(j)**2+v1(j)**2) … … 429 431 i = knindex(j) 430 432 rhod=p1lay(j)/RD/temp_air(j) 431 ! Nemoto and Nishimura 2004 show that steady-state saltation is achieved within a n intervalof about 10s433 ! Nemoto and Nishimura 2004 show that steady-state saltation is achieved within a time tau_eq_salt of about 10s 432 434 ! we thus prevent snowerosion (snow particle transfer from the saltation layer to the first model level) 433 ! to exceed (in intensity) integrated over 10sto exceed the total mass of snow particle in the saltation layer435 ! integrated over tau_eq_salt to exceed the total mass of snow particle in the saltation layer 434 436 ! (rho*qsalt*hsalt) 435 maxerosion=hsalt(j)*qsalt(j)*rhod/ 10.437 maxerosion=hsalt(j)*qsalt(j)*rhod/tau_eq_salt 436 438 fluxbs(j)=rhod*ws1(j)*cdragm(j)*(AcoefQBS(j)-qsalt(j)) & 437 439 / (1.-rhod*ws1(j)*cdragm(j)*BcoefQBS(j)*dtime) 438 440 !fluxbs(j)= rhod*ws1(j)*cdragm(j)*(qbs1(j)-qsalt(j)) 439 fluxbs(j)=max(-maxerosion,fluxbs(j)) 440 441 fluxbs(j)=min(0.,max(-maxerosion,fluxbs(j))) 441 442 ! for outputs 442 443
Note: See TracChangeset
for help on using the changeset viewer.