437subroutine micro_mg_tend ( &
438 mgncol, nlev, deltatin, &
447 relvar, accre_enhan_i, &
449 cldn, liqcldf, icecldf, qsatfac, &
450 qcsinksum_rate1ord, &
461 effc, effc_fn, effi, &
471 qgout, ngout, dgout, &
479 reff_rain, reff_snow, reff_grau, &
482 qcsevap, qisevap, qvres, &
483 cmeitot, vtrmc, vtrmi, &
488 qcsedten, qisedten, &
489 qrsedten, qssedten, &
491 mnuccctot, mnuccttot, msacwitot, &
492 psacwstot, bergstot, bergtot, &
494 qcrestot, prcitot, praitot, &
496 qirestot, mnuccrtot, mnuccritot, pracstot, &
498 meltsdttot, frzrdttot, mnuccdtot, &
500 pracgtot, psacwgtot, pgsacwtot, &
501 pgracstot, prdgtot, &
502 qmultgtot, qmultrgtot, psacrtot, &
503 npracgtot, nscngtot, ngracstot, &
504 nmultgtot, nmultrgtot, npsacwgtot, &
507 refl, arefl, areflz, &
508 frefl, csrfl, acsrfl, &
515 qgout2, ngout2, dgout2, freqg, &
519 prer_evap, xlat, xlon, lprnt, iccn, nlball)
571 integer,
intent(in) :: mgncol
572 integer,
intent(in) :: nlev
573 integer,
intent(in) :: nlball(mgncol)
574 real(r8),
intent(in) :: xlat,xlon
575 real(r8),
intent(in) :: deltatin
576 real(r8),
intent(in) :: t(mgncol,nlev)
577 real(r8),
intent(in) :: q(mgncol,nlev)
580 real(r8),
intent(in) :: qcn(mgncol,nlev)
581 real(r8),
intent(in) :: qin(mgncol,nlev)
582 real(r8),
intent(in) :: ncn(mgncol,nlev)
583 real(r8),
intent(in) :: nin(mgncol,nlev)
585 real(r8),
intent(in) :: qrn(mgncol,nlev)
586 real(r8),
intent(in) :: qsn(mgncol,nlev)
587 real(r8),
intent(in) :: nrn(mgncol,nlev)
588 real(r8),
intent(in) :: nsn(mgncol,nlev)
590 real(r8),
intent(in) :: qgr(mgncol,nlev)
591 real(r8),
intent(in) :: ngr(mgncol,nlev)
594 real(r8) :: relvar(mgncol,nlev)
595 real(r8) :: accre_enhan(mgncol,nlev)
597 real(r8),
intent(in) :: accre_enhan_i
600 real(r8),
intent(in) :: p(mgncol,nlev)
601 real(r8),
intent(in) :: pdel(mgncol,nlev)
603 real(r8),
intent(in) :: cldn(mgncol,nlev)
604 real(r8),
intent(in) :: liqcldf(mgncol,nlev)
605 real(r8),
intent(in) :: icecldf(mgncol,nlev)
606 real(r8),
intent(in) :: qsatfac(mgncol,nlev)
607 logical,
intent(in) :: lprnt
608 integer,
intent(in) :: iccn
613 real(r8),
intent(inout) :: naai(mgncol,nlev)
614 real(r8),
intent(in) :: npccnin(mgncol,nlev)
616 real(r8) :: npccn(mgncol,nlev)
620 real(r8),
intent(in) :: rndst(mgncol,nlev,10)
621 real(r8),
intent(in) :: nacon(mgncol,nlev,10)
625 real(r8),
intent(out) :: qcsinksum_rate1ord(mgncol,nlev)
627 real(r8),
intent(out) :: tlat(mgncol,nlev)
628 real(r8),
intent(out) :: qvlat(mgncol,nlev)
629 real(r8),
intent(out) :: qctend(mgncol,nlev)
630 real(r8),
intent(out) :: qitend(mgncol,nlev)
631 real(r8),
intent(out) :: nctend(mgncol,nlev)
632 real(r8),
intent(out) :: nitend(mgncol,nlev)
634 real(r8),
intent(out) :: qrtend(mgncol,nlev)
635 real(r8),
intent(out) :: qstend(mgncol,nlev)
636 real(r8),
intent(out) :: nrtend(mgncol,nlev)
637 real(r8),
intent(out) :: nstend(mgncol,nlev)
639 real(r8),
intent(out) :: qgtend(mgncol,nlev)
640 real(r8),
intent(out) :: ngtend(mgncol,nlev)
642 real(r8),
intent(out) :: effc(mgncol,nlev)
643 real(r8),
intent(out) :: effc_fn(mgncol,nlev)
644 real(r8),
intent(out) :: effi(mgncol,nlev)
645 real(r8),
intent(out) :: sadice(mgncol,nlev)
646 real(r8),
intent(out) :: sadsnow(mgncol,nlev)
647 real(r8),
intent(out) :: prect(mgncol)
648 real(r8),
intent(out) :: preci(mgncol)
649 real(r8),
intent(out) :: nevapr(mgncol,nlev)
650 real(r8),
intent(out) :: evapsnow(mgncol,nlev)
651 real(r8),
intent(out) :: am_evp_st(mgncol,nlev)
652 real(r8),
intent(out) :: prain(mgncol,nlev)
653 real(r8),
intent(out) :: prodsnow(mgncol,nlev)
654 real(r8),
intent(out) :: cmeout(mgncol,nlev)
655 real(r8),
intent(out) :: deffi(mgncol,nlev)
656 real(r8),
intent(out) :: pgamrad(mgncol,nlev)
657 real(r8),
intent(out) :: lamcrad(mgncol,nlev)
658 real(r8),
intent(out) :: qsout(mgncol,nlev)
659 real(r8),
intent(out) :: dsout(mgncol,nlev)
660 real(r8),
intent(out) :: lflx(mgncol,2:nlev+1)
661 real(r8),
intent(out) :: iflx(mgncol,2:nlev+1)
662 real(r8),
intent(out) :: rflx(mgncol,2:nlev+1)
663 real(r8),
intent(out) ::
sflx(mgncol,2:nlev+1)
665 real(r8),
intent(out) :: gflx(mgncol,2:nlev+1)
667 real(r8),
intent(out) :: qrout(mgncol,nlev)
668 real(r8),
intent(out) :: reff_rain(mgncol,nlev)
669 real(r8),
intent(out) :: reff_snow(mgncol,nlev)
671 real(r8),
intent(out) :: reff_grau(mgncol,nlev)
673 real(r8),
intent(out) :: qcsevap(mgncol,nlev)
674 real(r8),
intent(out) :: qisevap(mgncol,nlev)
675 real(r8),
intent(out) :: qvres(mgncol,nlev)
676 real(r8),
intent(out) :: cmeitot(mgncol,nlev)
677 real(r8),
intent(out) :: vtrmc(mgncol,nlev)
678 real(r8),
intent(out) :: vtrmi(mgncol,nlev)
679 real(r8),
intent(out) :: umr(mgncol,nlev)
680 real(r8),
intent(out) :: ums(mgncol,nlev)
682 real(r8),
intent(out) :: umg(mgncol,nlev)
683 real(r8),
intent(out) :: qgsedten(mgncol,nlev)
686 real(r8),
intent(out) :: qcsedten(mgncol,nlev)
687 real(r8),
intent(out) :: qisedten(mgncol,nlev)
688 real(r8),
intent(out) :: qrsedten(mgncol,nlev)
689 real(r8),
intent(out) :: qssedten(mgncol,nlev)
692 real(r8),
intent(out) :: pratot(mgncol,nlev)
693 real(r8),
intent(out) :: prctot(mgncol,nlev)
694 real(r8),
intent(out) :: mnuccctot(mgncol,nlev)
695 real(r8),
intent(out) :: mnuccttot(mgncol,nlev)
696 real(r8),
intent(out) :: msacwitot(mgncol,nlev)
697 real(r8),
intent(out) :: psacwstot(mgncol,nlev)
698 real(r8),
intent(out) :: bergstot(mgncol,nlev)
699 real(r8),
intent(out) :: bergtot(mgncol,nlev)
700 real(r8),
intent(out) :: melttot(mgncol,nlev)
701 real(r8),
intent(out) :: homotot(mgncol,nlev)
702 real(r8),
intent(out) :: qcrestot(mgncol,nlev)
703 real(r8),
intent(out) :: prcitot(mgncol,nlev)
704 real(r8),
intent(out) :: praitot(mgncol,nlev)
705 real(r8),
intent(out) :: qirestot(mgncol,nlev)
706 real(r8),
intent(out) :: mnuccrtot(mgncol,nlev)
707 real(r8),
intent(out) :: mnuccritot(mgncol,nlev)
708 real(r8),
intent(out) :: pracstot(mgncol,nlev)
709 real(r8),
intent(out) :: meltsdttot(mgncol,nlev)
710 real(r8),
intent(out) :: frzrdttot(mgncol,nlev)
711 real(r8),
intent(out) :: mnuccdtot(mgncol,nlev)
713 real(r8),
intent(out) :: pracgtot(mgncol,nlev)
714 real(r8),
intent(out) :: psacwgtot(mgncol,nlev)
715 real(r8),
intent(out) :: pgsacwtot(mgncol,nlev)
716 real(r8),
intent(out) :: pgracstot(mgncol,nlev)
717 real(r8),
intent(out) :: prdgtot(mgncol,nlev)
719 real(r8),
intent(out) :: qmultgtot(mgncol,nlev)
720 real(r8),
intent(out) :: qmultrgtot(mgncol,nlev)
721 real(r8),
intent(out) :: psacrtot(mgncol,nlev)
722 real(r8),
intent(out) :: npracgtot(mgncol,nlev)
723 real(r8),
intent(out) :: nscngtot(mgncol,nlev)
724 real(r8),
intent(out) :: ngracstot(mgncol,nlev)
725 real(r8),
intent(out) :: nmultgtot(mgncol,nlev)
726 real(r8),
intent(out) :: nmultrgtot(mgncol,nlev)
727 real(r8),
intent(out) :: npsacwgtot(mgncol,nlev)
729 real(r8),
intent(out) :: nrout(mgncol,nlev)
730 real(r8),
intent(out) :: nsout(mgncol,nlev)
731 real(r8),
intent(out) :: refl(mgncol,nlev)
732 real(r8),
intent(out) :: arefl(mgncol,nlev)
733 real(r8),
intent(out) :: areflz(mgncol,nlev)
734 real(r8),
intent(out) :: frefl(mgncol,nlev)
735 real(r8),
intent(out) :: csrfl(mgncol,nlev)
736 real(r8),
intent(out) :: acsrfl(mgncol,nlev)
737 real(r8),
intent(out) :: fcsrfl(mgncol,nlev)
738 real(r8),
intent(out) :: rercld(mgncol,nlev)
739 real(r8),
intent(out) :: ncai(mgncol,nlev)
740 real(r8),
intent(out) :: ncal(mgncol,nlev)
741 real(r8),
intent(out) :: qrout2(mgncol,nlev)
742 real(r8),
intent(out) :: qsout2(mgncol,nlev)
743 real(r8),
intent(out) :: nrout2(mgncol,nlev)
744 real(r8),
intent(out) :: nsout2(mgncol,nlev)
745 real(r8),
intent(out) :: drout2(mgncol,nlev)
746 real(r8),
intent(out) :: dsout2(mgncol,nlev)
747 real(r8),
intent(out) :: freqs(mgncol,nlev)
748 real(r8),
intent(out) :: freqr(mgncol,nlev)
749 real(r8),
intent(out) :: nfice(mgncol,nlev)
750 real(r8),
intent(out) :: qcrat(mgncol,nlev)
752 real(r8),
intent(out) :: qgout(mgncol,nlev)
753 real(r8),
intent(out) :: dgout(mgncol,nlev)
754 real(r8),
intent(out) :: ngout(mgncol,nlev)
756 real(r8),
intent(out) :: qgout2(mgncol,nlev)
757 real(r8),
intent(out) :: ngout2(mgncol,nlev)
758 real(r8),
intent(out) :: dgout2(mgncol,nlev)
759 real(r8),
intent(out) :: freqg(mgncol,nlev)
763 real(r8),
intent(out) :: prer_evap(mgncol,nlev)
784 real(r8) :: qc(mgncol,nlev)
785 real(r8) :: qi(mgncol,nlev)
786 real(r8) :: nc(mgncol,nlev)
787 real(r8) :: ni(mgncol,nlev)
788 real(r8) :: qr(mgncol,nlev)
789 real(r8) :: qs(mgncol,nlev)
790 real(r8) :: nr(mgncol,nlev)
791 real(r8) :: ns(mgncol,nlev)
793 real(r8) :: qg(mgncol,nlev)
794 real(r8) :: ng(mgncol,nlev)
805 real(r8) :: rho(mgncol,nlev)
806 real(r8) :: rhoinv(mgncol,nlev)
807 real(r8) :: dv(mgncol,nlev)
808 real(r8) :: mu(mgncol,nlev)
809 real(r8) :: sc(mgncol,nlev)
810 real(r8) :: rhof(mgncol,nlev)
813 real(r8) :: precip_frac(mgncol,nlev)
814 real(r8) :: cldm(mgncol,nlev)
815 real(r8) :: icldm(mgncol,nlev)
816 real(r8) :: lcldm(mgncol,nlev)
817 real(r8) :: qsfm(mgncol,nlev)
820 real(r8) :: qcic(mgncol,nlev)
821 real(r8) :: qiic(mgncol,nlev)
822 real(r8) :: qsic(mgncol,nlev)
823 real(r8) :: qric(mgncol,nlev)
825 real(r8) :: qgic(mgncol,nlev)
830 real(r8) :: ncic(mgncol,nlev)
831 real(r8) :: niic(mgncol,nlev)
832 real(r8) :: nsic(mgncol,nlev)
833 real(r8) :: nric(mgncol,nlev)
835 real(r8) :: ngic(mgncol,nlev)
839 real(r8) :: nimax(mgncol,nlev)
843 real(r8) :: lami(mgncol,nlev)
844 real(r8) :: n0i(mgncol,nlev)
846 real(r8) :: lamc(mgncol,nlev)
847 real(r8) :: pgam(mgncol,nlev)
849 real(r8) :: lams(mgncol,nlev)
850 real(r8) :: n0s(mgncol,nlev)
852 real(r8) :: lamr(mgncol,nlev)
853 real(r8) :: n0r(mgncol,nlev)
856 real(r8) :: lamg(mgncol,nlev)
857 real(r8) :: n0g(mgncol,nlev)
864 real(r8) :: minstsm(mgncol,nlev)
865 real(r8) :: ninstsm(mgncol,nlev)
868 real(r8) :: minstgm(mgncol,nlev)
869 real(r8) :: ninstgm(mgncol,nlev)
873 real(r8) :: minstrf(mgncol,nlev)
874 real(r8) :: ninstrf(mgncol,nlev)
877 real(r8) :: vap_dep(mgncol,nlev)
879 real(r8) :: ice_sublim(mgncol,nlev)
881 real(r8) :: nnuccd(mgncol,nlev)
882 real(r8) :: mnuccd(mgncol,nlev)
884 real(r8) :: mnuccc(mgncol,nlev)
885 real(r8) :: nnuccc(mgncol,nlev)
887 real(r8) :: mnucct(mgncol,nlev)
888 real(r8) :: nnucct(mgncol,nlev)
890 real(r8) :: mnudep(mgncol,nlev)
891 real(r8) :: nnudep(mgncol,nlev)
893 real(r8) :: msacwi(mgncol,nlev)
894 real(r8) :: nsacwi(mgncol,nlev)
896 real(r8) :: prc(mgncol,nlev)
897 real(r8) :: nprc(mgncol,nlev)
898 real(r8) :: nprc1(mgncol,nlev)
900 real(r8) :: nsagg(mgncol,nlev)
902 real(r8) :: nragg(mgncol,nlev)
904 real(r8) :: psacws(mgncol,nlev)
905 real(r8) :: npsacws(mgncol,nlev)
907 real(r8) :: pracs(mgncol,nlev)
908 real(r8) :: npracs(mgncol,nlev)
910 real(r8) :: mnuccr(mgncol,nlev)
911 real(r8) :: nnuccr(mgncol,nlev)
913 real(r8) :: mnuccri(mgncol,nlev)
914 real(r8) :: nnuccri(mgncol,nlev)
916 real(r8) :: pra(mgncol,nlev)
917 real(r8) :: npra(mgncol,nlev)
919 real(r8) :: prci(mgncol,nlev)
920 real(r8) :: nprci(mgncol,nlev)
922 real(r8) :: prai(mgncol,nlev)
923 real(r8) :: nprai(mgncol,nlev)
925 real(r8) :: pre(mgncol,nlev)
927 real(r8) :: prds(mgncol,nlev)
929 real(r8) :: nsubi(mgncol,nlev)
930 real(r8) :: nsubc(mgncol,nlev)
931 real(r8) :: nsubs(mgncol,nlev)
932 real(r8) :: nsubr(mgncol,nlev)
934 real(r8) :: berg(mgncol,nlev)
935 real(r8) :: bergs(mgncol,nlev)
939 real(r8) :: npracg(mgncol,nlev)
940 real(r8) :: nscng(mgncol,nlev)
941 real(r8) :: ngracs(mgncol,nlev)
942 real(r8) :: nmultg(mgncol,nlev)
943 real(r8) :: nmultrg(mgncol,nlev)
944 real(r8) :: npsacwg(mgncol,nlev)
946 real(r8) :: psacr(mgncol,nlev)
947 real(r8) :: pracg(mgncol,nlev)
948 real(r8) :: psacwg(mgncol,nlev)
949 real(r8) :: pgsacw(mgncol,nlev)
950 real(r8) :: pgracs(mgncol,nlev)
951 real(r8) :: prdg(mgncol,nlev)
953 real(r8) :: qmultg(mgncol,nlev)
954 real(r8) :: qmultrg(mgncol,nlev)
960 real(r8) :: uns(mgncol,nlev)
961 real(r8) :: unr(mgncol,nlev)
963 real(r8) :: ung(mgncol,nlev)
966 real(r8) :: arn(mgncol,nlev)
967 real(r8) :: asn(mgncol,nlev)
969 real(r8) :: agn(mgncol,nlev)
971 real(r8) :: acn(mgncol,nlev)
972 real(r8) :: ain(mgncol,nlev)
973 real(r8) :: ajn(mgncol,nlev)
976 real(r8) :: mi0l(mgncol)
979 real(r8) :: esl(mgncol,nlev)
980 real(r8) :: esi(mgncol,nlev)
984 real(r8) :: qvl(mgncol,nlev)
985 real(r8) :: qvi(mgncol,nlev)
989 real(r8) :: relhum(mgncol,nlev)
992 real(r8) :: fc(mgncol,nlev)
993 real(r8) :: fnc(mgncol,nlev)
994 real(r8) :: fi(mgncol,nlev)
995 real(r8) :: fni(mgncol,nlev)
998 real(r8) :: fg(mgncol,nlev)
999 real(r8) :: fng(mgncol,nlev)
1002 real(r8) :: fr(mgncol,nlev)
1003 real(r8) :: fnr(mgncol,nlev)
1004 real(r8) :: fs(mgncol,nlev)
1005 real(r8) :: fns(mgncol,nlev)
1007 real(r8) :: faloutc(nlev)
1008 real(r8) :: faloutnc(nlev)
1009 real(r8) :: falouti(nlev)
1010 real(r8) :: faloutni(nlev)
1012 real(r8) :: faloutr(nlev)
1013 real(r8) :: faloutnr(nlev)
1014 real(r8) :: falouts(nlev)
1015 real(r8) :: faloutns(nlev)
1018 real(r8) :: faltndnc
1020 real(r8) :: faltndni
1021 real(r8) :: faltndqie
1022 real(r8) :: faltndqce
1025 real(r8) :: faltndnr
1027 real(r8) :: faltndns
1030 real(r8) :: faloutg(nlev)
1031 real(r8) :: faloutng(nlev)
1033 real(r8) :: faltndng
1036 real(r8) :: rainrt(mgncol,nlev)
1047 real(r8) :: tx1, tx2, tx3, tx4, tx5, tx6, tx7, grho
1055 real(r8) :: dumc(mgncol,nlev)
1056 real(r8) :: dumnc(mgncol,nlev)
1057 real(r8) :: dumi(mgncol,nlev)
1058 real(r8) :: dumni(mgncol,nlev)
1059 real(r8) :: dumr(mgncol,nlev)
1060 real(r8) :: dumnr(mgncol,nlev)
1061 real(r8) :: dums(mgncol,nlev)
1062 real(r8) :: dumns(mgncol,nlev)
1064 real(r8) :: dumg(mgncol,nlev)
1065 real(r8) :: dumng(mgncol,nlev)
1069 real(r8) :: pdel_inv(mgncol,nlev)
1070 real(r8) :: ts_au_loc(mgncol)
1078 integer nstep, mdust, nlb, nstep_def
1082 real(r8) :: irad, ifrac
1088 real(r8),
parameter :: qimax=0.010_r8, qimin=0.005_r8, qiinv=one/(qimax-qimin)
1098 oneodt = one / deltat
1100 nstep_def = max(1, nint(deltat/5))
1124 if (microp_uniform)
then
1133 if (qc(i,k) >= qsmall)
then
1139 if (qi(i,k) >= qsmall)
then
1145 cldm(i,k) = max(icldm(i,k), lcldm(i,k))
1147 qsfm(i,k) = qsatfac(i,k)
1154 cldm(i,k) = max(cldn(i,k), mincld)
1155 lcldm(i,k) = max(liqcldf(i,k), mincld)
1156 icldm(i,k) = max(icecldf(i,k), mincld)
1157 qsfm(i,k) = qsatfac(i,k)
1178 rho(i,k) = p(i,k) / (r*t(i,k))
1179 rhoinv(i,k) = one / rho(i,k)
1180 dv(i,k) = 8.794e-5_r8 * t(i,k)**1.81_r8 / p(i,k)
1181 mu(i,k) = 1.496e-6_r8 * t(i,k)*sqrt(t(i,k)) / (t(i,k) + 120._r8)
1182 sc(i,k) = mu(i,k) / (rho(i,k)*dv(i,k))
1188 rhof(i,k) = (rhosu*rhoinv(i,k))**0.54_r8
1190 arn(i,k) = ar * rhof(i,k)
1191 asn(i,k) = as * rhof(i,k)
1193 agn(i,k) = agtmp * rhof(i,k)
1194 acn(i,k) = g*rhow/(18._r8*mu(i,k))
1195 tx1 = (rhosu*rhoinv(i,k))**0.35_r8
1205 accre_enhan(i,k) = accre_enhan_i
1207 esl(i,k) = min(
fpvsl(t(i,k)), p(i,k))
1208 qvl(i,k) = epsqs*esl(i,k) / (p(i,k)-omeps*esl(i,k))
1212 if (t(i,k) >= tmelt)
then
1217 esi(i,k) = min(
fpvsi(t(i,k)), p(i,k))
1218 qvi(i,k) = epsqs*esi(i,k) / (p(i,k)-omeps*esi(i,k))
1226 qvi(i,k) = qsfm(i,k) * qvi(i,k)
1228 qvl(i,k) = qsfm(i,k) * qvl(i,k)
1231 relhum(i,k) = max(zero, min(q(i,k)/max(qvl(i,k), qsmall), two))
1249 qcsedten(i,k) = zero
1250 qisedten(i,k) = zero
1251 qrsedten(i,k) = zero
1252 qssedten(i,k) = zero
1254 qgsedten(i,k) = zero
1260 mnuccctot(i,k) = zero
1261 mnuccttot(i,k) = zero
1262 msacwitot(i,k) = zero
1263 psacwstot(i,k) = zero
1264 bergstot(i,k) = zero
1268 qcrestot(i,k) = zero
1271 qirestot(i,k) = zero
1272 mnuccrtot(i,k) = zero
1274 mnuccritot(i,k) = zero
1277 pracstot(i,k) = zero
1278 meltsdttot(i,k) = zero
1279 frzrdttot(i,k) = zero
1280 mnuccdtot(i,k) = zero
1283 psacrtot(i,k) = zero
1284 pracgtot(i,k) = zero
1285 psacwgtot(i,k) = zero
1286 pgsacwtot(i,k) = zero
1287 pgracstot(i,k) = zero
1290 qmultgtot(i,k) = zero
1291 qmultrgtot(i,k) = zero
1292 npracgtot(i,k) = zero
1293 nscngtot(i,k) = zero
1294 ngracstot(i,k) = zero
1295 nmultgtot(i,k) = zero
1296 nmultrgtot(i,k) = zero
1297 npsacwgtot(i,k) = zero
1342 qcsinksum_rate1ord(i,k) = zero
1346 prer_evap(i,k) = zero
1347 evapsnow(i,k) = zero
1348 am_evp_st(i,k) = zero
1350 prodsnow(i,k) = zero
1353 precip_frac(i,k) = mincld
1432 reff_rain(i,k) = zero
1433 reff_snow(i,k) = zero
1435 reff_grau(i,k) = zero
1440 refl(i,k) = -9999._r8
1459 npccn(i,k) = npccnin(i,k)
1465 npccn(i,k) = max((npccnin(i,k)*lcldm(i,k)-nc(i,k))*oneodt, zero)
1500 if (qc(i,k) > qsmall .and. lcldm(i,k) >= mincld)
then
1501 npccn(i,k) = max((npccnin(i,k)*lcldm(i,k)-nc(i,k))*oneodt, zero)
1502 nc(i,k) = max(nc(i,k) + npccn(i,k)*deltat, zero)
1503 ncal(i,k) = nc(i,k) * rho(i,k) / lcldm(i,k)
1513 if (t(i,k) < icenuct)
then
1514 ncai(i,k) = 0.005_r8*exp(0.304_r8*(273.15_r8-t(i,k))) * 1000._r8
1516 ncai(i,k) = min(ncai(i,k), 355.0e3_r8)
1517 naai(i,k) = (ncai(i,k)*rhoinv(i,k) + naai(i,k)) * half
1518 ncai(i,k) = naai(i,k)*rho(i,k)
1525 elseif (iccn == 2)
then
1528 if (t(i,k) < icenuct)
then
1529 ncai(i,k) = naai(i,k)*rho(i,k)
1530 ncai(i,k) = min(ncai(i,k), 710.0e3_r8)
1531 naai(i,k) = ncai(i,k)*rhoinv(i,k)
1541 if (t(i,k) < icenuct)
then
1542 ncai(i,k) = 0.005_r8*exp(0.304_r8*(273.15_r8-t(i,k))) * 1000._r8
1543 ncai(i,k) = min(ncai(i,k), 355.0e3_r8)
1544 naai(i,k) = ncai(i,k)*rhoinv(i,k)
1567 where (naai > zero .and. t < icenuct .and. relhum*esl/esi > 1.05_r8)
1573 nnuccd = (naai-ni/icldm)/mtime*icldm
1574 nnuccd = max(nnuccd, zero)
1580 mnuccd = nnuccd * mi0
1600 if (t(i,k) > snowmelt)
then
1601 if (qs(i,k) > zero)
then
1604 dum = -(xlf/cpp) * qs(i,k)
1605 if (t(i,k)+dum < snowmelt)
then
1606 dum = min(one, max(zero, (cpp/xlf)*(t(i,k)-snowmelt)/qs(i,k)))
1611 minstsm(i,k) = dum*qs(i,k)
1612 ninstsm(i,k) = dum*ns(i,k)
1614 dum1 = - minstsm(i,k) * (xlf*oneodt)
1615 tlat(i,k) = tlat(i,k) + dum1
1616 meltsdttot(i,k) = meltsdttot(i,k) + dum1
1622 qs(i,k) = max(qs(i,k) - minstsm(i,k), zero)
1623 ns(i,k) = max(ns(i,k) - ninstsm(i,k), zero)
1624 qr(i,k) = max(qr(i,k) + minstsm(i,k), zero)
1625 nr(i,k) = max(nr(i,k) + ninstsm(i,k), zero)
1636 if (do_graupel .or. do_hail)
then
1642 if (t(i,k) > snowmelt)
then
1643 if (qg(i,k) > zero)
then
1646 dum = -(xlf/cpp) * qg(i,k)
1647 if (t(i,k)+dum < snowmelt)
then
1648 dum = max(zero, min(one, (cpp/xlf)*(t(i,k)-snowmelt)/qg(i,k)))
1653 minstgm(i,k) = dum*qg(i,k)
1654 ninstgm(i,k) = dum*ng(i,k)
1656 dum1 = - minstgm(i,k) * (xlf*oneodt)
1657 tlat(i,k) = dum1 + tlat(i,k)
1658 meltsdttot(i,k) = dum1 + meltsdttot(i,k)
1664 qg(i,k) = max(qg(i,k) - minstgm(i,k), zero)
1665 ng(i,k) = max(ng(i,k) - ninstgm(i,k), zero)
1666 qr(i,k) = max(qr(i,k) + minstgm(i,k), zero)
1667 nr(i,k) = max(nr(i,k) + ninstgm(i,k), zero)
1682 if (t(i,k) < rainfrze)
then
1684 if (qr(i,k) > zero)
then
1687 dum = (xlf/cpp) * qr(i,k)
1688 if (t(i,k)+dum > rainfrze)
then
1689 dum = -(t(i,k)-rainfrze) * (cpp/xlf)
1690 dum = min(one, max(zero, dum/qr(i,k)))
1695 minstrf(i,k) = dum*qr(i,k)
1696 ninstrf(i,k) = dum*nr(i,k)
1699 dum1 = minstrf(i,k) * (xlf*oneodt)
1700 tlat(i,k) = tlat(i,k) + dum1
1701 frzrdttot(i,k) = frzrdttot(i,k) + dum1
1703 qr(i,k) = max(qr(i,k) - minstrf(i,k), zero)
1704 nr(i,k) = max(nr(i,k) - ninstrf(i,k), zero)
1708 if(do_hail .or. do_graupel)
then
1709 qg(i,k) = max(qg(i,k) + minstrf(i,k), zero)
1710 ng(i,k) = max(ng(i,k) + ninstrf(i,k), zero)
1712 qs(i,k) = max(qs(i,k) + minstrf(i,k), zero)
1713 ns(i,k) = max(ns(i,k) + ninstrf(i,k), zero)
1738 if (qc(i,k) >= qsmall)
then
1740 dum = one / lcldm(i,k)
1742 qcic(i,k) = min(qc(i,k)*dum, 0.05_r8)
1743 ncic(i,k) = max(nc(i,k)*dum, zero)
1747 ncic(i,k) = ncnst * rhoinv(i,k)
1755 if (qi(i,k) >= qsmall)
then
1756 dum = one / icldm(i,k)
1758 qiic(i,k) = min(qi(i,k)*dum, 0.05_r8)
1759 niic(i,k) = max(ni(i,k)*dum, zero)
1763 niic(i,k) = ninst * rhoinv(i,k)
1781 micro_vert_loop:
do k=1,nlev
1783 if (trim(micro_mg_precip_frac_method) ==
'in_cloud')
then
1786 where (qc(:,k) < qsmall .and. qi(:,k) < qsmall)
1787 precip_frac(:,k) = precip_frac(:,k-1)
1791 else if (trim(micro_mg_precip_frac_method) ==
'max_overlap')
then
1802 where (qr(:,k-1) >= qsmall .or. qs(:,k-1) >= qsmall)
1803 precip_frac(:,k) = max(precip_frac(:,k-1), precip_frac(:,k))
1824 pgam(:,k), lamc(:,k), mgncol)
1837 if (.not. do_sb_physics)
then
1839 ncic(:,k), rho(:,k), relvar(:,k), prc(:,k), nprc(:,k), nprc1(:,k), mgncol)
1845 if (precip_frac(i,k) > mincld)
then
1846 dum = one / precip_frac(i,k)
1851 qric(i,k) = min(qr(i,k)*dum, 0.05_r8)
1852 nric(i,k) = nr(i,k) * dum
1858 if(qric(i,k) < qsmall)
then
1866 nric(i,k) = max(nric(i,k),zero)
1871 lami(:,k), mgncol, n0=n0i(:,k))
1874 if (do_sb_physics)
then
1875 if (do_liq_liu)
then
1877 qric(:,k),rho(:,k),relvar(:,k),prc(:,k),nprc(:,k),nprc1(:,k),mgncol)
1880 qric(:,k),rho(:,k),relvar(:,k),prc(:,k),nprc(:,k),nprc1(:,k), mgncol)
1890 if (qiic(i,k) >= qimax)
then
1892 ts_au_loc(i) = ts_au_min
1893 elseif (qiic(i,k) <= qimin)
then
1895 ts_au_loc(i) = ts_au
1898 ts_au_loc(i) = (ts_au*(qimax-qiic(i,k)) + ts_au_min*(qiic(i,k)-qimin)) * qiinv
1908 if(do_ice_gmao)
then
1910 n0i(:,k), dcs, ts_au_loc(:), prci(:,k), nprci(:,k), mgncol)
1913 dcs, ts_au_loc(:), prci(:,k), nprci(:,k), mgncol)
1926 if (precip_frac(i,k) > mincld)
then
1927 dum = one / precip_frac(i,k)
1932 qsic(i,k) = min(qs(i,k)*dum, 0.05_r8)
1933 nsic(i,k) = ns(i,k) * dum
1937 if(qsic(i,k) < qsmall)
then
1945 nsic(i,k) = max(nsic(i,k), zero)
1948 qgic(i,k) = min(qg(i,k)*dum, 0.01_r8)
1949 ngic(i,k) = ng(i,k) * dum
1952 if (qgic(i,k) < qsmall)
then
1960 ngic(i,k) = max(ngic(i,k), zero)
1970 lamr(:,k), mgncol, n0=n0r(:,k))
1973 if (lamr(i,k) >= qsmall)
then
1974 dum = arn(i,k) / lamr(i,k)**br
1975 dum1 = 9.1_r8*rhof(i,k)
1979 umr(i,k) = min(dum1, dum*gamma_br_plus4*oneo6)
1980 unr(i,k) = min(dum1, dum*gamma_br_plus1)
1992 lams(:,k), mgncol, n0=n0s(:,k))
1995 if (lams(i,k) >= qsmall)
then
1999 dum = asn(i,k) / lams(i,k)**bs
2000 dum1 = 1.2_r8*rhof(i,k)
2001 ums(i,k) = min(dum1, dum*gamma_bs_plus4*oneo6)
2002 uns(i,k) = min(dum1, dum*gamma_bs_plus1)
2010 if (do_graupel .or. do_hail)
then
2023 lamg(:,k), mgncol, n0=n0g(:,k))
2026 if (lamg(i,k) >= qsmall)
then
2030 dum = agn(i,k) / lamg(i,k)**bgtmp
2031 dum1 = 20._r8*rhof(i,k)
2032 umg(i,k) = min(dum1, dum*gamma_bg_plus4*oneo6)
2033 ung(i,k) = min(dum1, dum*gamma_bg_plus1)
2046 if (.not. use_hetfrz_classnuc)
then
2059 qcic(:,k), ncic(:,k), relvar(:,k), mnuccc(:,k), nnuccc(:,k), mgncol)
2071 where (qcic(:,k) >= qsmall .and. t(:,k) < 269.15_r8)
2072 where (nnuccc(:,k)*lcldm(:,k) > nnuccd(:,k))
2074 mnuccc(:,k) = mnuccc(:,k)*(nnuccd(:,k)/(nnuccc(:,k)*lcldm(:,k)))
2075 nnuccc(:,k) = nnuccd(:,k)/lcldm(:,k)
2083 mdust =
size(rndst,3)
2085 nacon(:,k,:), pgam(:,k), lamc(:,k), qcic(:,k), ncic(:,k), &
2086 relvar(:,k), mnucct(:,k), nnucct(:,k), mgncol, mdust)
2136 call snow_self_aggregation(t(:,k), rho(:,k), asn(:,k), rhosn, qsic(:,k), nsic(:,k), &
2139 call accrete_cloud_water_snow(t(:,k), rho(:,k), asn(:,k), uns(:,k), mu(:,k), &
2140 qcic(:,k), ncic(:,k), qsic(:,k), pgam(:,k), lamc(:,k), lams(:,k), n0s(:,k), &
2141 psacws(:,k), npsacws(:,k), mgncol)
2150 call accrete_rain_snow(t(:,k), rho(:,k), umr(:,k), ums(:,k), unr(:,k), uns(:,k), &
2151 qric(:,k), qsic(:,k), lamr(:,k), n0r(:,k), lams(:,k), n0s(:,k), &
2152 pracs(:,k), npracs(:,k), mgncol)
2155 mnuccr(:,k), nnuccr(:,k), mgncol)
2157 if (do_sb_physics)
then
2159 rho(:,k), relvar(:,k), pra(:,k), npra(:,k), mgncol)
2162 ncic(:,k), relvar(:,k), accre_enhan(:,k), pra(:,k), npra(:,k), mgncol)
2168 call accrete_cloud_ice_snow(t(:,k), rho(:,k), asn(:,k), qiic(:,k), niic(:,k), &
2169 qsic(:,k), lams(:,k), n0s(:,k), prai(:,k), nprai(:,k), mgncol)
2185 call bergeron_process_snow(t(:,k), rho(:,k), dv(:,k), mu(:,k), sc(:,k), &
2186 qvl(:,k), qvi(:,k), asn(:,k), qcic(:,k), qsic(:,k), lams(:,k), n0s(:,k), &
2193 bergs(:,k) = bergs(:,k) * micro_mg_berg_eff_factor
2199 icldm(:,k), rho(:,k), dv(:,k), qvl(:,k), qvi(:,k), &
2200 berg(:,k), vap_dep(:,k), ice_sublim(:,k), mgncol)
2209 ice_sublim(i,k) = max(ice_sublim(i,k), -qi(i,k)*oneodt)
2210 berg(i,k) = berg(i,k) * micro_mg_berg_eff_factor
2211 if (ice_sublim(i,k) < zero .and. qi(i,k) > qsmall .and. icldm(i,k) > mincld)
then
2212 nsubi(i,k) = sublim_factor * ice_sublim(i,k) * ni(i,k) / (qi(i,k) * icldm(i,k))
2230 if(do_hail .or. do_graupel)
then
2232 rho(:,k),lamr(:,k),n0r(:,k),lams(:,k),n0s(:,k), psacr(:,k), mgncol)
2235 n0g(:,k),lamg(:,k),bgtmp,agn(:,k), psacwg(:,k), npsacwg(:,k), mgncol)
2238 rho(:,k),rhosn,rhogtmp,asn(:,k),lams(:,k),n0s(:,k),deltat, &
2239 pgsacw(:,k),nscng(:,k),mgncol)
2248 umr(:,k),ung(:,k),unr(:,k),rho(:,k),n0r(:,k),lamr(:,k),n0g(:,k), &
2249 lamg(:,k), pracg(:,k),npracg(:,k),mgncol)
2262 qric(:,k),nric(:,k),nsic(:,k),n0s(:,k),lams(:,k),n0r(:,k),lamr(:,k), &
2263 deltat,pgracs(:,k),ngracs(:,k),mgncol)
2269 psacwg(:,k),pracg(:,k),qmultg(:,k),nmultg(:,k),qmultrg(:,k), &
2270 nmultrg(:,k),mgncol)
2277 dv(:,k), mu(:,k), sc(:,k), q(:,k), qvl(:,k), qvi(:,k), &
2278 lcldm(:,k), precip_frac(:,k), arn(:,k), asn(:,k), agn(:,k), bgtmp, &
2279 qcic(:,k), qiic(:,k), qric(:,k), qsic(:,k), qgic(:,k), &
2280 lamr(:,k), n0r(:,k), lams(:,k), n0s(:,k), lamg(:,k), n0g(:,k), &
2281 pre(:,k), prds(:,k), prdg(:,k), am_evp_st(:,k), mgncol)
2329 dv(:,k), mu(:,k), sc(:,k), q(:,k), qvl(:,k), qvi(:,k), &
2330 lcldm(:,k), precip_frac(:,k), arn(:,k), asn(:,k), qcic(:,k), qiic(:,k), &
2331 qric(:,k), qsic(:,k), lamr(:,k), n0r(:,k), lams(:,k), n0s(:,k), &
2332 pre(:,k), prds(:,k), am_evp_st(:,k), mgncol)
2353 dum = ( (prc(i,k) + pra(i,k) + mnuccc(i,k) + mnucct(i,k) + msacwi(i,k) &
2354 + psacws(i,k) + bergs(i,k) + qmultg(i,k) + psacwg(i,k) + pgsacw(i,k))*lcldm(i,k) &
2355 + berg(i,k) ) * deltat
2358 if (dum > qc(i,k) .and. abs(dum) > qsmall)
then
2363 ratio = qc(i,k) / dum * omsm
2365 qmultg(i,k) = ratio * qmultg(i,k)
2366 psacwg(i,k) = ratio * psacwg(i,k)
2367 pgsacw(i,k) = ratio * pgsacw(i,k)
2369 prc(i,k) = ratio * prc(i,k)
2370 pra(i,k) = ratio * pra(i,k)
2371 mnuccc(i,k) = ratio * mnuccc(i,k)
2372 mnucct(i,k) = ratio * mnucct(i,k)
2373 msacwi(i,k) = ratio * msacwi(i,k)
2374 psacws(i,k) = ratio * psacws(i,k)
2375 bergs(i,k) = ratio * bergs(i,k)
2376 berg(i,k) = ratio * berg(i,k)
2387 if (qc(i,k) >= qsmall)
then
2388 vap_dep(i,k) = vap_dep(i,k) * (one-qcrat(i,k))
2406 dum1 = vap_dep(i,k) + mnuccd(i,k)
2407 if (dum1 > 1.e-20_r8)
then
2408 dum = (q(i,k)-qvi(i,k))/(one + xxls_squared*qvi(i,k)/(cpp*rv*t(i,k)*t(i,k)))*oneodt
2409 dum = max(dum, zero)
2410 if (dum1 > dum)
then
2414 dum1 = mnuccd(i,k) / (vap_dep(i,k)+mnuccd(i,k))
2415 mnuccd(i,k) = dum*dum1
2416 vap_dep(i,k) = dum - mnuccd(i,k)
2430 dum = (nprc1(i,k) + npra(i,k) + nnuccc(i,k) + nnucct(i,k) &
2431 + npsacws(i,k) - nsubc(i,k) + npsacwg(i,k))*lcldm(i,k)*deltat
2434 if (dum > nc(i,k) .and. abs(dum) > qsmall)
then
2435 ratio = nc(i,k) / dum * omsm
2437 npsacwg(i,k) = ratio * npsacwg(i,k)
2440 nprc1(i,k) = ratio * nprc1(i,k)
2441 npra(i,k) = ratio * npra(i,k)
2442 nnuccc(i,k) = ratio * nnuccc(i,k)
2443 nnucct(i,k) = ratio * nnucct(i,k)
2444 npsacws(i,k) = ratio * npsacws(i,k)
2445 nsubc(i,k) = ratio * nsubc(i,k)
2454 if (lamr(i,k) > qsmall)
then
2455 if (one/lamr(i,k) < dcs)
then
2456 mnuccri(i,k) = mnuccr(i,k)
2457 nnuccri(i,k) = nnuccr(i,k)
2471 dum1 = - pre(i,k) + pracs(i,k) + mnuccr(i,k) + mnuccri(i,k) &
2472 + qmultrg(i,k) + pracg(i,k) + pgracs(i,k)
2473 dum3 = dum1 * precip_frac(i,k)
2474 dum2 = (pra(i,k)+prc(i,k))*lcldm(i,k)
2475 dum = (dum3 - dum2) * deltat
2479 if (dum > qr(i,k) .and. dum1 >= qsmall .and. abs(dum3) > qsmall)
then
2480 ratio = (qr(i,k)*oneodt + dum2) / dum3 * omsm
2482 qmultrg(i,k) = ratio * qmultrg(i,k)
2483 pracg(i,k) = ratio * pracg(i,k)
2484 pgracs(i,k) = ratio * pgracs(i,k)
2486 pre(i,k) = ratio * pre(i,k)
2487 pracs(i,k) = ratio * pracs(i,k)
2488 mnuccr(i,k) = ratio * mnuccr(i,k)
2489 mnuccri(i,k) = ratio * mnuccri(i,k)
2500 if (pre(i,k) < zero)
then
2501 dum = max(-one, pre(i,k)*deltat/qr(i,k))
2502 nsubr(i,k) = dum*nr(i,k) * oneodt
2515 dum1 = (-nsubr(i,k)+npracs(i,k)+nnuccr(i,k)+nnuccri(i,k)-nragg(i,k) &
2516 +npracg(i,k)+ngracs(i,k))*precip_frac(i,k)
2517 dum2 = nprc(i,k)*lcldm(i,k)
2518 dum = (dum1 - dum2) * deltat
2521 if (dum > nr(i,k) .and. abs(dum1) > qsmall)
then
2522 ratio = (nr(i,k)*oneodt + dum2) / dum1 * omsm
2525 npracg(i,k) = ratio * npracg(i,k)
2526 ngracs(i,k) = ratio * ngracs(i,k)
2528 nragg(i,k) = ratio * nragg(i,k)
2529 npracs(i,k) = ratio * npracs(i,k)
2530 nnuccr(i,k) = ratio * nnuccr(i,k)
2531 nsubr(i,k) = ratio * nsubr(i,k)
2532 nnuccri(i,k) = ratio * nnuccri(i,k)
2546 dum1 = (prci(i,k)+prai(i,k))*icldm(i,k)-ice_sublim(i,k)
2550 dum2 = vap_dep(i,k)+berg(i,k)+mnuccd(i,k) &
2551 + (mnuccc(i,k)+mnucct(i,k)+mnudep(i,k)+msacwi(i,k)+qmultg(i,k))*lcldm(i,k) &
2552 + (qmultrg(i,k)+mnuccri(i,k))*precip_frac(i,k)
2553 dum = (dum1 - dum2) * deltat
2556 if (dum > qi(i,k) .and. abs(dum1) > qsmall)
then
2557 ratio = (qi(i,k)*oneodt + dum2) / dum1 * omsm
2564 prci(i,k) = ratio * prci(i,k)
2565 prai(i,k) = ratio * prai(i,k)
2566 ice_sublim(i,k) = ratio * ice_sublim(i,k)
2579 if (use_hetfrz_classnuc)
then
2580 tmpfrz = nnuccc(i,k)
2585 dum1 = (nprci(i,k)+nprai(i,k)-nsubi(i,k))*icldm(i,k)
2588 dum2 = nnuccd(i,k)+(nnucct(i,k)+tmpfrz+nnudep(i,k)+nsacwi(i,k)+nmultg(i,k))*lcldm(i,k) &
2589 + (nmultrg(i,k)+nnuccri(i,k))*precip_frac(i,k)
2591 dum = (dum1 - dum2) * deltat
2593 if (dum > ni(i,k) .and. abs(dum1) > qsmall)
then
2594 ratio = (ni(i,k)*oneodt + dum2) / dum1 * omsm
2596 nprci(i,k) = ratio * nprci(i,k)
2597 nprai(i,k) = ratio * nprai(i,k)
2598 nsubi(i,k) = ratio * nsubi(i,k)
2610 if (do_hail .or. do_graupel)
then
2615 dum1 = (psacr(i,k) - prds(i,k)) * precip_frac(i,k)
2616 dum2 = pracs(i,k)*precip_frac(i,k) &
2617 + (prai(i,k)+prci(i,k))*icldm(i,k) + (bergs(i,k)+psacws(i,k))*lcldm(i,k)
2618 dum = (dum1 - dum2) * deltat
2619 if (dum > qs(i,k) .and. psacr(i,k)-prds(i,k) >= qsmall)
then
2620 ratio = (qs(i,k)*oneodt + dum2) / dum1 * omsm
2621 psacr(i,k) = ratio * psacr(i,k)
2622 prds(i,k) = ratio * prds(i,k)
2625 dum1 = - prds(i,k) * precip_frac(i,k)
2626 dum2 = (pracs(i,k)+mnuccr(i,k))*precip_frac(i,k) &
2627 + (prai(i,k)+prci(i,k))*icldm(i,k) + (bergs(i,k)+psacws(i,k))*lcldm(i,k)
2628 dum = (dum1 - dum2) * deltat
2629 if (dum > qs(i,k) .and. -prds(i,k) >= qsmall)
then
2630 ratio = (qs(i,k)*oneodt + dum2) / dum1 * omsm
2631 prds(i,k) = ratio * prds(i,k)
2661 if (do_hail .or. do_graupel)
then
2667 dum1 = precip_frac(i,k)* (-nsubs(i,k)-nsagg(i,k)+ngracs(i,k)) &
2668 - nscng(i,k)*lcldm(i,k)
2669 dum2 = nprci(i,k)*icldm(i,k)
2670 dum = (dum1 - dum2) * deltat
2672 if (dum > ns(i,k) .and. abs(dum1) > qsmall)
then
2673 ratio = (ns(i,k)*oneodt + dum2) / dum1 * omsm
2674 nscng(i,k) = ratio * nscng(i,k)
2675 ngracs(i,k) = ratio * ngracs(i,k)
2679 dum1 = precip_frac(i,k)* (-nsubs(i,k)-nsagg(i,k))
2680 dum2 = nnuccr(i,k)*precip_frac(i,k) + nprci(i,k)*icldm(i,k)
2681 dum = (dum1 - dum2) * deltat
2683 if (dum > ns(i,k) .and. abs(dum1) > qsmall)
then
2684 ratio = (ns(i,k)*oneodt + dum2) / dum1 * omsm
2687 nsubs(i,k) = ratio * nsubs(i,k)
2688 nsagg(i,k) = ratio * nsagg(i,k)
2694 if (do_hail .or. do_graupel)
then
2699 dum1 = -prdg(i,k) * precip_frac(i,k)
2700 dum2 = (pracg(i,k)+pgracs(i,k)+psacr(i,k)+mnuccr(i,k))*precip_frac(i,k) &
2701 + (psacwg(i,k)+pgsacw(i,k))*lcldm(i,k)
2702 dum = (dum1 - dum2) * deltat
2704 if (dum > qg(i,k) .and. abs(dum1) > qsmall)
then
2708 ratio = (qg(i,k)*oneodt + dum2) / dum1 * omsm
2710 prdg(i,k) = ratio * prdg(i,k)
2731 tx1 = pre(i,k) * precip_frac(i,k)
2732 tx2 = prds(i,k) * precip_frac(i,k)
2733 tx6 = prdg(i,k) * precip_frac(i,k)
2735 tx3 = tx1 + tx5 + ice_sublim(i,k)
2737 if (tx3 < -1.e-20_r8)
then
2739 tx4 = tx5 + ice_sublim(i,k) + vap_dep(i,k) + mnuccd(i,k)
2740 qtmp = q(i,k) - (tx1 + tx4) * deltat
2741 ttmp = t(i,k) + (tx1*xxlv + tx4*xxls) * (deltat/cpp)
2745 esn = min(
fpvsl(ttmp), p(i,k))
2746 qvn = epsqs*esn/(p(i,k)-omeps*esn) * qsfm(i,k)
2750 if (qtmp > qvn)
then
2759 tx5 = (vap_dep(i,k)+mnuccd(i,k)) * deltat
2761 ttmp = t(i,k) + tx5 * (xxls/cpp)
2765 esn = min(
fpvsl(ttmp), p(i,k))
2766 qvn = epsqs*esn / (p(i,k)-omeps*esn) * qsfm(i,k)
2769 dum = min(zero, (qtmp-qvn)/(one + xxlv_squared*qvn/(cpp*rv*ttmp*ttmp)))
2772 if (precip_frac(i,k) > mincld)
then
2773 tx4 = oneodt / precip_frac(i,k)
2777 pre(i,k) = dum*dum1*tx4
2781 esn = min(
fpvsi(ttmp), p(i,k))
2782 qvn = epsqs*esn / (p(i,k)-omeps*esn) * qsfm(i,k)
2786 dum = min(zero, (qtmp-qvn)/(one + xxls_squared*qvn/(cpp*rv*ttmp*ttmp)))
2789 prds(i,k) = dum*dum2*tx4
2791 prdg(i,k) = dum*dum3*tx4
2796 dum1 = one - dum1 - dum2 - dum3
2798 ice_sublim(i,k) = dum*dum1*oneodt
2823 qvlat(i,k) = qvlat(i,k)-(pre(i,k)+prds(i,k))*precip_frac(i,k) &
2824 -vap_dep(i,k)-ice_sublim(i,k)-mnuccd(i,k)-mnudep(i,k)*lcldm(i,k) &
2825 -prdg(i,k)*precip_frac(i,k)
2845 tlat(i,k) = tlat(i,k)+((pre(i,k)*precip_frac(i,k))*xxlv+ &
2846 ((prds(i,k)+prdg(i,k))*precip_frac(i,k)+vap_dep(i,k)+ice_sublim(i,k)+ &
2847 mnuccd(i,k)+mnudep(i,k)*lcldm(i,k))*xxls+ &
2848 ((bergs(i,k)+psacws(i,k)+mnuccc(i,k)+mnucct(i,k)+msacwi(i,k)+psacwg(i,k)+ &
2849 qmultg(i,k)+pgsacw(i,k))*lcldm(i,k)+ &
2850 (mnuccr(i,k)+pracs(i,k)+mnuccri(i,k)+pracg(i,k)+pgracs(i,k)+qmultrg(i,k))*precip_frac(i,k)+ &
2860 qctend(i,k) = qctend(i,k) + &
2861 (-pra(i,k)-prc(i,k)-mnuccc(i,k)-mnucct(i,k)-msacwi(i,k) - &
2862 psacws(i,k)-bergs(i,k)-qmultg(i,k)-psacwg(i,k)-pgsacw(i,k))*lcldm(i,k)-berg(i,k)
2870 qitend(i,k) = qitend(i,k)+ &
2871 (mnuccc(i,k)+mnucct(i,k)+mnudep(i,k)+msacwi(i,k)+qmultg(i,k)) * lcldm(i,k) &
2872 + (-prci(i,k)-prai(i,k)) * icldm(i,k) &
2873 + vap_dep(i,k)+berg(i,k)+ice_sublim(i,k)+mnuccd(i,k) &
2874 + (mnuccri(i,k)+qmultrg(i,k)) * precip_frac(i,k)
2881 qrtend(i,k) = qrtend(i,k) + (pra(i,k)+prc(i,k))*lcldm(i,k)+(pre(i,k)-pracs(i,k)- &
2882 mnuccr(i,k)-mnuccri(i,k)-qmultrg(i,k)-pracg(i,k)-pgracs(i,k))*precip_frac(i,k)
2888 if (do_hail.or.do_graupel)
then
2889 qgtend(i,k) = qgtend(i,k) + (pracg(i,k)+pgracs(i,k)+prdg(i,k)+psacr(i,k)+mnuccr(i,k))*precip_frac(i,k) &
2890 + (psacwg(i,k)+pgsacw(i,k))*lcldm(i,k)
2892 qstend(i,k) = qstend(i,k) + (prai(i,k)+prci(i,k))*icldm(i,k)+(psacws(i,k)+bergs(i,k))*lcldm(i,k) &
2893 + (prds(i,k)+pracs(i,k)-psacr(i,k))*precip_frac(i,k)
2897 qstend(i,k) = qstend(i,k) + (prai(i,k)+prci(i,k))*icldm(i,k)+(psacws(i,k)+bergs(i,k))*lcldm(i,k) &
2898 + (prds(i,k)+pracs(i,k)+mnuccr(i,k))*precip_frac(i,k)
2904 cmeout(i,k) = vap_dep(i,k) + ice_sublim(i,k) + mnuccd(i,k)
2907 cmeitot(i,k) = vap_dep(i,k) + ice_sublim(i,k) + mnuccd(i,k)
2918 evapsnow(i,k) = (-prds(i,k)-prdg(i,k)) * precip_frac(i,k)
2920 nevapr(i,k) = -pre(i,k) * precip_frac(i,k)
2921 prer_evap(i,k) = -pre(i,k) * precip_frac(i,k)
2927 prain(i,k) = (pra(i,k)+prc(i,k))*lcldm(i,k) &
2928 - (pracs(i,k)+mnuccr(i,k)+mnuccri(i,k))*precip_frac(i,k)
2929 if (do_hail .or. do_graupel)
then
2931 prodsnow(i,k) = (prai(i,k)+prci(i,k))*icldm(i,k)+(psacws(i,k)+bergs(i,k))*lcldm(i,k)+ &
2932 pracs(i,k)*precip_frac(i,k)
2934 prodsnow(i,k) = (prai(i,k)+prci(i,k))*icldm(i,k)+(psacws(i,k)+bergs(i,k))*lcldm(i,k)+ &
2935 (pracs(i,k)+mnuccr(i,k))*precip_frac(i,k)
2948 qcsinksum_rate1ord(i,k) = (pra(i,k)+prc(i,k)+psacws(i,k)+psacwg(i,k)+pgsacw(i,k))*lcldm(i,k)
2951 qcsinksum_rate1ord(i,k) = qcsinksum_rate1ord(i,k) / max(qc(i,k),1.0e-30_r8)
2955 pratot(i,k) = pra(i,k) * lcldm(i,k)
2956 prctot(i,k) = prc(i,k) * lcldm(i,k)
2957 mnuccctot(i,k) = mnuccc(i,k) * lcldm(i,k)
2958 mnuccttot(i,k) = mnucct(i,k) * lcldm(i,k)
2959 msacwitot(i,k) = msacwi(i,k) * lcldm(i,k)
2960 psacwstot(i,k) = psacws(i,k) * lcldm(i,k)
2961 bergstot(i,k) = bergs(i,k) * lcldm(i,k)
2962 bergtot(i,k) = berg(i,k)
2963 prcitot(i,k) = prci(i,k) * icldm(i,k)
2964 praitot(i,k) = prai(i,k) * icldm(i,k)
2965 mnuccdtot(i,k) = mnuccd(i,k) * icldm(i,k)
2967 pracstot(i,k) = pracs(i,k) * precip_frac(i,k)
2968 mnuccrtot(i,k) = mnuccr(i,k) * precip_frac(i,k)
2970 mnuccritot(i,k) = mnuccri(i,k) * precip_frac(i,k)
2974 psacrtot(i,k) = psacr(i,k) * precip_frac(i,k)
2975 pracgtot(i,k) = pracg(i,k) * precip_frac(i,k)
2976 psacwgtot(i,k) = psacwg(i,k) * lcldm(i,k)
2977 pgsacwtot(i,k) = pgsacw(i,k) * lcldm(i,k)
2978 pgracstot(i,k) = pgracs(i,k) * precip_frac(i,k)
2979 prdgtot(i,k) = prdg(i,k) * precip_frac(i,k)
2980 qmultgtot(i,k) = qmultg(i,k) * lcldm(i,k)
2981 qmultrgtot(i,k) = qmultrg(i,k) * precip_frac(i,k)
2982 npracgtot(i,k) = npracg(i,k) * precip_frac(i,k)
2983 nscngtot(i,k) = nscng(i,k) * lcldm(i,k)
2984 ngracstot(i,k) = ngracs(i,k) * precip_frac(i,k)
2985 nmultgtot(i,k) = nmultg(i,k) * lcldm(i,k)
2986 nmultrgtot(i,k) = nmultrg(i,k) * precip_frac(i,k)
2987 npsacwgtot(i,k) = npsacwg(i,k) * lcldm(i,k)
2994 nctend(i,k) = nctend(i,k) + (-nnuccc(i,k)-nnucct(i,k)-npsacws(i,k)+nsubc(i,k) &
2995 -npra(i,k)-nprc1(i,k)-npsacwg(i,k))*lcldm(i,k)
2998 if (use_hetfrz_classnuc)
then
2999 tmpfrz = nnuccc(i,k)
3007 nitend(i,k) = nitend(i,k) + nnuccd(i,k) &
3008 + (nnucct(i,k)+tmpfrz+nnudep(i,k)+nsacwi(i,k)+nmultg(i,k))*lcldm(i,k) &
3009 + (nsubi(i,k)-nprci(i,k)-nprai(i,k))*icldm(i,k) &
3010 + (nnuccri(i,k)+nmultrg(i,k))*precip_frac(i,k)
3013 if(do_graupel.or.do_hail)
then
3017 nstend(i,k) = nstend(i,k) + (nsubs(i,k)+nsagg(i,k)-ngracs(i,k))*precip_frac(i,k) &
3018 + nprci(i,k)*icldm(i,k)-nscng(i,k)*lcldm(i,k)
3020 ngtend(i,k) = ngtend(i,k) + nscng(i,k)*lcldm(i,k) &
3021 + (ngracs(i,k)+nnuccr(i,k))*precip_frac(i,k)
3025 nstend(i,k) = nstend(i,k) + (nsubs(i,k)+nsagg(i,k)+nnuccr(i,k))*precip_frac(i,k) &
3026 + nprci(i,k)*icldm(i,k)
3033 nrtend(i,k) = nrtend(i,k)+ nprc(i,k)*lcldm(i,k) &
3034 + (nsubr(i,k)-npracs(i,k)-nnuccr(i,k) &
3035 -nnuccri(i,k)+nragg(i,k)-npracg(i,k)-ngracs(i,k))*precip_frac(i,k)
3043 if (do_cldice .and. nitend(i,k) > zero .and. ni(i,k)+nitend(i,k)*deltat > nimax(i,k))
then
3044 nitend(i,k) = max(zero, (nimax(i,k)-ni(i,k))*oneodt)
3051 end do micro_vert_loop
3060 qrout(i,k) = qr(i,k)
3061 nrout(i,k) = nr(i,k) * rho(i,k)
3062 qsout(i,k) = qs(i,k)
3063 nsout(i,k) = ns(i,k) * rho(i,k)
3065 qgout(i,k) = qg(i,k)
3066 ngout(i,k) = ng(i,k) * rho(i,k)
3084 call calc_rercld(lamr, n0r, lamc, pgam, qric, qcic, ncic, rercld, mgncol, nlev)
3101 nctend(i,k) = nctend(i,k) + npccn(i,k)
3104 qstend(i,k) = qstend(i,k) + (qs(i,k)-qsn(i,k)) * oneodt
3107 nstend(i,k) = nstend(i,k) + (ns(i,k)-nsn(i,k)) * oneodt
3110 qrtend(i,k) = qrtend(i,k) + (qr(i,k)-qrn(i,k)) * oneodt
3113 nrtend(i,k) = nrtend(i,k) + (nr(i,k)-nrn(i,k)) * oneodt
3117 qgtend(i,k) = qgtend(i,k) + (qg(i,k)-qgr(i,k)) * oneodt
3130 ngtend(i,k) = ngtend(i,k) + (ng(i,k)-ngr(i,k)) * oneodt
3139 nevapr(i,k) = nevapr(i,k) + evapsnow(i,k)
3140 prain(i,k) = prain(i,k) + prodsnow(i,k)
3158 if (lcldm(i,k) > mincld)
then
3159 tx1 = one / lcldm(i,k)
3160 dumc(i,k) = (qc(i,k)+qctend(i,k)*deltat) * tx1
3161 dumnc(i,k) = max((nc(i,k)+nctend(i,k)*deltat)*tx1, zero)
3164 dumnc(i,k) = ncnst*rhoinv(i,k)
3170 if (icldm(i,k) > mincld)
then
3171 tx1 = one / icldm(i,k)
3172 dumi(i,k) = (qi(i,k)+qitend(i,k)*deltat) * tx1
3173 dumni(i,k) = max((ni(i,k)+nitend(i,k)*deltat)*tx1, zero)
3176 dumni(i,k) = ninst*rhoinv(i,k)
3182 if (precip_frac(i,k) > mincld)
then
3183 tx1 = one / precip_frac(i,k)
3184 dumr(i,k) = (qr(i,k)+qrtend(i,k)*deltat) * tx1
3185 dums(i,k) = (qs(i,k)+qstend(i,k)*deltat) * tx1
3187 dumnr(i,k) = max((nr(i,k)+nrtend(i,k)*deltat)*tx1, zero)
3188 dumns(i,k) = max((ns(i,k)+nstend(i,k)*deltat)*tx1, zero)
3191 dumg(i,k) = (qg(i,k)+qgtend(i,k)*deltat) * tx1
3193 if (dumg(i,k) > 0.01_r8)
then
3194 tx2 = dumg(i,k) - 0.01_r8
3196 dums(i,k) = dums(i,k) + tx2
3197 qstend(i,k) = (dums(i,k)*precip_frac(i,k) - qs(i,k)) * oneodt
3198 qgtend(i,k) = (dumg(i,k)*precip_frac(i,k) - qg(i,k)) * oneodt
3202 dumng(i,k) = max((ng(i,k)+ngtend(i,k)*deltat)*tx1, zero)
3205 dumng(i,k) = ngnst*rhoinv(i,k)
3229 pgam(:,k), lamc(:,k), mgncol)
3241 if (do_graupel .or. do_hail)
then
3255 if (dumc(i,k) >= qsmall)
then
3258 vtrmc(i,k) = acn(i,k)*gamma(pgam(i,k)+four+bc) &
3259 / (tx1*gamma(pgam(i,k)+four))
3261 fc(i,k) = grho * vtrmc(i,k)
3262 fnc(i,k) = grho * acn(i,k)*gamma(pgam(i,k)+one+bc) &
3263 / (tx1*gamma(pgam(i,k)+one))
3271 if (dumi(i,k) >= qsmall)
then
3273 tx3 = one / lami(i,k)
3274 tx1 = ain(i,k) * tx3**bi
3275 tx2 = 1.2_r8*rhof(i,k)
3276 vtrmi(i,k) = min(tx1*gamma_bi_plus4*oneo6, tx2)
3278 fi(i,k) = grho * vtrmi(i,k)
3279 fni(i,k) = grho * min(tx1*gamma_bi_plus1, tx2)
3283 irad = (1.5_r8 * 1e6_r8) * tx3
3284 ifrac = min(one, max(zero, (irad-18._r8)*half))
3286 if (ifrac < one)
then
3287 tx1 = ajn(i,k) / lami(i,k)**bj
3288 vtrmi(i,k) = ifrac*vtrmi(i,k) + (one-ifrac) * min(tx1*gamma_bj_plus4*oneo6, tx2)
3290 fi(i,k) = grho * vtrmi(i,k)
3291 fni(i,k) = ifrac * fni(i,k) + (one-ifrac) * grho * min(tx1*gamma_bj_plus1, tx2)
3301 if (dumr(i,k) >= qsmall)
then
3305 tx1 = arn(i,k) / lamr(i,k)**br
3306 tx2 = 9.1_r8*rhof(i,k)
3307 umr(i,k) = min(tx1*gamma_br_plus4*oneo6, tx2)
3308 unr(i,k) = min(tx1*gamma_br_plus1, tx2)
3310 fr(i,k) = grho * umr(i,k)
3311 fnr(i,k) = grho * unr(i,k)
3321 if (dums(i,k) >= qsmall)
then
3324 tx1 = asn(i,k) / lams(i,k)**bs
3325 tx2 = 1.2_r8*rhof(i,k)
3326 ums(i,k) = min(tx1*gamma_bs_plus4*oneo6, tx2)
3327 uns(i,k) = min(tx1*gamma_bs_plus1, tx2)
3329 fs(i,k) = grho * ums(i,k)
3330 fns(i,k) = grho * uns(i,k)
3337 if (do_graupel .or. do_hail)
then
3343 if (dumg(i,k) >= qsmall)
then
3346 tx1 = agn(i,k) / lamg(i,k)**bgtmp
3347 tx2 = 20._r8 * rhof(i,k)
3348 umg(i,k) = min(tx1*gamma_bg_plus4*oneo6, tx2)
3349 ung(i,k) = min(tx1*gamma_bg_plus1, tx2)
3351 fg(i,k) = grho * umg(i,k)
3352 fng(i,k) = grho * ung(i,k)
3363 dumc(i,k) = qc(i,k) + qctend(i,k)*deltat
3364 dumi(i,k) = qi(i,k) + qitend(i,k)*deltat
3365 dumr(i,k) = qr(i,k) + qrtend(i,k)*deltat
3366 dums(i,k) = qs(i,k) + qstend(i,k)*deltat
3368 dumnc(i,k) = nc(i,k) + nctend(i,k)*deltat
3369 dumni(i,k) = ni(i,k) + nitend(i,k)*deltat
3370 dumnr(i,k) = nr(i,k) + nrtend(i,k)*deltat
3371 dumns(i,k) = ns(i,k) + nstend(i,k)*deltat
3373 dumg(i,k) = qg(i,k) + qgtend(i,k)*deltat
3374 dumng(i,k) = ng(i,k) + ngtend(i,k)*deltat
3377 if (dumc(i,k) < qsmall) dumnc(i,k) = zero
3378 if (dumi(i,k) < qsmall) dumni(i,k) = zero
3379 if (dumr(i,k) < qsmall) dumnr(i,k) = zero
3380 if (dums(i,k) < qsmall) dumns(i,k) = zero
3381 if (dumg(i,k) < qsmall) dumng(i,k) = zero
3388 pdel_inv(i,k) = one / pdel(i,k)
3400 nstep = 1 + nint(max( maxval( fi(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
3401 maxval(fni(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
3402 nstep = min(nstep, nstep_def)
3421 tx7 = pdel_inv(i,k) * tx1
3422 dumi(i,k) = tx5 / (one + fi(i,k)*tx7)
3423 tx6 = (dumi(i,k)-tx5) * oneodt
3424 qitend(i,k) = qitend(i,k) + tx6
3426 dumni(i,k) = tx5 / (one + fni(i,k)*tx7)
3427 nitend(i,k) = nitend(i,k) + (dumni(i,k)-tx5) * oneodt
3430 qisedten(i,k) = qisedten(i,k) + tx6
3432 falouti(k) = fi(i,k) * dumi(i,k)
3433 faloutni(k) = fni(i,k) * dumni(i,k)
3435 iflx(i,k+1) = iflx(i,k+1) + falouti(k) * tx3
3446 if (icldm(i,k-1) > mincld)
then
3447 dum1 = max(zero, min(one, icldm(i,k)/icldm(i,k-1)))
3453 tx7 = pdel_inv(i,k) * tx1
3455 dumi(i,k) = (tx5 + falouti(k-1)*dum2) / (one + fi(i,k)*tx7)
3456 tx6 = (dumi(i,k)-tx5) * oneodt
3458 qitend(i,k) = qitend(i,k) + tx6
3460 dumni(i,k) = (tx5 + faloutni(k-1)*dum2) / (one + fni(i,k)*tx7)
3461 nitend(i,k) = nitend(i,k) + (dumni(i,k)-tx5) * oneodt
3464 qisedten(i,k) = qisedten(i,k) + tx6
3467 falouti(k) = fi(i,k) * dumi(i,k)
3468 faloutni(k) = fni(i,k) * dumni(i,k)
3470 dum2 = (one-dum1) * falouti(k-1) * pdel_inv(i,k) * tx2
3471 qvlat(i,k) = qvlat(i,k) + dum2
3472 qisevap(i,k) = qisevap(i,k) + dum2
3474 tlat(i,k) = tlat(i,k) - dum2 * xxls
3476 iflx(i,k+1) = iflx(i,k+1) + falouti(k) * tx3
3483 prect(i) = prect(i) + falouti(nlev) * (tx3*0.001_r8)
3484 preci(i) = preci(i) + falouti(nlev) * (tx3*0.001_r8)
3493 nstep = 1 + nint(max( maxval( fc(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
3494 maxval(fnc(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
3495 nstep = min(nstep, nstep_def)
3509 tx7 = pdel_inv(i,k) * tx1
3510 dumc(i,k) = tx5 / (one + fc(i,k)*tx7)
3511 tx6 = (dumc(i,k)-tx5) * oneodt
3512 qctend(i,k) = qctend(i,k) + tx6
3514 dumnc(i,k) = tx5 / (one + fnc(i,k)*tx7)
3515 nctend(i,k) = nctend(i,k) + (dumnc(i,k)-tx5) * oneodt
3519 qcsedten(i,k) = qcsedten(i,k) + tx6
3521 faloutc(k) = fc(i,k) * dumc(i,k)
3522 faloutnc(k) = fnc(i,k) * dumnc(i,k)
3524 lflx(i,k+1) = lflx(i,k+1) + faloutc(k) * tx3
3527 if (lcldm(i,k-1) > mincld)
then
3528 dum1 = max(zero, min(one, lcldm(i,k)/lcldm(i,k-1)))
3534 tx7 = pdel_inv(i,k) * tx1
3536 dumc(i,k) = (tx5 + faloutc(k-1)*dum2) / (one + fc(i,k)*tx7)
3537 tx6 = (dumc(i,k)-tx5) * oneodt
3538 qctend(i,k) = qctend(i,k) + tx6
3540 dumnc(i,k) = (tx5 + faloutnc(k-1)*dum2) / (one + fnc(i,k)*tx7)
3541 nctend(i,k) = nctend(i,k) + (dumnc(i,k)-tx5) * oneodt
3545 qcsedten(i,k) = qcsedten(i,k) + tx6
3547 faloutc(k) = fc(i,k) * dumc(i,k)
3548 faloutnc(k) = fnc(i,k) * dumnc(i,k)
3550 dum2 = (one-dum1) * faloutc(k-1) * pdel_inv(i,k) * tx2
3551 qvlat(i,k) = qvlat(i,k) + dum2
3552 qcsevap(i,k) = qcsevap(i,k) + dum2
3554 tlat(i,k) = tlat(i,k) - dum2 * xxlv
3556 lflx(i,k+1) = lflx(i,k+1) + faloutc(k) * tx3
3559 prect(i) = prect(i) + faloutc(nlev) * (tx3*0.001_r8)
3569 nstep = 1 + nint(max( maxval( fr(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
3570 maxval(fnr(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
3572 nstep = min(nstep, nstep_def)
3594 tx7 = pdel_inv(i,k) * tx1
3595 dumr(i,k) = tx5 / (one + fr(i,k)*tx7)
3596 tx6 = (dumr(i,k)-tx5) * oneodt
3597 qrtend(i,k) = qrtend(i,k) + tx6
3599 dumnr(i,k) = tx5 / (one + fnr(i,k)*tx7)
3600 nrtend(i,k) = nrtend(i,k) + (dumnr(i,k)-tx5) * oneodt
3603 qrsedten(i,k) = qrsedten(i,k) + tx6
3605 faloutr(k) = fr(i,k) * dumr(i,k)
3606 faloutnr(k) = fnr(i,k) * dumnr(i,k)
3608 rflx(i,k+1) = rflx(i,k+1) + faloutr(k) * tx3
3613 tx7 = pdel_inv(i,k) * tx1
3614 dumr(i,k) = (tx5 + faloutr(k-1)*tx7) / (one + fr(i,k)*tx7)
3615 tx6 = (dumr(i,k)-tx5) * oneodt
3616 qrtend(i,k) = qrtend(i,k) + tx6
3618 dumnr(i,k) = (tx5 + faloutnr(k-1)*tx7) / (one + fnr(i,k)*tx7)
3619 nrtend(i,k) = nrtend(i,k) + (dumnr(i,k)-tx5) * oneodt
3621 qrsedten(i,k) = qrsedten(i,k) + tx6
3623 faloutr(k) = fr(i,k) * dumr(i,k)
3624 faloutnr(k) = fnr(i,k) * dumnr(i,k)
3626 rflx(i,k+1) = rflx(i,k+1) + faloutr(k) * tx3
3629 prect(i) = prect(i) + faloutr(nlev) * (tx3*0.001_r8)
3638 nstep = 1 + nint(max( maxval( fs(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
3639 maxval(fns(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
3640 nstep = min(nstep, nstep_def)
3656 tx7 = pdel_inv(i,k) * tx1
3657 dums(i,k) = tx5 / (one + fs(i,k)*tx7)
3658 tx6 = (dums(i,k)-tx5) * oneodt
3659 qstend(i,k) = qstend(i,k) + tx6
3661 dumns(i,k) = tx5 / (one + fns(i,k)*tx7)
3662 nstend(i,k) = nstend(i,k) + (dumns(i,k)-tx5) * oneodt
3665 qssedten(i,k) = qssedten(i,k) + tx6
3667 falouts(k) = fs(i,k) * dums(i,k)
3668 faloutns(k) = fns(i,k) * dumns(i,k)
3670 sflx(i,k+1) =
sflx(i,k+1) + falouts(k) * tx3
3676 tx7 = pdel_inv(i,k) * tx1
3677 dums(i,k) = (tx5 + falouts(k-1)*tx7) / (one + fs(i,k)*tx7)
3678 tx6 = (dums(i,k)-tx5) * oneodt
3679 qstend(i,k) = qstend(i,k) + tx6
3681 dumns(i,k) = (tx5 + faloutns(k-1)*tx7) / (one + fns(i,k)*tx7)
3682 nstend(i,k) = nstend(i,k) + (dumns(i,k)-tx5) * oneodt
3685 qssedten(i,k) = qssedten(i,k) + tx6
3687 falouts(k) = fs(i,k) * dums(i,k)
3688 faloutns(k) = fns(i,k) * dumns(i,k)
3690 sflx(i,k+1) =
sflx(i,k+1) + falouts(k) * tx3
3693 prect(i) = prect(i) + falouts(nlev) * (tx3*0.001_r8)
3694 preci(i) = preci(i) + falouts(nlev) * (tx3*0.001_r8)
3701 if (do_graupel .or. do_hail)
then
3706 nstep = 1 + nint(max( maxval( fg(i,nlb:nlev)*pdel_inv(i,nlb:nlev)), &
3707 maxval(fng(i,nlb:nlev)*pdel_inv(i,nlb:nlev))) * deltat)
3708 nstep = min(nstep, nstep_def)
3723 tx7 = pdel_inv(i,k) * tx1
3724 dumg(i,k) = tx5 / (one + fg(i,k)*tx7)
3725 tx6 = (dumg(i,k)-tx5) * oneodt
3726 qgtend(i,k) = qgtend(i,k) + tx6
3728 dumng(i,k) = tx5 / (one + fng(i,k)*tx7)
3729 ngtend(i,k) = ngtend(i,k) + (dumng(i,k)-tx5) * oneodt
3732 qgsedten(i,k) = qgsedten(i,k) + tx6
3734 faloutg(k) = fg(i,k) * dumg(i,k)
3735 faloutng(k) = fng(i,k) * dumng(i,k)
3737 gflx(i,k+1) = gflx(i,k+1) + faloutg(k) * tx3
3742 tx7 = pdel_inv(i,k) * tx1
3743 dumg(i,k) = (tx5 + faloutg(k-1)*tx7) / (one + fg(i,k)*tx7)
3744 tx6 = (dumg(i,k)-tx5) * oneodt
3746 qgtend(i,k) = qgtend(i,k) + tx6
3748 dumng(i,k) = (tx5 + faloutng(k-1)*tx7) / (one + fng(i,k)*tx7)
3749 ngtend(i,k) = ngtend(i,k) + (dumng(i,k)-tx5) * oneodt
3752 qgsedten(i,k) = qgsedten(i,k) + tx6
3755 faloutg(k) = fg(i,k) * dumg(i,k)
3756 faloutng(k) = fng(i,k) * dumng(i,k)
3758 gflx(i,k+1) = gflx(i,k+1) + faloutg(k) * tx3
3765 prect(i) = prect(i) + faloutg(nlev) * (tx3*0.001_r8)
3766 preci(i) = preci(i) + faloutg(nlev) * (tx3*0.001_r8)
3784 dumc(i,k) = max(qc(i,k)+qctend(i,k)*deltat, zero)
3785 dumi(i,k) = max(qi(i,k)+qitend(i,k)*deltat, zero)
3786 dumnc(i,k) = max(nc(i,k)+nctend(i,k)*deltat, zero)
3787 dumni(i,k) = max(ni(i,k)+nitend(i,k)*deltat, zero)
3789 dumr(i,k) = max(qr(i,k)+qrtend(i,k)*deltat, zero)
3790 dumnr(i,k) = max(nr(i,k)+nrtend(i,k)*deltat, zero)
3791 dums(i,k) = max(qs(i,k)+qstend(i,k)*deltat, zero)
3792 dumns(i,k) = max(ns(i,k)+nstend(i,k)*deltat, zero)
3795 dumg(i,k) = max(qg(i,k)+qgtend(i,k)*deltat, zero)
3797 if (dumg(i,k) > 0.01_r8)
then
3798 tx2 = dumg(i,k) - 0.01_r8
3800 dums(i,k) = dums(i,k) + tx2
3801 qstend(i,k) = (dums(i,k) - qs(i,k)) * oneodt
3802 qgtend(i,k) = (dumg(i,k) - qg(i,k)) * oneodt
3805 dumng(i,k) = max(ng(i,k)+ngtend(i,k)*deltat, zero)
3810 dumnc(i,k) = ncnst*rhoinv(i,k)*lcldm(i,k)
3815 dumni(i,k) = ninst*rhoinv(i,k)*icldm(i,k)
3821 dumng(i,k) = ngnst*rhoinv(i,k)*precip_frac(i,k)
3825 if (dumc(i,k) < qsmall) dumnc(i,k) = zero
3826 if (dumi(i,k) < qsmall) dumni(i,k) = zero
3827 if (dumr(i,k) < qsmall) dumnr(i,k) = zero
3828 if (dums(i,k) < qsmall) dumns(i,k) = zero
3830 if (dumg(i,k) < qsmall) dumng(i,k) = zero
3845 tx1 = t(i,k) + tlat(i,k)*(deltat/cpp) - snowmelt
3846 if (tx1 > zero)
then
3847 if (dums(i,k) > zero)
then
3850 dum = -(xlf/cpp) * dums(i,k)
3851 if (tx1+dum < zero)
then
3852 dum = min(one, max(zero, -tx1/dum))
3858 qstend(i,k) = qstend(i,k) - tx2*dums(i,k)
3859 nstend(i,k) = nstend(i,k) - tx2*dumns(i,k)
3860 qrtend(i,k) = qrtend(i,k) + tx2*dums(i,k)
3861 nrtend(i,k) = nrtend(i,k) + tx2*dumns(i,k)
3863 dum1 = - xlf * tx2 * dums(i,k)
3864 tlat(i,k) = dum1 + tlat(i,k)
3865 meltsdttot(i,k) = dum1 + meltsdttot(i,k)
3871 if (do_graupel .or. do_hail)
then
3880 tx1 = t(i,k) + tlat(i,k)*(deltat/cpp) - snowmelt
3881 if (tx1 > zero)
then
3882 if (dumg(i,k) > zero)
then
3885 dum = -(xlf/cpp) * dumg(i,k)
3886 if (tx1+dum < zero)
then
3887 dum = max(zero, min(one, -tx1/dum))
3894 qgtend(i,k) = qgtend(i,k) - tx2*dumg(i,k)
3895 ngtend(i,k) = ngtend(i,k) - tx2*dumng(i,k)
3896 qrtend(i,k) = qrtend(i,k) + tx2*dumg(i,k)
3897 nrtend(i,k) = nrtend(i,k) + tx2*dumng(i,k)
3899 dum1 = - xlf*tx2*dumg(i,k)
3900 tlat(i,k) = dum1 + tlat(i,k)
3901 meltsdttot(i,k) = dum1 + meltsdttot(i,k)
3915 tx1 = t(i,k) + tlat(i,k) * (deltat/cpp) - rainfrze
3916 if (tx1 < zero)
then
3918 if (dumr(i,k) > zero)
then
3921 dum = (xlf/cpp) * dumr(i,k)
3922 if (tx1+dum > zero)
then
3923 dum = min(one, max(zero, -tx1/dum))
3928 qrtend(i,k) = qrtend(i,k) - tx2 * dumr(i,k)
3929 nrtend(i,k) = nrtend(i,k) - tx2 * dumnr(i,k)
3936 if (lamr(i,k) < one/dcs)
then
3938 if (do_hail .or. do_graupel)
then
3939 qgtend(i,k) = qgtend(i,k) + tx2 * dumr(i,k)
3940 ngtend(i,k) = ngtend(i,k) + tx2 * dumnr(i,k)
3942 qstend(i,k) = qstend(i,k) + tx2 * dumr(i,k)
3943 nstend(i,k) = nstend(i,k) + tx2 * dumnr(i,k)
3947 qitend(i,k) = qitend(i,k) + tx2 * dumr(i,k)
3948 nitend(i,k) = nitend(i,k) + tx2 * dumnr(i,k)
3951 dum1 = xlf*dum*dumr(i,k)*oneodt
3952 frzrdttot(i,k) = dum1 + frzrdttot(i,k)
3953 tlat(i,k) = dum1 + tlat(i,k)
3963 tx1 = t(i,k) + tlat(i,k) * (deltat/cpp) - tmelt
3964 if (tx1 > zero)
then
3965 if (dumi(i,k) > zero)
then
3969 dum = -dumi(i,k)*xlf/cpp
3970 if (tx1+dum < zero)
then
3971 dum = min(one, max(zero, -tx1/dum))
3977 qctend(i,k) = qctend(i,k) + tx2*dumi(i,k)
3980 melttot(i,k) = tx2*dumi(i,k)
3985 nctend(i,k) = nctend(i,k) + three*tx2*dumi(i,k)/(four*pi*5.12e-16_r8*rhow)
3987 qitend(i,k) = ((one-dum)*dumi(i,k)-qi(i,k)) * oneodt
3988 nitend(i,k) = ((one-dum)*dumni(i,k)-ni(i,k)) * oneodt
3989 tlat(i,k) = tlat(i,k) - xlf*tx2*dumi(i,k)
4004 tx1 = t(i,k) + tlat(i,k)*(deltat/cpp) - 233.15_r8
4005 if (tx1 < zero)
then
4006 if (dumc(i,k) > zero)
then
4009 dum = (xlf/cpp) * dumc(i,k)
4010 if (tx1+dum > zero)
then
4011 dum = min(one, max(zero, -tx1/dum))
4016 tx2 = dum * oneodt * dumc(i,k)
4017 qitend(i,k) = tx2 + qitend(i,k)
4023 nitend(i,k) = nitend(i,k) + tx2*(three/(four*pi*1.563e-14_r8* 500._r8))
4024 qctend(i,k) = ((one-dum)*dumc(i,k)-qc(i,k)) * oneodt
4025 nctend(i,k) = ((one-dum)*dumnc(i,k)-nc(i,k)) * oneodt
4026 tlat(i,k) = tlat(i,k) + xlf*tx2
4038 qtmp = q(i,k) + qvlat(i,k) * deltat
4039 ttmp = t(i,k) + tlat(i,k) * (deltat/cpp)
4043 esn = min(
fpvsl(ttmp), p(i,k))
4044 qvn = epsqs*esn/(p(i,k)-omeps*esn) * qsfm(i,k)
4048 if (qtmp > qvn .and. qvn > zero .and. allow_sed_supersat)
then
4050 dum = (qtmp-qvn)/(one+xxlv_squared*qvn/(cpp*rv*ttmp*ttmp)) * oneodt
4052 cmeout(i,k) = cmeout(i,k) + dum
4054 if (ttmp > 268.15_r8)
then
4058 else if (ttmp < 238.15_r8)
then
4061 dum1 = (268.15_r8-ttmp)/30._r8
4064 tx1 = xxls*dum1 + xxlv*(one-dum1)
4065 dum = (qtmp-qvn)/(one+tx1*tx1*qvn/(cpp*rv*ttmp*ttmp)) * oneodt
4066 tx2 = dum*(one-dum1)
4067 qctend(i,k) = qctend(i,k) + tx2
4069 qitend(i,k) = qitend(i,k) + dum*dum1
4070 qirestot(i,k) = dum*dum1
4071 qvlat(i,k) = qvlat(i,k) - dum
4074 tlat(i,k) = tlat(i,k) + dum*tx1
4091 if (lcldm(i,k) > mincld)
then
4092 tx1 = one / lcldm(i,k)
4096 if (icldm(i,k) > mincld)
then
4097 tx2 = one / icldm(i,k)
4101 if (precip_frac(i,k) > mincld)
then
4102 tx3 = one / precip_frac(i,k)
4106 dumc(i,k) = max(qc(i,k)+qctend(i,k)*deltat, zero) * tx1
4107 dumi(i,k) = max(qi(i,k)+qitend(i,k)*deltat, zero) * tx2
4108 dumnc(i,k) = max(nc(i,k)+nctend(i,k)*deltat, zero) * tx1
4109 dumni(i,k) = max(ni(i,k)+nitend(i,k)*deltat, zero) * tx2
4111 dumr(i,k) = max(qr(i,k)+qrtend(i,k)*deltat, zero) * tx3
4112 dumnr(i,k) = max(nr(i,k)+nrtend(i,k)*deltat, zero) * tx3
4113 dums(i,k) = max(qs(i,k)+qstend(i,k)*deltat, zero) * tx3
4114 dumns(i,k) = max(ns(i,k)+nstend(i,k)*deltat, zero) * tx3
4117 dumg(i,k) = max(qg(i,k)+qgtend(i,k)*deltat, zero) * tx3
4118 dumng(i,k) = max(ng(i,k)+ngtend(i,k)*deltat, zero) * tx3
4123 dumnc(i,k) = ncnst * rhoinv(i,k)
4128 dumni(i,k) = ninst * rhoinv(i,k)
4134 dumng(i,k) = ngnst*rhoinv(i,k)*precip_frac(i,k)
4141 dumc(i,k) = min(dumc(i,k), 10.e-3_r8)
4142 dumi(i,k) = min(dumi(i,k), 10.e-3_r8)
4144 dumr(i,k) = min(dumr(i,k), 10.e-3_r8)
4145 dums(i,k) = min(dums(i,k), 10.e-3_r8)
4147 dumg(i,k) = min(dumg(i,k), 10.e-3_r8)
4157 if (dumi(i,k) >= qsmall)
then
4163 if (dumni(i,k) /= tx1)
then
4165 nitend(i,k) = (dumni(i,k)*icldm(i,k)-ni(i,k)) * oneodt
4168 tx1 = one / lami(i,k)
4170 effi(i,k) = (three*1.e6_r8) * tx1
4171 sadice(i,k) = two*pi*(tx1*tx1*tx1)*dumni0*rho(i,k)*1.e-2_r8
4179 deffi(i,k) = effi(i,k) * (rhoi+rhoi)/rhows
4198 if (dumc(i,k) >= qsmall)
then
4208 nctend(i,k) = (ncnst*rhoinv(i,k)*lcldm(i,k)-nc(i,k)) * oneodt
4215 pgam(i,k), lamc(i,k))
4217 if (dum /= dumnc(i,k))
then
4219 nctend(i,k) = (dumnc(i,k)*lcldm(i,k)-nc(i,k)) * oneodt
4222 effc(i,k) = (half*1.e6_r8) * (pgam(i,k)+three) / lamc(i,k)
4224 lamcrad(i,k) = lamc(i,k)
4225 pgamrad(i,k) = pgam(i,k)
4233 dumnc(i,k) = 1.e8_r8
4238 pgam(i,k), lamc(i,k))
4240 effc_fn(i,k) = (half*1.e6_r8) * (pgam(i,k)+three)/lamc(i,k)
4255 if (dumr(i,k) >= qsmall)
then
4261 if (dum /= dumnr(i,k))
then
4263 nrtend(i,k) = (dumnr(i,k)*precip_frac(i,k)-nr(i,k)) *oneodt
4273 if (dums(i,k) >= qsmall)
then
4278 lams(i,k), n0=dumns0)
4280 if (dum /= dumns(i,k))
then
4282 nstend(i,k) = (dumns(i,k)*precip_frac(i,k)-ns(i,k)) * oneodt
4285 tx1 = (two*pi*1.e-2_r8) / (lams(i,k)*lams(i,k)*lams(i,k))
4286 sadsnow(i,k) = tx1*dumns0*rho(i,k)
4297 if (qc(i,k)+qctend(i,k)*deltat < qsmall) nctend(i,k) = -nc(i,k) * oneodt
4298 if (do_cldice .and. qi(i,k)+qitend(i,k)*deltat < qsmall) nitend(i,k) = -ni(i,k) * oneodt
4299 if (qr(i,k)+qrtend(i,k)*deltat < qsmall) nrtend(i,k) = -nr(i,k) * oneodt
4300 if (qs(i,k)+qstend(i,k)*deltat < qsmall) nstend(i,k) = -ns(i,k) * oneodt
4302 if (qg(i,k)+qgtend(i,k)*deltat < qsmall) ngtend(i,k) = -ng(i,k) * oneodt
4317 qc(i,k) = qc(i,k) + qctend(i,k)*deltat
4318 qi(i,k) = qi(i,k) + qitend(i,k)*deltat
4330 if (qrout(i,k) > 1.e-7_r8 .and. nrout(i,k) > zero)
then
4331 qrout2(i,k) = qrout(i,k) * precip_frac(i,k)
4332 nrout2(i,k) = nrout(i,k) * precip_frac(i,k)
4335 drout2(i,k) =
avg_diameter(qrout(i,k), nrout(i,k), rho(i,k), rhow)
4336 freqr(i,k) = precip_frac(i,k)
4338 reff_rain(i,k) = (1.e6_r8*1.5_r8) * drout2(i,k)
4344 reff_rain(i,k) = zero
4347 if (qsout(i,k) > 1.e-7_r8 .and. nsout(i,k) > zero)
then
4348 qsout2(i,k) = qsout(i,k) * precip_frac(i,k)
4349 nsout2(i,k) = nsout(i,k) * precip_frac(i,k)
4352 dsout2(i,k) =
avg_diameter(qsout(i,k), nsout(i,k), rho(i,k), rhosn)
4353 freqs(i,k) = precip_frac(i,k)
4355 dsout(i,k) = three*rhosn/rhows*dsout2(i,k)
4357 reff_snow(i,k) = (1.e6_r8*three) * dsout2(i,k)
4364 reff_snow(i,k) = zero
4379 if (qc(i,k) >= qsmall .and. (nc(i,k)+nctend(i,k)*deltat) > ten)
then
4380 tx1 = rho(i,k) / lcldm(i,k)
4381 tx2 = 1000._r8 * qc(i,k) * tx1
4382 dum = tx2 * tx2 * lcldm(i,k) &
4383 /(0.109_r8*(nc(i,k)+nctend(i,k)*deltat)*tx1*1.e-6_r8*precip_frac(i,k))
4390 if (qi(i,k) >= qsmall)
then
4392 dum1 = (qi(i,k)*rho(i,k)/icldm(i,k)*10000._r8)**(one/0.63_r8)*icldm(i,k)/precip_frac(i,k)
4397 if (qsout(i,k) >= qsmall)
then
4399 dum1 = dum1 + (qsout(i,k)*rho(i,k)*10000._r8)**(one/0.63_r8)
4402 refl(i,k) = dum + dum1
4409 if (rainrt(i,k) >= 0.001_r8)
then
4410 dum = rainrt(i,k) * rainrt(i,k)
4411 dum = log10(dum*dum*dum) + 16._r8
4415 dum = ten**(dum/ten)
4423 refl(i,k) = refl(i,k) + dum
4426 areflz(i,k) = refl(i,k) * precip_frac(i,k)
4430 if (refl(i,k) > minrefl)
then
4431 refl(i,k) = ten*log10(refl(i,k))
4433 refl(i,k) = -9999._r8
4437 if (refl(i,k) > mindbz)
then
4438 arefl(i,k) = refl(i,k) * precip_frac(i,k)
4439 frefl(i,k) = precip_frac(i,k)
4448 csrfl(i,k) = min(csmax,refl(i,k))
4451 if (csrfl(i,k) > csmin)
then
4452 acsrfl(i,k) = refl(i,k) * precip_frac(i,k)
4453 fcsrfl(i,k) = precip_frac(i,k)
4465 tx2 = qsout(i,k) + qi(i,k)
4466 tx1 = tx2 + qrout(i,k) + qc(i,k)
4467 if ( tx2 > qsmall .and. tx1 > qsmall)
then
4468 nfice(i,k) = min(tx2/tx1, one)