346subroutine micro_mg_tend ( &
347 mgncol, nlev, deltatin, &
353 relvar, accre_enhan_i, &
355 cldn, liqcldf, icecldf, qsatfac, &
356 qcsinksum_rate1ord, &
364 effc, effc_fn, effi, &
375 reff_rain, reff_snow, &
376 qcsevap, qisevap, qvres, &
377 cmeitot, vtrmc, vtrmi, &
379 qcsedten, qisedten, &
380 qrsedten, qssedten, &
382 mnuccctot, mnuccttot, msacwitot, &
383 psacwstot, bergstot, bergtot, &
385 qcrestot, prcitot, praitot, &
386 qirestot, mnuccrtot, pracstot, &
387 meltsdttot, frzrdttot, mnuccdtot, &
389 refl, arefl, areflz, &
390 frefl, csrfl, acsrfl, &
398 prer_evap, xlat, xlon, lprnt, iccn, nlball)
437 integer,
intent(in) :: mgncol
438 integer,
intent(in) :: nlev
439 integer,
intent(in) :: nlball(mgncol)
440 real(r8),
intent(in) :: xlat,xlon
441 real(r8),
intent(in) :: deltatin
442 real(r8),
intent(in) :: t(mgncol,nlev)
443 real(r8),
intent(in) :: q(mgncol,nlev)
446 real(r8),
intent(in) :: qcn(mgncol,nlev)
447 real(r8),
intent(in) :: qin(mgncol,nlev)
448 real(r8),
intent(in) :: ncn(mgncol,nlev)
449 real(r8),
intent(in) :: nin(mgncol,nlev)
451 real(r8),
intent(in) :: qrn(mgncol,nlev)
452 real(r8),
intent(in) :: qsn(mgncol,nlev)
453 real(r8),
intent(in) :: nrn(mgncol,nlev)
454 real(r8),
intent(in) :: nsn(mgncol,nlev)
456 real(r8) :: relvar(mgncol,nlev)
457 real(r8) :: accre_enhan(mgncol,nlev)
459 real(r8),
intent(in) :: accre_enhan_i
462 real(r8),
intent(in) :: p(mgncol,nlev)
463 real(r8),
intent(in) :: pdel(mgncol,nlev)
465 real(r8),
intent(in) :: cldn(mgncol,nlev)
466 real(r8),
intent(in) :: liqcldf(mgncol,nlev)
467 real(r8),
intent(in) :: icecldf(mgncol,nlev)
468 real(r8),
intent(in) :: qsatfac(mgncol,nlev)
469 logical,
intent(in) :: lprnt
470 integer,
intent(in) :: iccn
475 real(r8),
intent(inout) :: naai(mgncol,nlev)
476 real(r8),
intent(in) :: npccnin(mgncol,nlev)
478 real(r8) :: npccn(mgncol,nlev)
482 real(r8),
intent(in) :: rndst(mgncol,nlev,10)
483 real(r8),
intent(in) :: nacon(mgncol,nlev,10)
487 real(r8),
intent(out) :: qcsinksum_rate1ord(mgncol,nlev)
489 real(r8),
intent(out) :: tlat(mgncol,nlev)
490 real(r8),
intent(out) :: qvlat(mgncol,nlev)
491 real(r8),
intent(out) :: qctend(mgncol,nlev)
492 real(r8),
intent(out) :: qitend(mgncol,nlev)
493 real(r8),
intent(out) :: nctend(mgncol,nlev)
494 real(r8),
intent(out) :: nitend(mgncol,nlev)
496 real(r8),
intent(out) :: qrtend(mgncol,nlev)
497 real(r8),
intent(out) :: qstend(mgncol,nlev)
498 real(r8),
intent(out) :: nrtend(mgncol,nlev)
499 real(r8),
intent(out) :: nstend(mgncol,nlev)
500 real(r8),
intent(out) :: effc(mgncol,nlev)
501 real(r8),
intent(out) :: effc_fn(mgncol,nlev)
502 real(r8),
intent(out) :: effi(mgncol,nlev)
503 real(r8),
intent(out) :: sadice(mgncol,nlev)
504 real(r8),
intent(out) :: sadsnow(mgncol,nlev)
505 real(r8),
intent(out) :: prect(mgncol)
506 real(r8),
intent(out) :: preci(mgncol)
507 real(r8),
intent(out) :: nevapr(mgncol,nlev)
508 real(r8),
intent(out) :: evapsnow(mgncol,nlev)
509 real(r8),
intent(out) :: am_evp_st(mgncol,nlev)
510 real(r8),
intent(out) :: prain(mgncol,nlev)
511 real(r8),
intent(out) :: prodsnow(mgncol,nlev)
512 real(r8),
intent(out) :: cmeout(mgncol,nlev)
513 real(r8),
intent(out) :: deffi(mgncol,nlev)
514 real(r8),
intent(out) :: pgamrad(mgncol,nlev)
515 real(r8),
intent(out) :: lamcrad(mgncol,nlev)
516 real(r8),
intent(out) :: qsout(mgncol,nlev)
517 real(r8),
intent(out) :: dsout(mgncol,nlev)
518 real(r8),
intent(out) :: lflx(mgncol,2:nlev+1)
519 real(r8),
intent(out) :: iflx(mgncol,2:nlev+1)
520 real(r8),
intent(out) :: rflx(mgncol,2:nlev+1)
521 real(r8),
intent(out) ::
sflx(mgncol,2:nlev+1)
522 real(r8),
intent(out) :: qrout(mgncol,nlev)
523 real(r8),
intent(out) :: reff_rain(mgncol,nlev)
524 real(r8),
intent(out) :: reff_snow(mgncol,nlev)
525 real(r8),
intent(out) :: qcsevap(mgncol,nlev)
526 real(r8),
intent(out) :: qisevap(mgncol,nlev)
527 real(r8),
intent(out) :: qvres(mgncol,nlev)
528 real(r8),
intent(out) :: cmeitot(mgncol,nlev)
529 real(r8),
intent(out) :: vtrmc(mgncol,nlev)
530 real(r8),
intent(out) :: vtrmi(mgncol,nlev)
531 real(r8),
intent(out) :: umr(mgncol,nlev)
532 real(r8),
intent(out) :: ums(mgncol,nlev)
533 real(r8),
intent(out) :: qcsedten(mgncol,nlev)
534 real(r8),
intent(out) :: qisedten(mgncol,nlev)
535 real(r8),
intent(out) :: qrsedten(mgncol,nlev)
536 real(r8),
intent(out) :: qssedten(mgncol,nlev)
539 real(r8),
intent(out) :: pratot(mgncol,nlev)
540 real(r8),
intent(out) :: prctot(mgncol,nlev)
541 real(r8),
intent(out) :: mnuccctot(mgncol,nlev)
542 real(r8),
intent(out) :: mnuccttot(mgncol,nlev)
543 real(r8),
intent(out) :: msacwitot(mgncol,nlev)
544 real(r8),
intent(out) :: psacwstot(mgncol,nlev)
545 real(r8),
intent(out) :: bergstot(mgncol,nlev)
546 real(r8),
intent(out) :: bergtot(mgncol,nlev)
547 real(r8),
intent(out) :: melttot(mgncol,nlev)
548 real(r8),
intent(out) :: homotot(mgncol,nlev)
549 real(r8),
intent(out) :: qcrestot(mgncol,nlev)
550 real(r8),
intent(out) :: prcitot(mgncol,nlev)
551 real(r8),
intent(out) :: praitot(mgncol,nlev)
552 real(r8),
intent(out) :: qirestot(mgncol,nlev)
553 real(r8),
intent(out) :: mnuccrtot(mgncol,nlev)
554 real(r8),
intent(out) :: pracstot(mgncol,nlev)
555 real(r8),
intent(out) :: meltsdttot(mgncol,nlev)
556 real(r8),
intent(out) :: frzrdttot(mgncol,nlev)
557 real(r8),
intent(out) :: mnuccdtot(mgncol,nlev)
558 real(r8),
intent(out) :: nrout(mgncol,nlev)
559 real(r8),
intent(out) :: nsout(mgncol,nlev)
560 real(r8),
intent(out) :: refl(mgncol,nlev)
561 real(r8),
intent(out) :: arefl(mgncol,nlev)
562 real(r8),
intent(out) :: areflz(mgncol,nlev)
563 real(r8),
intent(out) :: frefl(mgncol,nlev)
564 real(r8),
intent(out) :: csrfl(mgncol,nlev)
565 real(r8),
intent(out) :: acsrfl(mgncol,nlev)
566 real(r8),
intent(out) :: fcsrfl(mgncol,nlev)
567 real(r8),
intent(out) :: rercld(mgncol,nlev)
568 real(r8),
intent(out) :: ncai(mgncol,nlev)
569 real(r8),
intent(out) :: ncal(mgncol,nlev)
570 real(r8),
intent(out) :: qrout2(mgncol,nlev)
571 real(r8),
intent(out) :: qsout2(mgncol,nlev)
572 real(r8),
intent(out) :: nrout2(mgncol,nlev)
573 real(r8),
intent(out) :: nsout2(mgncol,nlev)
574 real(r8),
intent(out) :: drout2(mgncol,nlev)
575 real(r8),
intent(out) :: dsout2(mgncol,nlev)
576 real(r8),
intent(out) :: freqs(mgncol,nlev)
577 real(r8),
intent(out) :: freqr(mgncol,nlev)
578 real(r8),
intent(out) :: nfice(mgncol,nlev)
579 real(r8),
intent(out) :: qcrat(mgncol,nlev)
581 real(r8),
intent(out) :: prer_evap(mgncol,nlev)
602 real(r8) :: qc(mgncol,nlev)
603 real(r8) :: qi(mgncol,nlev)
604 real(r8) :: nc(mgncol,nlev)
605 real(r8) :: ni(mgncol,nlev)
606 real(r8) :: qr(mgncol,nlev)
607 real(r8) :: qs(mgncol,nlev)
608 real(r8) :: nr(mgncol,nlev)
609 real(r8) :: ns(mgncol,nlev)
617 real(r8) :: rho(mgncol,nlev)
618 real(r8) :: rhoinv(mgncol,nlev)
619 real(r8) :: dv(mgncol,nlev)
620 real(r8) :: mu(mgncol,nlev)
621 real(r8) :: sc(mgncol,nlev)
622 real(r8) :: rhof(mgncol,nlev)
625 real(r8) :: precip_frac(mgncol,nlev)
626 real(r8) :: cldm(mgncol,nlev)
627 real(r8) :: icldm(mgncol,nlev)
628 real(r8) :: lcldm(mgncol,nlev)
629 real(r8) :: qsfm(mgncol,nlev)
632 real(r8) :: qcic(mgncol,nlev)
633 real(r8) :: qiic(mgncol,nlev)
634 real(r8) :: qsic(mgncol,nlev)
635 real(r8) :: qric(mgncol,nlev)
638 real(r8) :: ncic(mgncol,nlev)
639 real(r8) :: niic(mgncol,nlev)
640 real(r8) :: nsic(mgncol,nlev)
641 real(r8) :: nric(mgncol,nlev)
643 real(r8) :: nimax(mgncol,nlev)
647 real(r8) :: lami(mgncol,nlev)
648 real(r8) :: n0i(mgncol,nlev)
650 real(r8) :: lamc(mgncol,nlev)
651 real(r8) :: pgam(mgncol,nlev)
653 real(r8) :: lams(mgncol,nlev)
654 real(r8) :: n0s(mgncol,nlev)
656 real(r8) :: lamr(mgncol,nlev)
657 real(r8) :: n0r(mgncol,nlev)
662 real(r8) :: minstsm(mgncol,nlev)
663 real(r8) :: ninstsm(mgncol,nlev)
665 real(r8) :: minstrf(mgncol,nlev)
666 real(r8) :: ninstrf(mgncol,nlev)
669 real(r8) :: vap_dep(mgncol,nlev)
671 real(r8) :: ice_sublim(mgncol,nlev)
673 real(r8) :: nnuccd(mgncol,nlev)
674 real(r8) :: mnuccd(mgncol,nlev)
676 real(r8) :: mnuccc(mgncol,nlev)
677 real(r8) :: nnuccc(mgncol,nlev)
679 real(r8) :: mnucct(mgncol,nlev)
680 real(r8) :: nnucct(mgncol,nlev)
682 real(r8) :: mnudep(mgncol,nlev)
683 real(r8) :: nnudep(mgncol,nlev)
685 real(r8) :: msacwi(mgncol,nlev)
686 real(r8) :: nsacwi(mgncol,nlev)
688 real(r8) :: prc(mgncol,nlev)
689 real(r8) :: nprc(mgncol,nlev)
690 real(r8) :: nprc1(mgncol,nlev)
692 real(r8) :: nsagg(mgncol,nlev)
694 real(r8) :: nragg(mgncol,nlev)
696 real(r8) :: psacws(mgncol,nlev)
697 real(r8) :: npsacws(mgncol,nlev)
699 real(r8) :: pracs(mgncol,nlev)
700 real(r8) :: npracs(mgncol,nlev)
702 real(r8) :: mnuccr(mgncol,nlev)
703 real(r8) :: nnuccr(mgncol,nlev)
705 real(r8) :: mnuccri(mgncol,nlev)
706 real(r8) :: nnuccri(mgncol,nlev)
708 real(r8) :: pra(mgncol,nlev)
709 real(r8) :: npra(mgncol,nlev)
711 real(r8) :: prci(mgncol,nlev)
712 real(r8) :: nprci(mgncol,nlev)
714 real(r8) :: prai(mgncol,nlev)
715 real(r8) :: nprai(mgncol,nlev)
717 real(r8) :: pre(mgncol,nlev)
719 real(r8) :: prds(mgncol,nlev)
721 real(r8) :: nsubi(mgncol,nlev)
722 real(r8) :: nsubc(mgncol,nlev)
723 real(r8) :: nsubs(mgncol,nlev)
724 real(r8) :: nsubr(mgncol,nlev)
726 real(r8) :: berg(mgncol,nlev)
727 real(r8) :: bergs(mgncol,nlev)
732 real(r8) :: uns(mgncol,nlev)
733 real(r8) :: unr(mgncol,nlev)
735 real(r8) :: arn(mgncol,nlev)
736 real(r8) :: asn(mgncol,nlev)
737 real(r8) :: acn(mgncol,nlev)
738 real(r8) :: ain(mgncol,nlev)
739 real(r8) :: ajn(mgncol,nlev)
742 real(r8) :: mi0l(mgncol)
745 real(r8) :: esl(mgncol,nlev)
746 real(r8) :: esi(mgncol,nlev)
750 real(r8) :: qvl(mgncol,nlev)
751 real(r8) :: qvi(mgncol,nlev)
755 real(r8) :: relhum(mgncol,nlev)
758 real(r8) :: fc(mgncol,nlev)
759 real(r8) :: fnc(mgncol,nlev)
760 real(r8) :: fi(mgncol,nlev)
761 real(r8) :: fni(mgncol,nlev)
763 real(r8) :: fr(mgncol,nlev)
764 real(r8) :: fnr(mgncol,nlev)
765 real(r8) :: fs(mgncol,nlev)
766 real(r8) :: fns(mgncol,nlev)
768 real(r8) :: faloutc(nlev)
769 real(r8) :: faloutnc(nlev)
770 real(r8) :: falouti(nlev)
771 real(r8) :: faloutni(nlev)
773 real(r8) :: faloutr(nlev)
774 real(r8) :: faloutnr(nlev)
775 real(r8) :: falouts(nlev)
776 real(r8) :: faloutns(nlev)
782 real(r8) :: faltndqie
783 real(r8) :: faltndqce
790 real(r8) :: rainrt(mgncol,nlev)
799 real(r8) :: tx1, tx2, tx3, tx4, tx5, tx6, tx7, grho
807 real(r8) :: dumc(mgncol,nlev)
808 real(r8) :: dumnc(mgncol,nlev)
809 real(r8) :: dumi(mgncol,nlev)
810 real(r8) :: dumni(mgncol,nlev)
811 real(r8) :: dumr(mgncol,nlev)
812 real(r8) :: dumnr(mgncol,nlev)
813 real(r8) :: dums(mgncol,nlev)
814 real(r8) :: dumns(mgncol,nlev)
817 real(r8) :: pdel_inv(mgncol,nlev)
818 real(r8) :: ts_au_loc(mgncol)
826 integer nstep, mdust, nlb, nstep_def
829 real(r8) :: irad, ifrac, tsfac
833 real(r8),
parameter :: qimax=0.010, qimin=0.005, qiinv=one/(qimax-qimin)
843 oneodt = one / deltat
845 nstep_def = max(1, nint(deltat/5))
866 if (microp_uniform)
then
875 if (qc(i,k) >= qsmall)
then
881 if (qi(i,k) >= qsmall)
then
887 cldm(i,k) = max(icldm(i,k), lcldm(i,k))
889 qsfm(i,k) = qsatfac(i,k)
896 cldm(i,k) = max(cldn(i,k), mincld)
897 lcldm(i,k) = max(liqcldf(i,k), mincld)
898 icldm(i,k) = max(icecldf(i,k), mincld)
899 qsfm(i,k) = qsatfac(i,k)
917 rho(i,k) = p(i,k) / (r*t(i,k))
918 rhoinv(i,k) = one / rho(i,k)
919 dv(i,k) = 8.794e-5_r8 * t(i,k)**1.81_r8 / p(i,k)
920 mu(i,k) = 1.496e-6_r8 * t(i,k)*sqrt(t(i,k)) / (t(i,k) + 120._r8)
921 sc(i,k) = mu(i,k) / (rho(i,k)*dv(i,k))
927 rhof(i,k) = (rhosu*rhoinv(i,k))**0.54_r8
929 arn(i,k) = ar*rhof(i,k)
930 asn(i,k) = as*rhof(i,k)
931 acn(i,k) = g*rhow/(18._r8*mu(i,k))
932 tx1 = (rhosu*rhoinv(i,k))**0.35_r8
942 accre_enhan(i,k) = accre_enhan_i
944 esl(i,k) = min(
fpvsl(t(i,k)), p(i,k))
945 qvl(i,k) = epsqs*esl(i,k) / (p(i,k)-omeps*esl(i,k))
949 if (t(i,k) >= tmelt)
then
954 esi(i,k) = min(
fpvsi(t(i,k)), p(i,k))
955 qvi(i,k) = epsqs*esi(i,k) / (p(i,k)-omeps*esi(i,k))
963 qvi(i,k) = qsfm(i,k) * qvi(i,k)
965 qvl(i,k) = qsfm(i,k) * qvl(i,k)
968 relhum(i,k) = max(zero, min(q(i,k)/max(qvl(i,k), qsmall), two))
994 mnuccctot(i,k) = zero
995 mnuccttot(i,k) = zero
996 msacwitot(i,k) = zero
997 psacwstot(i,k) = zero
1002 qcrestot(i,k) = zero
1005 qirestot(i,k) = zero
1006 mnuccrtot(i,k) = zero
1007 pracstot(i,k) = zero
1008 meltsdttot(i,k) = zero
1009 frzrdttot(i,k) = zero
1010 mnuccdtot(i,k) = zero
1030 qcsinksum_rate1ord(i,k) = zero
1034 prer_evap(i,k) = zero
1035 evapsnow(i,k) = zero
1036 am_evp_st(i,k) = zero
1038 prodsnow(i,k) = zero
1041 precip_frac(i,k) = mincld
1100 reff_rain(i,k) = zero
1101 reff_snow(i,k) = zero
1103 refl(i,k) = -9999._r8
1122 npccn(i,k) = npccnin(i,k)
1128 npccn(i,k) = max((npccnin(i,k)*lcldm(i,k)-nc(i,k))*oneodt, zero)
1150 where (qc >= qsmall)
1151 npccn = max((npccnin*lcldm-nc)*oneodt, zero)
1152 nc = max(nc + npccn*deltat, zero)
1161 if (t(i,k) < icenuct)
then
1162 ncai(i,k) = 0.005_r8*exp(0.304_r8*(273.15_r8-t(i,k))) * 1000._r8
1164 ncai(i,k) = min(ncai(i,k), 355.0e3_r8)
1165 naai(i,k) = (ncai(i,k)*rhoinv(i,k) + naai(i,k)) * half
1166 ncai(i,k) = naai(i,k)*rho(i,k)
1173 elseif (iccn == 2)
then
1176 if (t(i,k) < icenuct)
then
1177 ncai(i,k) = naai(i,k)*rho(i,k)
1187 if (t(i,k) < icenuct)
then
1188 ncai(i,k) = 0.005_r8*exp(0.304_r8*(273.15_r8-t(i,k))) * 1000._r8
1189 ncai(i,k) = min(ncai(i,k), 355.0e3_r8)
1190 naai(i,k) = ncai(i,k)*rhoinv(i,k)
1218 where (naai > zero .and. t < icenuct .and. relhum*esl/esi > 1.05_r8)
1222 nnuccd = (naai-ni/icldm)/mtime*icldm
1223 nnuccd = max(nnuccd, zero)
1229 mnuccd = nnuccd * mi0
1249 if (t(i,k) > snowmelt)
then
1250 if (qs(i,k) > zero)
then
1253 dum = -(xlf/cpp) * qs(i,k)
1254 if (t(i,k)+dum < snowmelt)
then
1255 dum = min(one, max(zero, (cpp/xlf)*(t(i,k)-snowmelt)/qs(i,k)))
1260 minstsm(i,k) = dum*qs(i,k)
1261 ninstsm(i,k) = dum*ns(i,k)
1263 dum1 = - minstsm(i,k) * (xlf*oneodt)
1264 tlat(i,k) = tlat(i,k) + dum1
1265 meltsdttot(i,k) = meltsdttot(i,k) + dum1
1267 qs(i,k) = max(qs(i,k) - minstsm(i,k), zero)
1268 ns(i,k) = max(ns(i,k) - ninstsm(i,k), zero)
1269 qr(i,k) = max(qr(i,k) + minstsm(i,k), zero)
1270 nr(i,k) = max(nr(i,k) + ninstsm(i,k), zero)
1282 if (t(i,k) < rainfrze)
then
1284 if (qr(i,k) > zero)
then
1287 dum = (xlf/cpp) * qr(i,k)
1288 if (t(i,k)+dum > rainfrze)
then
1289 dum = -(t(i,k)-rainfrze) * (cpp/xlf)
1290 dum = min(one, max(zero, dum/qr(i,k)))
1295 minstrf(i,k) = dum*qr(i,k)
1296 ninstrf(i,k) = dum*nr(i,k)
1299 dum1 = minstrf(i,k) * (xlf*oneodt)
1300 tlat(i,k) = tlat(i,k) + dum1
1301 frzrdttot(i,k) = frzrdttot(i,k) + dum1
1303 qr(i,k) = max(qr(i,k) - minstrf(i,k), zero)
1304 nr(i,k) = max(nr(i,k) - ninstrf(i,k), zero)
1305 qs(i,k) = max(qs(i,k) + minstrf(i,k), zero)
1306 ns(i,k) = max(ns(i,k) + ninstrf(i,k), zero)
1321 if (qc(i,k) >= qsmall)
then
1322 dum = one / lcldm(i,k)
1324 qcic(i,k) = min(qc(i,k)*dum, 0.05_r8)
1325 ncic(i,k) = max(nc(i,k)*dum, zero)
1329 ncic(i,k) = ncnst * rhoinv(i,k)
1336 if (qi(i,k) >= qsmall)
then
1337 dum = one / icldm(i,k)
1339 qiic(i,k) = min(qi(i,k)*dum, 0.05_r8)
1340 niic(i,k) = max(ni(i,k)*dum, zero)
1344 niic(i,k) = ninst * rhoinv(i,k)
1362 micro_vert_loop:
do k=1,nlev
1364 if (trim(micro_mg_precip_frac_method) ==
'in_cloud')
then
1367 where (qc(:,k) < qsmall .and. qi(:,k) < qsmall)
1368 precip_frac(:,k) = precip_frac(:,k-1)
1372 else if (trim(micro_mg_precip_frac_method) ==
'max_overlap')
then
1379 where (qr(:,k-1) >= qsmall .or. qs(:,k-1) >= qsmall)
1380 precip_frac(:,k) = max(precip_frac(:,k-1),precip_frac(:,k))
1396 pgam(:,k), lamc(:,k), mgncol)
1404 if (.not. do_sb_physics)
then
1406 ncic(:,k), rho(:,k), relvar(:,k), prc(:,k), nprc(:,k), nprc1(:,k), mgncol)
1412 if (precip_frac(i,k) > mincld)
then
1413 dum = one / precip_frac(i,k)
1418 qric(i,k) = min(qr(i,k)*dum, 0.05_r8)
1419 nric(i,k) = nr(i,k) * dum
1425 if(qric(i,k) < qsmall)
then
1433 nric(i,k) = max(nric(i,k),zero)
1438 lami(:,k), mgncol, n0=n0i(:,k))
1444 if (do_sb_physics)
then
1445 if (do_liq_liu)
then
1447 qric(:,k),rho(:,k),relvar(:,k),prc(:,k),nprc(:,k),nprc1(:,k),mgncol)
1450 qric(:,k),rho(:,k),relvar(:,k),prc(:,k),nprc(:,k),nprc1(:,k), mgncol)
1460 if (qiic(i,k) >= qimax)
then
1462 ts_au_loc(i) = ts_au_min
1463 elseif (qiic(i,k) <= qimin)
then
1465 ts_au_loc(i) = ts_au
1468 ts_au_loc(i) = (ts_au*(qimax-qiic(i,k)) + ts_au_min*(qiic(i,k)-qimin)) * qiinv
1473 if(do_ice_gmao)
then
1475 n0i(:,k), dcs, ts_au_loc(:), prci(:,k), nprci(:,k), mgncol)
1478 dcs, ts_au_loc(:), prci(:,k), nprci(:,k), mgncol)
1491 if (precip_frac(i,k) > mincld)
then
1492 dum = one / precip_frac(i,k)
1498 qsic(i,k) = min(qs(i,k)*dum, 0.10_r8)
1499 nsic(i,k) = ns(i,k) * dum
1503 if(qsic(i,k) < qsmall)
then
1511 nsic(i,k) = max(nsic(i,k), zero)
1520 lamr(:,k), mgncol, n0=n0r(:,k))
1523 if (lamr(i,k) >= qsmall)
then
1524 dum = arn(i,k) / lamr(i,k)**br
1525 dum1 = 9.1_r8*rhof(i,k)
1529 umr(i,k) = min(dum1, dum*gamma_br_plus4*oneo6)
1530 unr(i,k) = min(dum1, dum*gamma_br_plus1)
1542 lams(:,k), mgncol, n0=n0s(:,k))
1545 if (lams(i,k) >= qsmall)
then
1549 dum = asn(i,k) / lams(i,k)**bs
1550 dum1 = 1.2_r8*rhof(i,k)
1551 ums(i,k) = min(dum1, dum*gamma_bs_plus4*oneo6)
1552 uns(i,k) = min(dum1, dum*gamma_bs_plus1)
1561 if (.not. use_hetfrz_classnuc)
then
1567 qcic(:,k), ncic(:,k), relvar(:,k), mnuccc(:,k), nnuccc(:,k), mgncol)
1572 where (qcic(:,k) >= qsmall .and. t(:,k) < 269.15_r8)
1573 where (nnuccc(:,k)*lcldm(:,k) > nnuccd(:,k))
1575 mnuccc(:,k) = mnuccc(:,k)*(nnuccd(:,k)/(nnuccc(:,k)*lcldm(:,k)))
1576 nnuccc(:,k) = nnuccd(:,k)/lcldm(:,k)
1580 mdust =
size(rndst,3)
1582 nacon(:,k,:), pgam(:,k), lamc(:,k), qcic(:,k), ncic(:,k), &
1583 relvar(:,k), mnucct(:,k), nnucct(:,k), mgncol, mdust)
1629 call snow_self_aggregation(t(:,k), rho(:,k), asn(:,k), rhosn, qsic(:,k), nsic(:,k), &
1632 call accrete_cloud_water_snow(t(:,k), rho(:,k), asn(:,k), uns(:,k), mu(:,k), &
1633 qcic(:,k), ncic(:,k), qsic(:,k), pgam(:,k), lamc(:,k), lams(:,k), n0s(:,k), &
1634 psacws(:,k), npsacws(:,k), mgncol)
1643 call accrete_rain_snow(t(:,k), rho(:,k), umr(:,k), ums(:,k), unr(:,k), uns(:,k), &
1644 qric(:,k), qsic(:,k), lamr(:,k), n0r(:,k), lams(:,k), n0s(:,k), &
1645 pracs(:,k), npracs(:,k), mgncol)
1648 mnuccr(:,k), nnuccr(:,k), mgncol)
1650 if (do_sb_physics)
then
1652 rho(:,k), relvar(:,k), pra(:,k), npra(:,k), mgncol)
1655 ncic(:,k), relvar(:,k), accre_enhan(:,k), pra(:,k), npra(:,k), mgncol)
1661 call accrete_cloud_ice_snow(t(:,k), rho(:,k), asn(:,k), qiic(:,k), niic(:,k), &
1662 qsic(:,k), lams(:,k), n0s(:,k), prai(:,k), nprai(:,k), mgncol)
1669 dv(:,k), mu(:,k), sc(:,k), q(:,k), qvl(:,k), qvi(:,k), &
1670 lcldm(:,k), precip_frac(:,k), arn(:,k), asn(:,k), qcic(:,k), qiic(:,k), &
1671 qric(:,k), qsic(:,k), lamr(:,k), n0r(:,k), lams(:,k), n0s(:,k), &
1672 pre(:,k), prds(:,k), am_evp_st(:,k), mgncol)
1674 call bergeron_process_snow(t(:,k), rho(:,k), dv(:,k), mu(:,k), sc(:,k), &
1675 qvl(:,k), qvi(:,k), asn(:,k), qcic(:,k), qsic(:,k), lams(:,k), n0s(:,k), &
1678 bergs(:,k)=bergs(:,k)*micro_mg_berg_eff_factor
1684 icldm(:,k), rho(:,k), dv(:,k), qvl(:,k), qvi(:,k), &
1685 berg(:,k), vap_dep(:,k), ice_sublim(:,k), mgncol)
1689 ice_sublim(i,k) = max(ice_sublim(i,k), -qi(i,k)*oneodt)
1691 berg(i,k) = berg(i,k)*micro_mg_berg_eff_factor
1693 if (vap_dep(i,k) < zero .and. qi(i,k) > qsmall .and. icldm(i,k) > mincld)
then
1694 nsubi(i,k) = vap_dep(i,k) * ni(i,k) / (qi(i,k) * icldm(i,k))
1721 dum = ((prc(i,k)+pra(i,k)+mnuccc(i,k)+mnucct(i,k)+msacwi(i,k)+ &
1722 psacws(i,k)+bergs(i,k))*lcldm(i,k)+berg(i,k))*deltat
1724 if (dum > qc(i,k))
then
1725 ratio = qc(i,k) / dum * omsm
1727 prc(i,k) = ratio * prc(i,k)
1728 pra(i,k) = ratio * pra(i,k)
1729 mnuccc(i,k) = ratio * mnuccc(i,k)
1730 mnucct(i,k) = ratio * mnucct(i,k)
1731 msacwi(i,k) = ratio * msacwi(i,k)
1732 psacws(i,k) = ratio * psacws(i,k)
1733 bergs(i,k) = ratio * bergs(i,k)
1734 berg(i,k) = ratio * berg(i,k)
1743 if (qc(i,k) >= qsmall)
then
1744 vap_dep(i,k) = vap_dep(i,k)*(1._r8-qcrat(i,k))
1762 dum1 = vap_dep(i,k) + mnuccd(i,k)
1763 if (dum1 > 1.e-20_r8)
then
1764 dum = (q(i,k)-qvi(i,k))/(one + xxls_squared*qvi(i,k)/(cpp*rv*t(i,k)*t(i,k)))*oneodt
1765 dum = max(dum, zero)
1766 if (dum1 > dum)
then
1770 dum1 = mnuccd(i,k) / (vap_dep(i,k)+mnuccd(i,k))
1771 mnuccd(i,k) = dum*dum1
1772 vap_dep(i,k) = dum - mnuccd(i,k)
1783 dum = (nprc1(i,k)+npra(i,k)+nnuccc(i,k)+nnucct(i,k)+ &
1784 npsacws(i,k)-nsubc(i,k))*lcldm(i,k) * deltat
1786 if (dum > nc(i,k))
then
1787 ratio = nc(i,k) / dum * omsm
1789 nprc1(i,k) = ratio * nprc1(i,k)
1790 npra(i,k) = ratio * npra(i,k)
1791 nnuccc(i,k) = ratio * nnuccc(i,k)
1792 nnucct(i,k) = ratio * nnucct(i,k)
1793 npsacws(i,k) = ratio * npsacws(i,k)
1794 nsubc(i,k) = ratio * nsubc(i,k)
1803 if (lamr(i,k) > qsmall)
then
1804 if (one/lamr(i,k) < dcs)
then
1805 mnuccri(i,k) = mnuccr(i,k)
1806 nnuccri(i,k) = nnuccr(i,k)
1819 dum1 = -pre(i,k) + pracs(i,k) + mnuccr(i,k) + mnuccri(i,k)
1820 dum3 = dum1 * precip_frac(i,k)
1821 dum2 = (pra(i,k)+prc(i,k))*lcldm(i,k)
1822 dum = (dum3 - dum2) * deltat
1825 if (dum > qr(i,k) .and. dum1 >= qsmall)
then
1826 ratio = (qr(i,k)*oneodt + dum2) / dum3 * omsm
1827 pre(i,k) = ratio * pre(i,k)
1828 pracs(i,k) = ratio * pracs(i,k)
1829 mnuccr(i,k) = ratio * mnuccr(i,k)
1830 mnuccri(i,k) = ratio * mnuccri(i,k)
1841 if (pre(i,k) < zero)
then
1842 dum = max(-one, pre(i,k)*deltat/qr(i,k))
1843 nsubr(i,k) = dum*nr(i,k) * oneodt
1852 dum1 = (-nsubr(i,k)+npracs(i,k)+nnuccr(i,k)+nnuccri(i,k)-nragg(i,k))*precip_frac(i,k)
1853 dum2 = nprc(i,k)*lcldm(i,k)
1854 dum = (dum1 - dum2) * deltat
1856 if (dum > nr(i,k))
then
1857 ratio = (nr(i,k)*oneodt + dum2) / dum1 * omsm
1859 nragg(i,k) = ratio * nragg(i,k)
1860 npracs(i,k) = ratio * npracs(i,k)
1861 nnuccr(i,k) = ratio * nnuccr(i,k)
1862 nsubr(i,k) = ratio * nsubr(i,k)
1863 nnuccri(i,k) = ratio * nnuccri(i,k)
1875 dum1 = (prci(i,k)+prai(i,k))*icldm(i,k)-ice_sublim(i,k)
1876 dum2 = vap_dep(i,k)+berg(i,k)+mnuccd(i,k) &
1877 + (mnuccc(i,k)+mnucct(i,k)+mnudep(i,k)+msacwi(i,k))*lcldm(i,k) &
1878 + mnuccri(i,k)*precip_frac(i,k)
1879 dum = (dum1 - dum2) * deltat
1881 if (dum > qi(i,k))
then
1882 ratio = (qi(i,k)*oneodt + dum2) / dum1 * omsm
1884 prci(i,k) = ratio * prci(i,k)
1885 prai(i,k) = ratio * prai(i,k)
1886 ice_sublim(i,k) = ratio * ice_sublim(i,k)
1899 if (use_hetfrz_classnuc)
then
1900 tmpfrz = nnuccc(i,k)
1904 dum1 = (nprci(i,k)+nprai(i,k)-nsubi(i,k))*icldm(i,k)
1905 dum2 = nnuccd(i,k)+(nnucct(i,k)+tmpfrz+nnudep(i,k)+nsacwi(i,k))*lcldm(i,k) &
1906 + nnuccri(i,k)*precip_frac(i,k)
1907 dum = (dum1 - dum2) * deltat
1909 if (dum > ni(i,k))
then
1910 ratio = (ni(i,k)*oneodt + dum2) / dum1 * omsm
1912 nprci(i,k) = ratio * nprci(i,k)
1913 nprai(i,k) = ratio * nprai(i,k)
1914 nsubi(i,k) = ratio * nsubi(i,k)
1925 dum1 = - prds(i,k) * precip_frac(i,k)
1926 dum2 = (pracs(i,k)+mnuccr(i,k))*precip_frac(i,k) &
1927 + (prai(i,k)+prci(i,k))*icldm(i,k) + (bergs(i,k)+psacws(i,k))*lcldm(i,k)
1928 dum = (dum1 - dum2) * deltat
1930 if (dum > qs(i,k) .and. -prds(i,k) >= qsmall)
then
1931 ratio = (qs(i,k)*oneodt + dum2) / dum1 * omsm
1933 prds(i,k) = ratio * prds(i,k)
1946 dum1 = precip_frac(i,k)* (-nsubs(i,k)-nsagg(i,k))
1947 dum2 = nnuccr(i,k)*precip_frac(i,k) + nprci(i,k)*icldm(i,k)
1948 dum = (dum1 - dum2) * deltat
1950 if (dum > ns(i,k))
then
1951 ratio = (ns(i,k)*oneodt + dum2) / dum1 * omsm
1953 nsubs(i,k) = ratio * nsubs(i,k)
1954 nsagg(i,k) = ratio * nsagg(i,k)
1966 tx1 = pre(i,k) * precip_frac(i,k)
1967 tx2 = prds(i,k) * precip_frac(i,k)
1968 tx3 = tx1 + tx2 + ice_sublim(i,k)
1969 if (tx3 < -1.e-20_r8)
then
1971 tx4 = tx2 + ice_sublim(i,k) + vap_dep(i,k) + mnuccd(i,k)
1972 qtmp = q(i,k) - (tx1 + tx4) * deltat
1973 ttmp = t(i,k) + (tx1*xxlv + tx4*xxls) * (deltat/cpp)
1977 esn = min(
fpvsl(ttmp), p(i,k))
1978 qvn = epsqs*esn/(p(i,k)-omeps*esn) * qsfm(i,k)
1981 if (qtmp > qvn)
then
1987 tx5 = (vap_dep(i,k)+mnuccd(i,k)) * deltat
1989 ttmp = t(i,k) + tx5 * (xxls/cpp)
1993 esn = min(
fpvsl(ttmp), p(i,k))
1994 qvn = epsqs*esn / (p(i,k)-omeps*esn) * qsfm(i,k)
1996 dum = (qtmp-qvn) / (one + xxlv_squared*qvn/(cpp*rv*ttmp*ttmp))
1997 dum = min(dum, zero)
2000 if (precip_frac(i,k) > mincld)
then
2001 tx4 = oneodt / precip_frac(i,k)
2005 pre(i,k) = dum*dum1*tx4
2009 esn = min(
fpvsi(ttmp), p(i,k))
2010 qvn = epsqs*esn / (p(i,k)-omeps*esn) * qsfm(i,k)
2013 dum = (qtmp-qvn) / (one + xxls_squared*qvn/(cpp*rv*ttmp*ttmp))
2014 dum = min(dum, zero)
2017 prds(i,k) = dum*dum2*tx4
2020 dum1 = one - dum1 - dum2
2021 ice_sublim(i,k) = dum*dum1*oneodt
2043 qvlat(i,k) = qvlat(i,k) - (pre(i,k)+prds(i,k))*precip_frac(i,k)-&
2044 vap_dep(i,k)-ice_sublim(i,k)-mnuccd(i,k)-mnudep(i,k)*lcldm(i,k)
2055 tlat(i,k) = tlat(i,k) + ((pre(i,k)*precip_frac(i,k)) &
2056 *xxlv+(prds(i,k)*precip_frac(i,k)+vap_dep(i,k)+ice_sublim(i,k)+mnuccd(i,k)+mnudep(i,k)*lcldm(i,k))*xxls+ &
2057 ((bergs(i,k)+psacws(i,k)+mnuccc(i,k)+mnucct(i,k)+msacwi(i,k))*lcldm(i,k)+(mnuccr(i,k)+ &
2058 pracs(i,k)+mnuccri(i,k))*precip_frac(i,k)+berg(i,k))*xlf)
2062 qctend(i,k) = qctend(i,k) + (-pra(i,k)-prc(i,k)-mnuccc(i,k)-mnucct(i,k)-msacwi(i,k)- &
2063 psacws(i,k)-bergs(i,k))*lcldm(i,k)-berg(i,k)
2066 qitend(i,k) = qitend(i,k) + &
2067 (mnuccc(i,k)+mnucct(i,k)+mnudep(i,k)+msacwi(i,k))*lcldm(i,k)+(-prci(i,k)- &
2068 prai(i,k))*icldm(i,k)+vap_dep(i,k)+berg(i,k)+ice_sublim(i,k)+ &
2069 mnuccd(i,k)+mnuccri(i,k)*precip_frac(i,k)
2072 qrtend(i,k) = qrtend(i,k) + (pra(i,k)+prc(i,k))*lcldm(i,k)+(pre(i,k)-pracs(i,k)- &
2073 mnuccr(i,k)-mnuccri(i,k))*precip_frac(i,k)
2075 qstend(i,k) = qstend(i,k) + (prai(i,k)+prci(i,k))*icldm(i,k)+(psacws(i,k)+bergs(i,k))*lcldm(i,k) &
2076 + (prds(i,k)+pracs(i,k)+mnuccr(i,k))*precip_frac(i,k)
2079 cmeout(i,k) = vap_dep(i,k) + ice_sublim(i,k) + mnuccd(i,k)
2082 cmeitot(i,k) = vap_dep(i,k) + ice_sublim(i,k) + mnuccd(i,k)
2088 evapsnow(i,k) = -prds(i,k) * precip_frac(i,k)
2089 nevapr(i,k) = -pre(i,k) * precip_frac(i,k)
2090 prer_evap(i,k) = -pre(i,k) * precip_frac(i,k)
2094 prain(i,k) = (pra(i,k)+prc(i,k))*lcldm(i,k)+(-pracs(i,k)- &
2095 mnuccr(i,k)-mnuccri(i,k))*precip_frac(i,k)
2096 prodsnow(i,k) = (prai(i,k)+prci(i,k))*icldm(i,k)+(psacws(i,k)+bergs(i,k))*lcldm(i,k)+(&
2097 pracs(i,k)+mnuccr(i,k))*precip_frac(i,k)
2105 qcsinksum_rate1ord(i,k) = (pra(i,k)+prc(i,k)+psacws(i,k))*lcldm(i,k)
2107 qcsinksum_rate1ord(i,k) = qcsinksum_rate1ord(i,k) / max(qc(i,k),1.0e-30_r8)
2111 pratot(i,k) = pra(i,k) * lcldm(i,k)
2112 prctot(i,k) = prc(i,k) * lcldm(i,k)
2113 mnuccctot(i,k) = mnuccc(i,k) * lcldm(i,k)
2114 mnuccttot(i,k) = mnucct(i,k) * lcldm(i,k)
2115 msacwitot(i,k) = msacwi(i,k) * lcldm(i,k)
2116 psacwstot(i,k) = psacws(i,k) * lcldm(i,k)
2117 bergstot(i,k) = bergs(i,k) * lcldm(i,k)
2118 bergtot(i,k) = berg(i,k)
2119 prcitot(i,k) = prci(i,k) * icldm(i,k)
2120 praitot(i,k) = prai(i,k) * icldm(i,k)
2121 mnuccdtot(i,k) = mnuccd(i,k) * icldm(i,k)
2123 pracstot(i,k) = pracs(i,k) * precip_frac(i,k)
2124 mnuccrtot(i,k) = mnuccr(i,k) * precip_frac(i,k)
2127 nctend(i,k) = nctend(i,k) + (-nnuccc(i,k)-nnucct(i,k)-npsacws(i,k)+nsubc(i,k) &
2128 - npra(i,k)-nprc1(i,k))*lcldm(i,k)
2131 if (use_hetfrz_classnuc)
then
2132 tmpfrz = nnuccc(i,k)
2136 nitend(i,k) = nitend(i,k) + nnuccd(i,k)+ &
2137 (nnucct(i,k)+tmpfrz+nnudep(i,k)+nsacwi(i,k))*lcldm(i,k)+(nsubi(i,k)-nprci(i,k)- &
2138 nprai(i,k))*icldm(i,k)+nnuccri(i,k)*precip_frac(i,k)
2141 nstend(i,k) = nstend(i,k) + (nsubs(i,k)+nsagg(i,k)+nnuccr(i,k))*precip_frac(i,k) &
2142 + nprci(i,k)*icldm(i,k)
2144 nrtend(i,k) = nrtend(i,k) + nprc(i,k)*lcldm(i,k)+(nsubr(i,k)-npracs(i,k)-nnuccr(i,k) &
2145 - nnuccri(i,k)+nragg(i,k))*precip_frac(i,k)
2152 if (do_cldice .and. nitend(i,k) > zero .and. ni(i,k)+nitend(i,k)*deltat > nimax(i,k))
then
2153 nitend(i,k) = max(zero, (nimax(i,k)-ni(i,k))*oneodt)
2160 end do micro_vert_loop
2169 qrout(i,k) = qr(i,k)
2170 nrout(i,k) = nr(i,k) * rho(i,k)
2171 qsout(i,k) = qs(i,k)
2172 nsout(i,k) = ns(i,k) * rho(i,k)
2189 call calc_rercld(lamr, n0r, lamc, pgam, qric, qcic, ncic, rercld, mgncol, nlev)
2206 nctend(i,k) = nctend(i,k) + npccn(i,k)
2209 qstend(i,k) = qstend(i,k) + (qs(i,k)-qsn(i,k)) * oneodt
2212 nstend(i,k) = nstend(i,k) + (ns(i,k)-nsn(i,k)) * oneodt
2215 qrtend(i,k) = qrtend(i,k) + (qr(i,k)-qrn(i,k)) * oneodt
2218 nrtend(i,k) = nrtend(i,k) + (nr(i,k)-nrn(i,k)) * oneodt
2226 nevapr(i,k) = nevapr(i,k) + evapsnow(i,k)
2227 prain(i,k) = prain(i,k) + prodsnow(i,k)
2243 if (lcldm(i,k) > mincld)
then
2244 tx1 = one / lcldm(i,k)
2245 dumc(i,k) = (qc(i,k)+qctend(i,k)*deltat) * tx1
2246 dumnc(i,k) = max((nc(i,k)+nctend(i,k)*deltat)*tx1, zero)
2251 if (icldm(i,k) > mincld)
then
2252 tx1 = one / icldm(i,k)
2253 dumi(i,k) = (qi(i,k)+qitend(i,k)*deltat) * tx1
2254 dumni(i,k) = max((ni(i,k)+nitend(i,k)*deltat)*tx1, zero)
2259 if (precip_frac(i,k) > mincld)
then
2260 tx1 = one / precip_frac(i,k)
2261 dumr(i,k) = (qr(i,k)+qrtend(i,k)*deltat) * tx1
2262 dums(i,k) = (qs(i,k)+qstend(i,k)*deltat) * tx1
2264 dumnr(i,k) = max((nr(i,k)+nrtend(i,k)*deltat)*tx1, zero)
2265 dumns(i,k) = max((ns(i,k)+nstend(i,k)*deltat)*tx1, zero)
2275 dumnc(i,k) = ncnst*rhoinv(i,k)
2280 dumni(i,k) = ninst*rhoinv(i,k)
2294 pgam(:,k), lamc(:,k), mgncol)
2316 if (dumc(i,k) >= qsmall)
then
2319 vtrmc(i,k) = acn(i,k)*gamma(four+bc+pgam(i,k)) &
2320 / (tx1*gamma(pgam(i,k)+four))
2322 fc(i,k) = grho*vtrmc(i,k)
2324 fnc(i,k) = grho* acn(i,k)*gamma(pgam(i,k)+one+bc) &
2325 / (tx1*gamma(pgam(i,k)+one))
2333 if (dumi(i,k) >= qsmall)
then
2335 tx3 = one / lami(i,k)
2336 tx1 = ain(i,k) * tx3**bi
2337 tx2 = 1.2_r8*rhof(i,k)
2338 vtrmi(i,k) = min(tx1*gamma_bi_plus4*oneo6, tx2)
2340 fi(i,k) = grho * vtrmi(i,k)
2341 fni(i,k) = grho * min(tx1*gamma_bi_plus1, tx2)
2345 irad = (1.5_r8 * 1e6_r8) * tx3
2346 ifrac = min(one, max(zero, (irad-18._r8)*half))
2348 if (ifrac < one)
then
2349 tx1 = ajn(i,k) / lami(i,k)**bj
2350 vtrmi(i,k) = ifrac*vtrmi(i,k) + (one-ifrac) * min(tx1*gamma_bj_plus4*oneo6, tx2)
2352 fi(i,k) = grho*vtrmi(i,k)
2353 fni(i,k) = ifrac * fni(i,k) + (one-ifrac) * grho * min(tx1*gamma_bj_plus1, tx2)
2364 if (dumr(i,k) >= qsmall)
then
2368 tx1 = arn(i,k) / lamr(i,k)**br
2369 tx2 = 9.1_r8*rhof(i,k)
2370 umr(i,k) = min(tx1*gamma_br_plus4*oneo6, tx2)
2371 unr(i,k) = min(tx1*gamma_br_plus1, tx2)
2373 fr(i,k) = grho * umr(i,k)
2374 fnr(i,k) = grho * unr(i,k)
2385 if (dums(i,k) >= qsmall)
then
2388 tx1 = asn(i,k) / lams(i,k)**bs
2389 tx2 = 1.2_r8*rhof(i,k)
2390 ums(i,k) = min(tx1*gamma_bs_plus4*oneo6, tx2)
2391 uns(i,k) = min(tx1*gamma_bs_plus1, tx2)
2393 fs(i,k) = grho * ums(i,k)
2394 fns(i,k) = grho * uns(i,k)
2404 dumc(i,k) = qc(i,k) + qctend(i,k)*deltat
2405 dumi(i,k) = qi(i,k) + qitend(i,k)*deltat
2406 dumr(i,k) = qr(i,k) + qrtend(i,k)*deltat
2407 dums(i,k) = qs(i,k) + qstend(i,k)*deltat
2409 dumnc(i,k) = nc(i,k) + nctend(i,k)*deltat
2410 dumni(i,k) = ni(i,k) + nitend(i,k)*deltat
2411 dumnr(i,k) = nr(i,k) + nrtend(i,k)*deltat
2412 dumns(i,k) = ns(i,k) + nstend(i,k)*deltat
2414 if (dumc(i,k) < qsmall) dumnc(i,k) = zero
2415 if (dumi(i,k) < qsmall) dumni(i,k) = zero
2416 if (dumr(i,k) < qsmall) dumnr(i,k) = zero
2417 if (dums(i,k) < qsmall) dumns(i,k) = zero
2424 pdel_inv(i,k) = one / pdel(i,k)
2436 nstep = 1 + nint(max( maxval( fi(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
2437 maxval(fni(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
2439 nstep = min(nstep, nstep_def)
2458 tx7 = pdel_inv(i,k) * tx1
2459 dumi(i,k) = tx5 / (one + fi(i,k)*tx7)
2460 tx6 = (dumi(i,k)-tx5) * oneodt
2461 qitend(i,k) = qitend(i,k) + tx6
2463 dumni(i,k) = tx5 / (one + fni(i,k)*tx7)
2464 nitend(i,k) = nitend(i,k) + (dumni(i,k)-tx5) * oneodt
2467 qisedten(i,k) = qisedten(i,k) + tx6
2469 falouti(k) = fi(i,k) * dumi(i,k)
2470 faloutni(k) = fni(i,k) * dumni(i,k)
2472 iflx(i,k+1) = iflx(i,k+1) + falouti(k) * tx3
2483 if (icldm(i,k-1) > mincld)
then
2484 dum1 = max(zero, min(one, icldm(i,k)/icldm(i,k-1)))
2490 tx7 = pdel_inv(i,k) * tx1
2492 dumi(i,k) = (tx5 + falouti(k-1)*dum2) / (one + fi(i,k)*tx7)
2493 tx6 = (dumi(i,k)-tx5) * oneodt
2495 qitend(i,k) = qitend(i,k) + tx6
2497 dumni(i,k) = (tx5 + faloutni(k-1)*dum2) / (one + fni(i,k)*tx7)
2498 nitend(i,k) = nitend(i,k) + (dumni(i,k)-tx5) * oneodt
2501 qisedten(i,k) = qisedten(i,k) + tx6
2504 falouti(k) = fi(i,k) * dumi(i,k)
2505 faloutni(k) = fni(i,k) * dumni(i,k)
2507 dum2 = (one-dum1) * falouti(k-1) * pdel_inv(i,k) * tx2
2508 qvlat(i,k) = qvlat(i,k) + dum2
2509 qisevap(i,k) = qisevap(i,k) + dum2
2511 tlat(i,k) = tlat(i,k) - dum2 * xxls
2513 iflx(i,k+1) = iflx(i,k+1) + falouti(k) * tx3
2520 prect(i) = prect(i) + falouti(nlev) * (tx3*0.001_r8)
2521 preci(i) = preci(i) + falouti(nlev) * (tx3*0.001_r8)
2530 nstep = 1 + nint(max( maxval( fc(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
2531 maxval(fnc(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
2533 nstep = min(nstep, nstep_def)
2547 tx7 = pdel_inv(i,k) * tx1
2548 dumc(i,k) = tx5 / (one + fc(i,k)*tx7)
2549 tx6 = (dumc(i,k)-tx5) * oneodt
2550 qctend(i,k) = qctend(i,k) + tx6
2552 dumnc(i,k) = tx5 / (one + fnc(i,k)*tx7)
2553 nctend(i,k) = nctend(i,k) + (dumnc(i,k)-tx5) * oneodt
2557 qcsedten(i,k) = qcsedten(i,k) + tx6
2559 faloutc(k) = fc(i,k) * dumc(i,k)
2560 faloutnc(k) = fnc(i,k) * dumnc(i,k)
2562 lflx(i,k+1) = lflx(i,k+1) + faloutc(k) * tx3
2565 if (lcldm(i,k-1) > mincld)
then
2566 dum1 = max(zero, min(one, lcldm(i,k)/lcldm(i,k-1)))
2572 tx7 = pdel_inv(i,k) * tx1
2574 dumc(i,k) = (tx5 + faloutc(k-1)*dum2) / (one + fc(i,k)*tx7)
2575 tx6 = (dumc(i,k)-tx5) * oneodt
2576 qctend(i,k) = qctend(i,k) + tx6
2578 dumnc(i,k) = (tx5 + faloutnc(k-1)*dum2) / (one + fnc(i,k)*tx7)
2579 nctend(i,k) = nctend(i,k) + (dumnc(i,k)-tx5) * oneodt
2583 qcsedten(i,k) = qcsedten(i,k) + tx6
2585 faloutc(k) = fc(i,k) * dumc(i,k)
2586 faloutnc(k) = fnc(i,k) * dumnc(i,k)
2588 dum2 = (one-dum1) * faloutc(k-1) * pdel_inv(i,k) * tx2
2589 qvlat(i,k) = qvlat(i,k) + dum2
2590 qcsevap(i,k) = qcsevap(i,k) + dum2
2592 tlat(i,k) = tlat(i,k) - dum2 * xxlv
2594 lflx(i,k+1) = lflx(i,k+1) + faloutc(k) * tx3
2597 prect(i) = prect(i) + faloutc(nlev) * (tx3*0.001_r8)
2607 nstep = 1 + nint(max( maxval( fr(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
2608 maxval(fnr(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
2610 nstep = min(nstep, nstep_def)
2632 tx7 = pdel_inv(i,k) * tx1
2633 dumr(i,k) = tx5 / (one + fr(i,k)*tx7)
2634 tx6 = (dumr(i,k)-tx5) * oneodt
2635 qrtend(i,k) = qrtend(i,k) + tx6
2637 dumnr(i,k) = tx5 / (one + fnr(i,k)*tx7)
2638 nrtend(i,k) = nrtend(i,k) + (dumnr(i,k)-tx5) * oneodt
2641 qrsedten(i,k) = qrsedten(i,k) + tx6
2643 faloutr(k) = fr(i,k) * dumr(i,k)
2644 faloutnr(k) = fnr(i,k) * dumnr(i,k)
2646 rflx(i,k+1) = rflx(i,k+1) + faloutr(k) * tx3
2651 tx7 = pdel_inv(i,k) * tx1
2652 dumr(i,k) = (tx5 + faloutr(k-1)*tx7) / (one + fr(i,k)*tx7)
2653 tx6 = (dumr(i,k)-tx5) * oneodt
2654 qrtend(i,k) = qrtend(i,k) + tx6
2656 dumnr(i,k) = (tx5 + faloutnr(k-1)*tx7) / (one + fnr(i,k)*tx7)
2657 nrtend(i,k) = nrtend(i,k) + (dumnr(i,k)-tx5) * oneodt
2661 qrsedten(i,k) = qrsedten(i,k) + tx6
2663 faloutr(k) = fr(i,k) * dumr(i,k)
2664 faloutnr(k) = fnr(i,k) * dumnr(i,k)
2666 rflx(i,k+1) = rflx(i,k+1) + faloutr(k) * tx3
2669 prect(i) = prect(i) + faloutr(nlev) * (tx3*0.001_r8)
2677 nstep = 1 + nint(max( maxval( fs(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
2678 maxval(fns(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
2679 nstep = min(nstep, nstep_def)
2694 tx7 = pdel_inv(i,k) * tx1
2695 dums(i,k) = tx5 / (one + fs(i,k)*tx7)
2696 tx6 = (dums(i,k)-tx5) * oneodt
2697 qstend(i,k) = qstend(i,k) + tx6
2699 dumns(i,k) = tx5 / (one + fns(i,k)*tx7)
2700 nstend(i,k) = nstend(i,k) + (dumns(i,k)-tx5) * oneodt
2703 qssedten(i,k) = qssedten(i,k) + tx6
2705 falouts(k) = fs(i,k) * dums(i,k)
2706 faloutns(k) = fns(i,k) * dumns(i,k)
2708 sflx(i,k+1) =
sflx(i,k+1) + falouts(k) * tx3
2714 tx7 = pdel_inv(i,k) * tx1
2715 dums(i,k) = (tx5 + falouts(k-1)*tx7) / (one + fs(i,k)*tx7)
2716 tx6 = (dums(i,k)-tx5) * oneodt
2717 qstend(i,k) = qstend(i,k) + tx6
2719 dumns(i,k) = (tx5 + faloutns(k-1)*tx7) / (one + fns(i,k)*tx7)
2720 nstend(i,k) = nstend(i,k) + (dumns(i,k)-tx5) * oneodt
2723 qssedten(i,k) = qssedten(i,k) + tx6
2725 falouts(k) = fs(i,k) * dums(i,k)
2726 faloutns(k) = fns(i,k) * dumns(i,k)
2728 sflx(i,k+1) =
sflx(i,k+1) + falouts(k) * tx3
2731 prect(i) = prect(i) + falouts(nlev) * (tx3*0.001_r8)
2732 preci(i) = preci(i) + falouts(nlev) * (tx3*0.001_r8)
2748 dumc(i,k) = max(qc(i,k)+qctend(i,k)*deltat, zero)
2749 dumi(i,k) = max(qi(i,k)+qitend(i,k)*deltat, zero)
2750 dumnc(i,k) = max(nc(i,k)+nctend(i,k)*deltat, zero)
2751 dumni(i,k) = max(ni(i,k)+nitend(i,k)*deltat, zero)
2753 dumr(i,k) = max(qr(i,k)+qrtend(i,k)*deltat, zero)
2754 dumnr(i,k) = max(nr(i,k)+nrtend(i,k)*deltat, zero)
2755 dums(i,k) = max(qs(i,k)+qstend(i,k)*deltat, zero)
2756 dumns(i,k) = max(ns(i,k)+nstend(i,k)*deltat, zero)
2760 dumnc(i,k) = ncnst*rhoinv(i,k)*lcldm(i,k)
2765 dumni(i,k) = ninst*rhoinv(i,k)*icldm(i,k)
2768 if (dumc(i,k) < qsmall) dumnc(i,k) = zero
2769 if (dumi(i,k) < qsmall) dumni(i,k) = zero
2770 if (dumr(i,k) < qsmall) dumnr(i,k) = zero
2771 if (dums(i,k) < qsmall) dumns(i,k) = zero
2785 tx1 = t(i,k) + tlat(i,k)*(deltat/cpp) - snowmelt
2786 if (tx1 > zero)
then
2787 if (dums(i,k) > zero)
then
2790 dum = -(xlf/cpp) * dums(i,k)
2791 if (tx1+dum < zero)
then
2792 dum = min(one, max(zero, -tx1/dum))
2798 qstend(i,k) = qstend(i,k) - tx1*dums(i,k)
2799 nstend(i,k) = nstend(i,k) - tx1*dumns(i,k)
2800 qrtend(i,k) = qrtend(i,k) + tx1*dums(i,k)
2801 nrtend(i,k) = nrtend(i,k) + tx1*dumns(i,k)
2803 dum1 = - xlf * tx1 * dums(i,k)
2804 tlat(i,k) = tlat(i,k) + dum1
2805 meltsdttot(i,k) = meltsdttot(i,k) + dum1
2815 tx1 = t(i,k) + tlat(i,k) * (deltat/cpp) - rainfrze
2816 if (tx1 < zero)
then
2818 if (dumr(i,k) > zero)
then
2821 dum = (xlf/cpp) * dumr(i,k)
2822 if (tx1+dum > zero)
then
2823 dum = min(one, max(zero, -tx1/dum))
2828 qrtend(i,k) = qrtend(i,k) - tx2 * dumr(i,k)
2829 nrtend(i,k) = nrtend(i,k) - tx2 * dumnr(i,k)
2836 if (lamr(i,k) < one/dcs)
then
2837 qstend(i,k) = qstend(i,k) + tx2 * dumr(i,k)
2838 nstend(i,k) = nstend(i,k) + tx2 * dumnr(i,k)
2840 qitend(i,k) = qitend(i,k) + tx2 * dumr(i,k)
2841 nitend(i,k) = nitend(i,k) + tx2 * dumnr(i,k)
2844 dum1 = xlf*dum*dumr(i,k)*oneodt
2845 frzrdttot(i,k) = dum1 + frzrdttot(i,k)
2846 tlat(i,k) = dum1 + tlat(i,k)
2856 tx1 = t(i,k) + tlat(i,k) * (deltat/cpp) - tmelt
2857 if (tx1 > zero)
then
2858 if (dumi(i,k) > zero)
then
2862 dum = -dumi(i,k)*xlf/cpp
2863 if (tx1+dum < zero)
then
2864 dum = min(one, max(zero, -tx1/dum))
2870 qctend(i,k) = qctend(i,k) + tx2*dumi(i,k)
2873 melttot(i,k) = tx2*dumi(i,k)
2878 nctend(i,k) = nctend(i,k) + three*tx2*dumi(i,k)/(four*pi*5.12e-16_r8*rhow)
2880 qitend(i,k) = ((one-dum)*dumi(i,k)-qi(i,k)) * oneodt
2881 nitend(i,k) = ((one-dum)*dumni(i,k)-ni(i,k)) * oneodt
2882 tlat(i,k) = tlat(i,k) - xlf*tx2*dumi(i,k)
2896 tx1 = t(i,k) + tlat(i,k)*(deltat/cpp) - 233.15_r8
2897 if (tx1 < zero)
then
2898 if (dumc(i,k) > zero)
then
2901 dum = (xlf/cpp) * dumc(i,k)
2902 if (tx1+dum > zero)
then
2903 dum = min(one, max(zero, -tx1/dum))
2908 tx2 = dum * oneodt * dumc(i,k)
2909 qitend(i,k) = tx2 + qitend(i,k)
2915 nitend(i,k) = nitend(i,k) + tx2*(three/(four*pi*1.563e-14_r8* 500._r8))
2916 qctend(i,k) = ((one-dum)*dumc(i,k)-qc(i,k)) * oneodt
2917 nctend(i,k) = ((one-dum)*dumnc(i,k)-nc(i,k)) * oneodt
2918 tlat(i,k) = tlat(i,k) + xlf*tx2
2930 qtmp = q(i,k) + qvlat(i,k) * deltat
2931 ttmp = t(i,k) + tlat(i,k) * (deltat/cpp)
2935 esn = min(
fpvsl(ttmp), p(i,k))
2936 qvn = epsqs*esn/(p(i,k)-omeps*esn) * qsfm(i,k)
2939 if (qtmp > qvn .and. qvn > 0 .and. allow_sed_supersat)
then
2941 dum = (qtmp-qvn)/(one+xxlv_squared*qvn/(cpp*rv*ttmp*ttmp)) * oneodt
2943 cmeout(i,k) = cmeout(i,k) + dum
2945 if (ttmp > 268.15_r8)
then
2949 else if (ttmp < 238.15_r8)
then
2952 dum1 = (268.15_r8-ttmp)/30._r8
2955 tx1 = xxls*dum1 + xxlv*(one-dum1)
2956 dum = (qtmp-qvn)/(one+tx1*tx1*qvn/(cpp*rv*ttmp*ttmp)) * oneodt
2957 tx2 = dum*(one-dum1)
2958 qctend(i,k) = qctend(i,k) + tx2
2960 qitend(i,k) = qitend(i,k) + dum*dum1
2961 qirestot(i,k) = dum*dum1
2962 qvlat(i,k) = qvlat(i,k) - dum
2965 tlat(i,k) = tlat(i,k) + dum*tx1
2981 if (lcldm(i,k) > mincld)
then
2982 tx1 = one / lcldm(i,k)
2986 if (icldm(i,k) > mincld)
then
2987 tx2 = one / icldm(i,k)
2991 if (precip_frac(i,k) > mincld)
then
2992 tx3 = one / precip_frac(i,k)
2996 dumc(i,k) = max(qc(i,k)+qctend(i,k)*deltat, zero) * tx1
2997 dumi(i,k) = max(qi(i,k)+qitend(i,k)*deltat, zero) * tx2
2998 dumnc(i,k) = max(nc(i,k)+nctend(i,k)*deltat, zero) * tx1
2999 dumni(i,k) = max(ni(i,k)+nitend(i,k)*deltat, zero) * tx2
3001 dumr(i,k) = max(qr(i,k)+qrtend(i,k)*deltat, zero) * tx3
3002 dumnr(i,k) = max(nr(i,k)+nrtend(i,k)*deltat, zero) * tx3
3003 dums(i,k) = max(qs(i,k)+qstend(i,k)*deltat, zero) * tx3
3004 dumns(i,k) = max(ns(i,k)+nstend(i,k)*deltat, zero) * tx3
3008 dumnc(i,k) = ncnst * rhoinv(i,k)
3013 dumni(i,k) = ninst * rhoinv(i,k)
3019 dumc(i,k) = min(dumc(i,k), 10.e-3_r8)
3020 dumi(i,k) = min(dumi(i,k), 10.e-3_r8)
3022 dumr(i,k) = min(dumr(i,k), 10.e-3_r8)
3023 dums(i,k) = min(dums(i,k), 10.e-3_r8)
3032 if (dumi(i,k) >= qsmall)
then
3038 if (dumni(i,k) /= tx1)
then
3040 nitend(i,k) = (dumni(i,k)*icldm(i,k)-ni(i,k)) * oneodt
3043 tx1 = one / lami(i,k)
3045 effi(i,k) = (three*1.e6_r8) * tx1
3046 sadice(i,k) = two*pi*(tx1*tx1*tx1)*dumni0*rho(i,k)*1.e-2_r8
3055 deffi(i,k) = effi(i,k) * (rhoi+rhoi)/rhows
3074 if (dumc(i,k) >= qsmall)
then
3084 nctend(i,k) = (ncnst*rhoinv(i,k)*lcldm(i,k)-nc(i,k)) * oneodt
3091 pgam(i,k), lamc(i,k))
3093 if (dum /= dumnc(i,k))
then
3095 nctend(i,k) = (dumnc(i,k)*lcldm(i,k)-nc(i,k)) * oneodt
3098 effc(i,k) = (half*1.e6_r8) * (pgam(i,k)+three) / lamc(i,k)
3100 lamcrad(i,k) = lamc(i,k)
3101 pgamrad(i,k) = pgam(i,k)
3109 dumnc(i,k) = 1.e8_r8
3114 pgam(i,k), lamc(i,k))
3116 effc_fn(i,k) = (half*1.e6_r8) * (pgam(i,k)+three)/lamc(i,k)
3131 if (dumr(i,k) >= qsmall)
then
3137 if (dum /= dumnr(i,k))
then
3139 nrtend(i,k) = (dumnr(i,k)*precip_frac(i,k)-nr(i,k)) *oneodt
3149 if (dums(i,k) >= qsmall)
then
3154 lams(i,k), n0=dumns0)
3156 if (dum /= dumns(i,k))
then
3158 nstend(i,k) = (dumns(i,k)*precip_frac(i,k)-ns(i,k)) * oneodt
3161 tx1 = (two*pi*1.e-2_r8) / (lams(i,k)*lams(i,k)*lams(i,k))
3162 sadsnow(i,k) = tx1*dumns0*rho(i,k)
3173 if (qc(i,k)+qctend(i,k)*deltat < qsmall) nctend(i,k) = -nc(i,k) * oneodt
3174 if (do_cldice .and. qi(i,k)+qitend(i,k)*deltat < qsmall) nitend(i,k) = -ni(i,k) * oneodt
3175 if (qr(i,k)+qrtend(i,k)*deltat < qsmall) nrtend(i,k) = -nr(i,k) * oneodt
3176 if (qs(i,k)+qstend(i,k)*deltat < qsmall) nstend(i,k) = -ns(i,k) * oneodt
3190 qc(i,k) = qc(i,k) + qctend(i,k)*deltat
3191 qi(i,k) = qi(i,k) + qitend(i,k)*deltat
3203 if (qrout(i,k) > 1.e-7_r8 .and. nrout(i,k) > zero)
then
3204 qrout2(i,k) = qrout(i,k) * precip_frac(i,k)
3205 nrout2(i,k) = nrout(i,k) * precip_frac(i,k)
3208 drout2(i,k) =
avg_diameter(qrout(i,k), nrout(i,k), rho(i,k), rhow)
3209 freqr(i,k) = precip_frac(i,k)
3211 reff_rain(i,k) = (1.e6_r8*three) * drout2(i,k)
3217 reff_rain(i,k) = zero
3220 if (qsout(i,k) > 1.e-7_r8 .and. nsout(i,k) > zero)
then
3221 qsout2(i,k) = qsout(i,k) * precip_frac(i,k)
3222 nsout2(i,k) = nsout(i,k) * precip_frac(i,k)
3225 dsout2(i,k) =
avg_diameter(qsout(i,k), nsout(i,k), rho(i,k), rhosn)
3226 freqs(i,k) = precip_frac(i,k)
3228 dsout(i,k) = three*rhosn/rhows*dsout2(i,k)
3230 reff_snow(i,k) = (1.e6_r8*three) * dsout2(i,k)
3237 reff_snow(i,k) = zero
3251 if (qc(i,k) >= qsmall .and. (nc(i,k)+nctend(i,k)*deltat) > ten)
then
3252 tx1 = rho(i,k) / lcldm(i,k)
3253 tx2 = 1000._r8 * qc(i,k) * tx1
3254 dum = tx2 * tx2 * lcldm(i,k) &
3255 /(0.109_r8*(nc(i,k)+nctend(i,k)*deltat)*tx1*1.e-6_r8*precip_frac(i,k))
3261 if (qi(i,k) >= qsmall)
then
3263 dum1 = (qi(i,k)*rho(i,k)/icldm(i,k)*10000._r8)**(one/0.63_r8)*icldm(i,k)/precip_frac(i,k)
3268 if (qsout(i,k) >= qsmall)
then
3270 dum1 = dum1 + (qsout(i,k)*rho(i,k)*10000._r8)**(one/0.63_r8)
3273 refl(i,k) = dum + dum1
3280 if (rainrt(i,k) >= 0.001_r8)
then
3281 dum = rainrt(i,k) * rainrt(i,k)
3282 dum = log10(dum*dum*dum) + 16._r8
3286 dum = ten**(dum/ten)
3294 refl(i,k) = refl(i,k) + dum
3297 areflz(i,k) = refl(i,k) * precip_frac(i,k)
3301 if (refl(i,k) > minrefl)
then
3302 refl(i,k) = ten*log10(refl(i,k))
3304 refl(i,k) = -9999._r8
3308 if (refl(i,k) > mindbz)
then
3309 arefl(i,k) = refl(i,k) * precip_frac(i,k)
3310 frefl(i,k) = precip_frac(i,k)
3319 csrfl(i,k) = min(csmax,refl(i,k))
3322 if (csrfl(i,k) > csmin)
then
3323 acsrfl(i,k) = refl(i,k) * precip_frac(i,k)
3324 fcsrfl(i,k) = precip_frac(i,k)
3336 tx2 = qsout(i,k) + qi(i,k)
3337 tx1 = tx2 + qrout(i,k) + qc(i,k)
3338 if ( tx2 > qsmall .and. tx1 > qsmall)
then
3339 nfice(i,k) = min(tx2/tx1, one)