52 real,
parameter::
r3 = 1./3.
75 SUBROUTINE update_dz_c_fwd(is, ie, js, je, km, ng, dt, dp0, zs, area, &
76 & ut, vt, gz, ws, npx, npy, sw_corner, se_corner, ne_corner, nw_corner&
81 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km, npx, npy,
grid_type 82 LOGICAL,
INTENT(IN) :: sw_corner, se_corner, ne_corner, nw_corner
83 REAL,
INTENT(IN) :: dt
84 REAL,
INTENT(IN) :: dp0(km)
85 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: ut, vt
86 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng),
INTENT(IN) :: area
87 REAL,
INTENT(INOUT) :: gz(is-ng:ie+ng, js-ng:je+ng, km+1)
88 REAL,
INTENT(IN) :: zs(is-ng:ie+ng, js-ng:je+ng)
89 REAL :: ws(is-ng:ie+ng, js-ng:je+ng)
91 REAL :: gz2(is-ng:ie+ng, js-ng:je+ng)
92 REAL,
DIMENSION(is-1:ie+2, js-1:je+1) :: xfx, fx
93 REAL,
DIMENSION(is-1:ie+1, js-1:je+2) :: yfx, fy
94 REAL,
PARAMETER :: r14=1./14.
96 INTEGER :: is1, ie1, js1, je1
98 REAL :: rdt, top_ratio, bot_ratio, int_ratio
119 top_ratio = dp0(1)/(dp0(1)+dp0(2))
120 bot_ratio = dp0(km)/(dp0(km-1)+dp0(km))
137 xfx(i, j) = ut(i, j, 1) + (ut(i, j, 1)-ut(i, j, 2))*&
144 yfx(i, j) = vt(i, j, 1) + (vt(i, j, 1)-vt(i, j, 2))*&
149 ELSE IF (k .EQ. km + 1)
THEN 154 xfx(i, j) = ut(i, j, km) + (ut(i, j, km)-ut(i, j, km-1))*&
163 yfx(i, j) = vt(i, j, km) + (vt(i, j, km)-vt(i, j, km-1))*&
172 int_ratio = 1./(dp0(k-1)+dp0(k))
176 xfx(i, j) = (dp0(k)*ut(i, j, k-1)+dp0(k-1)*ut(i, j, k))*&
183 yfx(i, j) = (dp0(k)*vt(i, j, k-1)+dp0(k-1)*vt(i, j, k))*&
192 gz2(i, j) = gz(i, j, k)
197 & se_corner, ne_corner, nw_corner)
204 IF (xfx(i, j) .GT. 0.)
THEN 206 fx(i, j) = gz2(i-1, j)
214 fx(i, j) = xfx(i, j)*fx(i, j)
219 & se_corner, ne_corner, nw_corner)
226 IF (yfx(i, j) .GT. 0.)
THEN 228 fy(i, j) = gz2(i, j-1)
236 fy(i, j) = yfx(i, j)*fy(i, j)
242 gz(i, j, k) = (gz2(i, j)*area(i, j)+(fx(i, j)-fx(i+1, j))+(fy(&
243 & i, j)-fy(i, j+1)))/(area(i, j)+(xfx(i, j)-xfx(i+1, j))+(yfx(&
244 & i, j)-yfx(i, j+1)))
253 ws(i, j) = (zs(i, j)-gz(i, j, km+1))*rdt
257 IF (gz(i, j, k) .LT. gz(i, j, k+1) +
dz_min)
THEN 259 gz(i, j, k) = gz(i, j, k+1) +
dz_min 263 gz(i, j, k) = gz(i, j, k)
305 SUBROUTINE update_dz_c_bwd(is, ie, js, je, km, ng, dt, dp0, zs, area, &
306 & ut, ut_ad, vt, vt_ad, gz, gz_ad, ws, ws_ad, npx, npy, sw_corner, &
307 & se_corner, ne_corner, nw_corner, bd, grid_type)
310 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km, npx, npy,
grid_type 311 LOGICAL,
INTENT(IN) :: sw_corner, se_corner, ne_corner, nw_corner
312 REAL,
INTENT(IN) :: dt
313 REAL,
INTENT(IN) :: dp0(km)
314 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: ut, vt
315 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km) :: ut_ad, vt_ad
316 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng),
INTENT(IN) :: area
317 REAL,
INTENT(INOUT) :: gz(is-ng:ie+ng, js-ng:je+ng, km+1)
318 REAL,
INTENT(INOUT) :: gz_ad(is-ng:ie+ng, js-ng:je+ng, km+1)
319 REAL,
INTENT(IN) :: zs(is-ng:ie+ng, js-ng:je+ng)
320 REAL :: ws(is-ng:ie+ng, js-ng:je+ng)
321 REAL :: ws_ad(is-ng:ie+ng, js-ng:je+ng)
322 REAL :: gz2(is-ng:ie+ng, js-ng:je+ng)
323 REAL :: gz2_ad(is-ng:ie+ng, js-ng:je+ng)
324 REAL,
DIMENSION(is-1:ie+2, js-1:je+1) :: xfx, fx
325 REAL,
DIMENSION(is-1:ie+2, js-1:je+1) :: xfx_ad, fx_ad
326 REAL,
DIMENSION(is-1:ie+1, js-1:je+2) :: yfx, fy
327 REAL,
DIMENSION(is-1:ie+1, js-1:je+2) :: yfx_ad, fy_ad
328 REAL,
PARAMETER :: r14=1./14.
330 INTEGER :: is1, ie1, js1, je1
332 REAL :: rdt, top_ratio, bot_ratio, int_ratio
375 IF (branch .EQ. 0)
THEN 377 gz_ad(i, j, k+1) = gz_ad(i, j, k+1) + gz_ad(i, j, k)
386 gz_ad(i, j, km+1) = gz_ad(i, j, km+1) - rdt*ws_ad(i, j)
399 temp = area(i, j) + xfx(i, j) - xfx(i+1, j) + yfx(i, j) - yfx(&
401 temp_ad = gz_ad(i, j, k)/temp
402 temp_ad0 = -((area(i, j)*gz2(i, j)+fx(i, j)-fx(i+1, j)+fy(i, j&
403 & )-fy(i, j+1))*temp_ad/temp)
404 gz2_ad(i, j) = gz2_ad(i, j) + area(i, j)*temp_ad
405 fx_ad(i, j) = fx_ad(i, j) + temp_ad
406 fx_ad(i+1, j) = fx_ad(i+1, j) - temp_ad
407 fy_ad(i, j) = fy_ad(i, j) + temp_ad
408 fy_ad(i, j+1) = fy_ad(i, j+1) - temp_ad
409 xfx_ad(i, j) = xfx_ad(i, j) + temp_ad0
410 xfx_ad(i+1, j) = xfx_ad(i+1, j) - temp_ad0
411 yfx_ad(i, j) = yfx_ad(i, j) + temp_ad0
412 yfx_ad(i, j+1) = yfx_ad(i, j+1) - temp_ad0
419 yfx_ad(i, j) = yfx_ad(i, j) + fy(i, j)*fy_ad(i, j)
420 fy_ad(i, j) = yfx(i, j)*fy_ad(i, j)
422 IF (branch .EQ. 0)
THEN 424 gz2_ad(i, j-1) = gz2_ad(i, j-1) + fy_ad(i, j)
428 gz2_ad(i, j) = gz2_ad(i, j) + fy_ad(i, j)
435 & , npy, sw_corner, se_corner, &
436 & ne_corner, nw_corner)
440 xfx_ad(i, j) = xfx_ad(i, j) + fx(i, j)*fx_ad(i, j)
441 fx_ad(i, j) = xfx(i, j)*fx_ad(i, j)
443 IF (branch .EQ. 0)
THEN 445 gz2_ad(i-1, j) = gz2_ad(i-1, j) + fx_ad(i, j)
449 gz2_ad(i, j) = gz2_ad(i, j) + fx_ad(i, j)
456 & , npy, sw_corner, se_corner, &
457 & ne_corner, nw_corner)
461 gz_ad(i, j, k) = gz_ad(i, j, k) + gz2_ad(i, j)
466 IF (branch .EQ. 0)
THEN 470 vt_ad(i, j, k-1) = vt_ad(i, j, k-1) + int_ratio*dp0(k)*&
472 vt_ad(i, j, k) = vt_ad(i, j, k) + int_ratio*dp0(k-1)*yfx_ad(&
480 ut_ad(i, j, k-1) = ut_ad(i, j, k-1) + int_ratio*dp0(k)*&
482 ut_ad(i, j, k) = ut_ad(i, j, k) + int_ratio*dp0(k-1)*xfx_ad(&
488 ELSE IF (branch .EQ. 1)
THEN 492 vt_ad(i, j, km) = vt_ad(i, j, km) + (bot_ratio+1.0)*yfx_ad(i&
494 vt_ad(i, j, km-1) = vt_ad(i, j, km-1) - bot_ratio*yfx_ad(i, &
502 ut_ad(i, j, km) = ut_ad(i, j, km) + (bot_ratio+1.0)*xfx_ad(i&
504 ut_ad(i, j, km-1) = ut_ad(i, j, km-1) - bot_ratio*xfx_ad(i, &
513 vt_ad(i, j, 1) = vt_ad(i, j, 1) + (top_ratio+1.0)*yfx_ad(i, &
515 vt_ad(i, j, 2) = vt_ad(i, j, 2) - top_ratio*yfx_ad(i, j)
522 ut_ad(i, j, 1) = ut_ad(i, j, 1) + (top_ratio+1.0)*xfx_ad(i, &
524 ut_ad(i, j, 2) = ut_ad(i, j, 2) - top_ratio*xfx_ad(i, j)
531 SUBROUTINE update_dz_c(is, ie, js, je, km, ng, dt, dp0, zs, area, ut, &
532 & vt, gz, ws, npx, npy, sw_corner, se_corner, ne_corner, nw_corner, bd&
537 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km, npx, npy,
grid_type 538 LOGICAL,
INTENT(IN) :: sw_corner, se_corner, ne_corner, nw_corner
539 REAL,
INTENT(IN) :: dt
540 REAL,
INTENT(IN) :: dp0(km)
541 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: ut, vt
542 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng),
INTENT(IN) :: area
543 REAL,
INTENT(INOUT) :: gz(is-ng:ie+ng, js-ng:je+ng, km+1)
544 REAL,
INTENT(IN) :: zs(is-ng:ie+ng, js-ng:je+ng)
545 REAL,
INTENT(OUT) :: ws(is-ng:ie+ng, js-ng:je+ng)
547 REAL :: gz2(is-ng:ie+ng, js-ng:je+ng)
548 REAL,
DIMENSION(is-1:ie+2, js-1:je+1) :: xfx, fx
549 REAL,
DIMENSION(is-1:ie+1, js-1:je+2) :: yfx, fy
550 REAL,
PARAMETER :: r14=1./14.
552 INTEGER :: is1, ie1, js1, je1
554 REAL :: rdt, top_ratio, bot_ratio, int_ratio
558 top_ratio = dp0(1)/(dp0(1)+dp0(2))
559 bot_ratio = dp0(km)/(dp0(km-1)+dp0(km))
575 xfx(i, j) = ut(i, j, 1) + (ut(i, j, 1)-ut(i, j, 2))*&
581 yfx(i, j) = vt(i, j, 1) + (vt(i, j, 1)-vt(i, j, 2))*&
585 ELSE IF (k .EQ. km + 1)
THEN 589 xfx(i, j) = ut(i, j, km) + (ut(i, j, km)-ut(i, j, km-1))*&
597 yfx(i, j) = vt(i, j, km) + (vt(i, j, km)-vt(i, j, km-1))*&
604 int_ratio = 1./(dp0(k-1)+dp0(k))
607 xfx(i, j) = (dp0(k)*ut(i, j, k-1)+dp0(k-1)*ut(i, j, k))*&
613 yfx(i, j) = (dp0(k)*vt(i, j, k-1)+dp0(k-1)*vt(i, j, k))*&
620 gz2(i, j) = gz(i, j, k)
623 IF (
grid_type .LT. 3)
CALL fill_4corners(gz2, 1, bd, npx, npy, &
624 & sw_corner, se_corner, ne_corner&
628 IF (xfx(i, j) .GT. 0.)
THEN 629 fx(i, j) = gz2(i-1, j)
633 fx(i, j) = xfx(i, j)*fx(i, j)
636 IF (
grid_type .LT. 3)
CALL fill_4corners(gz2, 2, bd, npx, npy, &
637 & sw_corner, se_corner, ne_corner&
641 IF (yfx(i, j) .GT. 0.)
THEN 642 fy(i, j) = gz2(i, j-1)
646 fy(i, j) = yfx(i, j)*fy(i, j)
651 gz(i, j, k) = (gz2(i, j)*area(i, j)+(fx(i, j)-fx(i+1, j))+(fy(&
652 & i, j)-fy(i, j+1)))/(area(i, j)+(xfx(i, j)-xfx(i+1, j))+(yfx(&
653 & i, j)-yfx(i, j+1)))
661 ws(i, j) = (zs(i, j)-gz(i, j, km+1))*rdt
665 IF (gz(i, j, k) .LT. gz(i, j, k+1) +
dz_min)
THEN 666 gz(i, j, k) = gz(i, j, k+1) +
dz_min 668 gz(i, j, k) = gz(i, j, k)
695 & npx, npy, area, rarea, dp0, zs, zh, crx, cry, xfx, yfx, delz, ws, &
696 & rdt, gridstruct, bd, hord_pert)
699 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km, npx, npy
700 INTEGER,
INTENT(IN) :: hord, hord_pert
701 REAL,
INTENT(IN) :: rdt
702 REAL,
INTENT(IN) :: dp0(km)
703 REAL,
INTENT(IN) :: area(is-ng:ie+ng, js-ng:je+ng)
704 REAL,
INTENT(IN) :: rarea(is-ng:ie+ng, js-ng:je+ng)
705 REAL,
INTENT(INOUT) :: damp(km+1)
706 INTEGER,
INTENT(INOUT) :: ndif(km+1)
707 REAL,
INTENT(IN) :: zs(is-ng:ie+ng, js-ng:je+ng)
708 REAL,
INTENT(INOUT) :: zh(is-ng:ie+ng, js-ng:je+ng, km+1)
709 REAL,
INTENT(OUT) :: delz(is-ng:ie+ng, js-ng:je+ng, km)
710 REAL,
DIMENSION(is:ie+1, js-ng:je+ng, km),
INTENT(INOUT) :: crx, xfx
711 REAL,
DIMENSION(is-ng:ie+ng, js:je+1, km),
INTENT(INOUT) :: cry, yfx
712 REAL :: ws(is:ie, js:je)
716 REAL,
DIMENSION(is:ie+1, js-ng:je+ng, km+1) :: crx_adv, xfx_adv
717 REAL,
DIMENSION(is-ng:ie+ng, js:je+1, km+1) :: cry_adv, yfx_adv
718 REAL,
DIMENSION(is:ie+1, js:je) :: fx
719 REAL,
DIMENSION(is:ie, js:je+1) :: fy
720 REAL,
DIMENSION(is-ng:ie+ng+1, js-ng:je+ng) :: fx2
721 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng+1) :: fy2
722 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng) :: wk2, z2
723 REAL :: ra_x(is:ie, js-ng:je+ng)
724 REAL :: ra_y(is-ng:ie+ng, js:je)
726 INTEGER :: i, j, k, isd, ied, jsd, jed
727 LOGICAL :: uniform_grid
749 uniform_grid = .false.
751 damp(km+1) = damp(km)
753 ndif(km+1) = ndif(km)
763 & jed, j, km, dp0, uniform_grid, 0)
764 IF (j .LE. je + 1 .AND. j .GE. js)
THEN 767 & arg1, j, km, dp0, uniform_grid, 0)
781 ra_x(i, j) = area(i, j) + (xfx_adv(i, j, k)-xfx_adv(i+1, j, k)&
788 ra_y(i, j) = area(i, j) + (yfx_adv(i, j, k)-yfx_adv(i, j+1, k)&
792 IF (damp(k) .GT. 1.e-5)
THEN 796 z2(i, j) = zh(i, j, k)
799 IF (hord .EQ. hord_pert)
THEN 800 CALL fv_tp_2d_fwd(z2, crx_adv(is:ie+1, jsd:jed, k), cry_adv&
801 & (isd:ied, js:je+1, k), npx, npy, hord, fx, fy, &
802 & xfx_adv(is:ie+1, jsd:jed, k), yfx_adv(isd:ied, &
803 & js:je+1, k), gridstruct, bd, ra_x, ra_y)
809 CALL fv_tp_2d(z2, crx_adv(is:ie+1, jsd:jed, k), cry_adv(isd:&
810 & ied, js:je+1, k), npx, npy, hord_pert, fx, fy, xfx_adv&
811 & (is:ie+1, jsd:jed, k), yfx_adv(isd:ied, js:je+1, k), &
812 & gridstruct, bd, ra_x, ra_y)
815 CALL del6_vt_flux(ndif(k), npx, npy, damp(k), z2, wk2, fx2, fy2&
819 zh(i, j, k) = (z2(i, j)*area(i, j)+(fx(i, j)-fx(i+1, j))+(fy&
820 & (i, j)-fy(i, j+1)))/(ra_x(i, j)+ra_y(i, j)-area(i, j)) + (&
821 & fx2(i, j)-fx2(i+1, j)+(fy2(i, j)-fy2(i, j+1)))*rarea(i, j)
826 IF (hord .EQ. hord_pert)
THEN 827 CALL fv_tp_2d_fwd(zh(isd:ied, jsd:jed, k), crx_adv(is:ie+1&
828 & , jsd:jed, k), cry_adv(isd:ied, js:je+1, k), &
829 & npx, npy, hord, fx, fy, xfx_adv(is:ie+1, jsd:&
830 & jed, k), yfx_adv(isd:ied, js:je+1, k), &
831 & gridstruct, bd, ra_x, ra_y)
836 CALL pushrealarray(zh(isd:ied, jsd:jed, k), (ied-isd+1)*(jed-&
838 CALL fv_tp_2d(zh(isd:ied, jsd:jed, k), crx_adv(is:ie+1, jsd:&
839 & jed, k), cry_adv(isd:ied, js:je+1, k), npx, npy, &
840 & hord_pert, fx, fy, xfx_adv(is:ie+1, jsd:jed, k), &
841 & yfx_adv(isd:ied, js:je+1, k), gridstruct, bd, ra_x, &
848 zh(i, j, k) = (zh(i, j, k)*area(i, j)+(fx(i, j)-fx(i+1, j))+&
849 & (fy(i, j)-fy(i, j+1)))/(ra_x(i, j)+ra_y(i, j)-area(i, j))
861 ws(i, j) = (zs(i, j)-zh(i, j, km+1))*rdt
865 IF (zh(i, j, k) .LT. zh(i, j, k+1) +
dz_min)
THEN 867 zh(i, j, k) = zh(i, j, k+1) +
dz_min 871 zh(i, j, k) = zh(i, j, k)
912 & npx, npy, area, rarea, dp0, zs, zh, zh_ad, crx, crx_ad, cry, cry_ad&
913 & , xfx, xfx_ad, yfx, yfx_ad, delz, ws, ws_ad, rdt, gridstruct, bd, &
917 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km, npx, npy
918 INTEGER,
INTENT(IN) :: hord, hord_pert
919 REAL,
INTENT(IN) :: rdt
920 REAL,
INTENT(IN) :: dp0(km)
921 REAL,
INTENT(IN) :: area(is-ng:ie+ng, js-ng:je+ng)
922 REAL,
INTENT(IN) :: rarea(is-ng:ie+ng, js-ng:je+ng)
923 REAL,
INTENT(INOUT) :: damp(km+1)
924 INTEGER,
INTENT(INOUT) :: ndif(km+1)
925 REAL,
INTENT(IN) :: zs(is-ng:ie+ng, js-ng:je+ng)
926 REAL,
INTENT(INOUT) :: zh(is-ng:ie+ng, js-ng:je+ng, km+1)
927 REAL,
INTENT(INOUT) :: zh_ad(is-ng:ie+ng, js-ng:je+ng, km+1)
928 REAL,
INTENT(OUT) :: delz(is-ng:ie+ng, js-ng:je+ng, km)
929 REAL,
DIMENSION(is:ie+1, js-ng:je+ng, km),
INTENT(INOUT) :: crx, xfx
930 REAL,
DIMENSION(is:ie+1, js-ng:je+ng, km),
INTENT(INOUT) :: crx_ad, &
932 REAL,
DIMENSION(is-ng:ie+ng, js:je+1, km),
INTENT(INOUT) :: cry, yfx
933 REAL,
DIMENSION(is-ng:ie+ng, js:je+1, km),
INTENT(INOUT) :: cry_ad, &
935 REAL :: ws(is:ie, js:je)
936 REAL :: ws_ad(is:ie, js:je)
938 REAL,
DIMENSION(is:ie+1, js-ng:je+ng, km+1) :: crx_adv, xfx_adv
939 REAL,
DIMENSION(is:ie+1, js-ng:je+ng, km+1) :: crx_adv_ad, &
941 REAL,
DIMENSION(is-ng:ie+ng, js:je+1, km+1) :: cry_adv, yfx_adv
942 REAL,
DIMENSION(is-ng:ie+ng, js:je+1, km+1) :: cry_adv_ad, &
944 REAL,
DIMENSION(is:ie+1, js:je) :: fx
945 REAL,
DIMENSION(is:ie+1, js:je) :: fx_ad
946 REAL,
DIMENSION(is:ie, js:je+1) :: fy
947 REAL,
DIMENSION(is:ie, js:je+1) :: fy_ad
948 REAL,
DIMENSION(is-ng:ie+ng+1, js-ng:je+ng) :: fx2
949 REAL,
DIMENSION(is-ng:ie+ng+1, js-ng:je+ng) :: fx2_ad
950 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng+1) :: fy2
951 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng+1) :: fy2_ad
952 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng) :: wk2, z2
953 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng) :: wk2_ad, z2_ad
954 REAL :: ra_x(is:ie, js-ng:je+ng)
955 REAL :: ra_x_ad(is:ie, js-ng:je+ng)
956 REAL :: ra_y(is-ng:ie+ng, js:je)
957 REAL :: ra_y_ad(is-ng:ie+ng, js:je)
958 INTEGER :: i, j, k, isd, ied, jsd, jed
959 LOGICAL :: uniform_grid
990 CALL poprealarray(yfx_adv, (ie+2*ng-is+1)*(je-js+2)*(km+1))
991 CALL poprealarray(crx_adv, (ie-is+2)*(je+2*ng-js+1)*(km+1))
992 CALL poprealarray(xfx_adv, (ie-is+2)*(je+2*ng-js+1)*(km+1))
996 CALL poprealarray(cry_adv, (ie+2*ng-is+1)*(je-js+2)*(km+1))
1007 IF (branch .EQ. 0)
THEN 1009 zh_ad(i, j, k+1) = zh_ad(i, j, k+1) + zh_ad(i, j, k)
1010 zh_ad(i, j, k) = 0.0
1018 zh_ad(i, j, km+1) = zh_ad(i, j, km+1) - rdt*ws_ad(i, j)
1036 IF (branch .EQ. 0)
THEN 1040 temp0 = ra_x(i, j) - area(i, j) + ra_y(i, j)
1041 temp_ad2 = zh_ad(i, j, k)/temp0
1042 temp_ad3 = -((area(i, j)*zh(i, j, k)+fx(i, j)-fx(i+1, j)+fy(&
1043 & i, j)-fy(i, j+1))*temp_ad2/temp0)
1044 fx_ad(i, j) = fx_ad(i, j) + temp_ad2
1045 fx_ad(i+1, j) = fx_ad(i+1, j) - temp_ad2
1046 fy_ad(i, j) = fy_ad(i, j) + temp_ad2
1047 fy_ad(i, j+1) = fy_ad(i, j+1) - temp_ad2
1048 ra_x_ad(i, j) = ra_x_ad(i, j) + temp_ad3
1049 ra_y_ad(i, j) = ra_y_ad(i, j) + temp_ad3
1050 zh_ad(i, j, k) = area(i, j)*temp_ad2
1054 IF (branch .EQ. 0)
THEN 1055 CALL poprealarray(zh(isd:ied, jsd:jed, k), (ied-isd+1)*(jed-&
1059 CALL fv_tp_2d_adm(zh(isd:ied, jsd:jed, k), zh_ad(isd:ied, jsd:&
1060 & jed, k), crx_adv(is:ie+1, jsd:jed, k), crx_adv_ad(&
1061 & is:ie+1, jsd:jed, k), cry_adv(isd:ied, js:je+1, k)&
1062 & , cry_adv_ad(isd:ied, js:je+1, k), npx, npy, &
1063 & hord_pert, fx, fx_ad, fy, fy_ad, xfx_adv(is:ie+1, &
1064 & jsd:jed, k), xfx_adv_ad(is:ie+1, jsd:jed, k), &
1065 & yfx_adv(isd:ied, js:je+1, k), yfx_adv_ad(isd:ied, &
1066 & js:je+1, k), gridstruct, bd, ra_x, ra_x_ad, ra_y, &
1069 CALL fv_tp_2d_bwd(zh(isd:ied, jsd:jed, k), zh_ad(isd:ied, &
1070 & jsd:jed, k), crx_adv(is:ie+1, jsd:jed, k), &
1071 & crx_adv_ad(is:ie+1, jsd:jed, k), cry_adv(isd:&
1072 & ied, js:je+1, k), cry_adv_ad(isd:ied, js:je+1, &
1073 & k), npx, npy, hord, fx, fx_ad, fy, fy_ad, &
1074 & xfx_adv(is:ie+1, jsd:jed, k), xfx_adv_ad(is:ie+&
1075 & 1, jsd:jed, k), yfx_adv(isd:ied, js:je+1, k), &
1076 & yfx_adv_ad(isd:ied, js:je+1, k), gridstruct, bd&
1077 & , ra_x, ra_x_ad, ra_y, ra_y_ad)
1082 temp = ra_x(i, j) - area(i, j) + ra_y(i, j)
1083 temp_ad = zh_ad(i, j, k)/temp
1084 temp_ad0 = -((area(i, j)*z2(i, j)+fx(i, j)-fx(i+1, j)+fy(i, &
1085 & j)-fy(i, j+1))*temp_ad/temp)
1086 temp_ad1 = rarea(i, j)*zh_ad(i, j, k)
1087 z2_ad(i, j) = z2_ad(i, j) + area(i, j)*temp_ad
1088 fx_ad(i, j) = fx_ad(i, j) + temp_ad
1089 fx_ad(i+1, j) = fx_ad(i+1, j) - temp_ad
1090 fy_ad(i, j) = fy_ad(i, j) + temp_ad
1091 fy_ad(i, j+1) = fy_ad(i, j+1) - temp_ad
1092 ra_x_ad(i, j) = ra_x_ad(i, j) + temp_ad0
1093 ra_y_ad(i, j) = ra_y_ad(i, j) + temp_ad0
1094 fx2_ad(i, j) = fx2_ad(i, j) + temp_ad1
1095 fx2_ad(i+1, j) = fx2_ad(i+1, j) - temp_ad1
1096 fy2_ad(i, j) = fy2_ad(i, j) + temp_ad1
1097 fy2_ad(i, j+1) = fy2_ad(i, j+1) - temp_ad1
1098 zh_ad(i, j, k) = 0.0
1102 & , wk2_ad, fx2, fx2_ad, fy2, fy2_ad, gridstruct, &
1105 IF (branch .EQ. 0)
THEN 1106 CALL fv_tp_2d_bwd(z2, z2_ad, crx_adv(is:ie+1, jsd:jed, k), &
1107 & crx_adv_ad(is:ie+1, jsd:jed, k), cry_adv(isd:&
1108 & ied, js:je+1, k), cry_adv_ad(isd:ied, js:je+1, &
1109 & k), npx, npy, hord, fx, fx_ad, fy, fy_ad, &
1110 & xfx_adv(is:ie+1, jsd:jed, k), xfx_adv_ad(is:ie+&
1111 & 1, jsd:jed, k), yfx_adv(isd:ied, js:je+1, k), &
1112 & yfx_adv_ad(isd:ied, js:je+1, k), gridstruct, bd&
1113 & , ra_x, ra_x_ad, ra_y, ra_y_ad)
1118 CALL fv_tp_2d_adm(z2, z2_ad, crx_adv(is:ie+1, jsd:jed, k), &
1119 & crx_adv_ad(is:ie+1, jsd:jed, k), cry_adv(isd:ied, &
1120 & js:je+1, k), cry_adv_ad(isd:ied, js:je+1, k), npx&
1121 & , npy, hord_pert, fx, fx_ad, fy, fy_ad, xfx_adv(is&
1122 & :ie+1, jsd:jed, k), xfx_adv_ad(is:ie+1, jsd:jed, k&
1123 & ), yfx_adv(isd:ied, js:je+1, k), yfx_adv_ad(isd:&
1124 & ied, js:je+1, k), gridstruct, bd, ra_x, ra_x_ad, &
1130 zh_ad(i, j, k) = zh_ad(i, j, k) + z2_ad(i, j)
1138 yfx_adv_ad(i, j, k) = yfx_adv_ad(i, j, k) + ra_y_ad(i, j)
1139 yfx_adv_ad(i, j+1, k) = yfx_adv_ad(i, j+1, k) - ra_y_ad(i, j)
1146 xfx_adv_ad(i, j, k) = xfx_adv_ad(i, j, k) + ra_x_ad(i, j)
1147 xfx_adv_ad(i+1, j, k) = xfx_adv_ad(i+1, j, k) - ra_x_ad(i, j)
1155 & , cry_adv, cry_adv_ad, yfx_adv&
1156 & , yfx_adv_ad, isd, ied, js, &
1157 & arg1, j, km, dp0, uniform_grid&
1161 & crx_adv_ad, xfx_adv, xfx_adv_ad, is, arg1, jsd, &
1162 & jed, j, km, dp0, uniform_grid, 0)
1167 SUBROUTINE update_dz_d(ndif, damp, hord, is, ie, js, je, km, ng, npx, &
1168 & npy, area, rarea, dp0, zs, zh, crx, cry, xfx, yfx, delz, ws, rdt, &
1169 & gridstruct, bd, hord_pert)
1172 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km, npx, npy
1173 INTEGER,
INTENT(IN) :: hord, hord_pert
1174 REAL,
INTENT(IN) :: rdt
1175 REAL,
INTENT(IN) :: dp0(km)
1176 REAL,
INTENT(IN) :: area(is-ng:ie+ng, js-ng:je+ng)
1177 REAL,
INTENT(IN) :: rarea(is-ng:ie+ng, js-ng:je+ng)
1178 REAL,
INTENT(INOUT) :: damp(km+1)
1179 INTEGER,
INTENT(INOUT) :: ndif(km+1)
1180 REAL,
INTENT(IN) :: zs(is-ng:ie+ng, js-ng:je+ng)
1181 REAL,
INTENT(INOUT) :: zh(is-ng:ie+ng, js-ng:je+ng, km+1)
1182 REAL,
INTENT(OUT) :: delz(is-ng:ie+ng, js-ng:je+ng, km)
1183 REAL,
DIMENSION(is:ie+1, js-ng:je+ng, km),
INTENT(INOUT) :: crx, xfx
1184 REAL,
DIMENSION(is-ng:ie+ng, js:je+1, km),
INTENT(INOUT) :: cry, yfx
1185 REAL,
INTENT(OUT) :: ws(is:ie, js:je)
1189 REAL,
DIMENSION(is:ie+1, js-ng:je+ng, km+1) :: crx_adv, xfx_adv
1190 REAL,
DIMENSION(is-ng:ie+ng, js:je+1, km+1) :: cry_adv, yfx_adv
1191 REAL,
DIMENSION(is:ie+1, js:je) :: fx
1192 REAL,
DIMENSION(is:ie, js:je+1) :: fy
1193 REAL,
DIMENSION(is-ng:ie+ng+1, js-ng:je+ng) :: fx2
1194 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng+1) :: fy2
1195 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng) :: wk2, z2
1196 REAL :: ra_x(is:ie, js-ng:je+ng)
1197 REAL :: ra_y(is-ng:ie+ng, js:je)
1199 INTEGER :: i, j, k, isd, ied, jsd, jed
1200 LOGICAL :: uniform_grid
1203 uniform_grid = .false.
1204 damp(km+1) = damp(km)
1205 ndif(km+1) = ndif(km)
1214 CALL edge_profile(crx, xfx, crx_adv, xfx_adv, is, arg1, jsd, jed, &
1215 & j, km, dp0, uniform_grid, 0)
1216 IF (j .LE. je + 1 .AND. j .GE. js)
THEN 1218 CALL edge_profile(cry, yfx, cry_adv, yfx_adv, isd, ied, js, arg1&
1219 & , j, km, dp0, uniform_grid, 0)
1229 ra_x(i, j) = area(i, j) + (xfx_adv(i, j, k)-xfx_adv(i+1, j, k)&
1235 ra_y(i, j) = area(i, j) + (yfx_adv(i, j, k)-yfx_adv(i, j+1, k)&
1239 IF (damp(k) .GT. 1.e-5)
THEN 1242 z2(i, j) = zh(i, j, k)
1245 IF (hord .EQ. hord_pert)
THEN 1246 CALL fv_tp_2d(z2, crx_adv(is:ie+1, jsd:jed, k), cry_adv(isd&
1247 & :ied, js:je+1, k), npx, npy, hord, fx, fy, xfx_adv(&
1248 & is:ie+1, jsd:jed, k), yfx_adv(isd:ied, js:je+1, k)&
1249 & , gridstruct, bd, ra_x, ra_y)
1251 CALL fv_tp_2d(z2, crx_adv(is:ie+1, jsd:jed, k), cry_adv(isd:&
1252 & ied, js:je+1, k), npx, npy, hord_pert, fx, fy, xfx_adv&
1253 & (is:ie+1, jsd:jed, k), yfx_adv(isd:ied, js:je+1, k), &
1254 & gridstruct, bd, ra_x, ra_y)
1256 CALL del6_vt_flux(ndif(k), npx, npy, damp(k), z2, wk2, fx2, fy2&
1260 zh(i, j, k) = (z2(i, j)*area(i, j)+(fx(i, j)-fx(i+1, j))+(fy&
1261 & (i, j)-fy(i, j+1)))/(ra_x(i, j)+ra_y(i, j)-area(i, j)) + (&
1262 & fx2(i, j)-fx2(i+1, j)+(fy2(i, j)-fy2(i, j+1)))*rarea(i, j)
1266 IF (hord .EQ. hord_pert)
THEN 1267 CALL fv_tp_2d(zh(isd:ied, jsd:jed, k), crx_adv(is:ie+1, jsd&
1268 & :jed, k), cry_adv(isd:ied, js:je+1, k), npx, npy, &
1269 & hord, fx, fy, xfx_adv(is:ie+1, jsd:jed, k), yfx_adv&
1270 & (isd:ied, js:je+1, k), gridstruct, bd, ra_x, ra_y)
1272 CALL fv_tp_2d(zh(isd:ied, jsd:jed, k), crx_adv(is:ie+1, jsd:&
1273 & jed, k), cry_adv(isd:ied, js:je+1, k), npx, npy, &
1274 & hord_pert, fx, fy, xfx_adv(is:ie+1, jsd:jed, k), &
1275 & yfx_adv(isd:ied, js:je+1, k), gridstruct, bd, ra_x, &
1280 zh(i, j, k) = (zh(i, j, k)*area(i, j)+(fx(i, j)-fx(i+1, j))+&
1281 & (fy(i, j)-fy(i, j+1)))/(ra_x(i, j)+ra_y(i, j)-area(i, j))
1291 ws(i, j) = (zs(i, j)-zh(i, j, km+1))*rdt
1295 IF (zh(i, j, k) .LT. zh(i, j, k+1) +
dz_min)
THEN 1296 zh(i, j, k) = zh(i, j, k+1) +
dz_min 1298 zh(i, j, k) = zh(i, j, k)
1325 & cappa, cp, ptop, hs, w3, pt, q_con, delp, gz, pef, ws, p_fac, a_imp&
1328 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km
1329 INTEGER,
INTENT(IN) :: ms
1330 REAL,
INTENT(IN) :: dt, akap, cp, ptop, p_fac, a_imp, scale_m
1331 REAL,
INTENT(IN) :: ws(is-ng:ie+ng, js-ng:je+ng)
1332 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: pt, &
1334 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: q_con, &
1336 REAL,
INTENT(IN) :: hs(is-ng:ie+ng, js-ng:je+ng)
1337 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: w3
1339 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km+1),
INTENT(INOUT) :: gz
1340 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km+1) :: pef
1342 REAL,
DIMENSION(is-1:ie+1, km) :: dm, dz2, w2, pm2, gm2, cp2
1343 REAL,
DIMENSION(is-1:ie+1, km+1) :: pem, pe2, peg
1374 dm(i, k) = delp(i, j, k)
1387 pem(i, k) = pem(i, k-1) + dm(i, k-1)
1393 dz2(i, k) = gz(i, j, k+1) - gz(i, j, k)
1395 pm2(i, k) = dm(i, k)/log(pem(i, k+1)/pem(i, k))
1397 dm(i, k) = dm(i, k)*rgrav
1399 w2(i, k) = w3(i, j, k)
1402 IF (a_imp .LT. -0.01)
THEN 1404 & , dm, pem, w2, dz2, pt(is1:ie1, j, 1:km), ws(&
1405 & is1:ie1, j), p_fac, scale_m)
1407 ELSE IF (a_imp .LE. 0.5)
THEN 1409 & , pm2, w2, dz2, pt(is1:ie1, j, 1:km), ws(is1:ie1, j), &
1414 & akap, pe2, dm, pm2, pem, w2, dz2, pt(is1:ie1, j, &
1415 & 1:km), ws(is1:ie1, j), p_fac)
1422 pef(i, j, k) = pe2(i, k) + pem(i, k)
1428 gz(i, j, km+1) = hs(i, j)
1433 gz(i, j, k) = gz(i, j, k+1) - dz2(i, k)*
grav 1469 & cappa, cp, ptop, hs, w3, w3_ad, pt, pt_ad, q_con, delp, delp_ad, gz&
1470 & , gz_ad, pef, pef_ad, ws, ws_ad, p_fac, a_imp, scale_m)
1472 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km
1473 INTEGER,
INTENT(IN) :: ms
1474 REAL,
INTENT(IN) :: dt, akap, cp, ptop, p_fac, a_imp, scale_m
1475 REAL,
INTENT(IN) :: ws(is-ng:ie+ng, js-ng:je+ng)
1476 REAL :: ws_ad(is-ng:ie+ng, js-ng:je+ng)
1477 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: pt, &
1479 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km) :: pt_ad, delp_ad
1480 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: q_con, &
1482 REAL,
INTENT(IN) :: hs(is-ng:ie+ng, js-ng:je+ng)
1483 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: w3
1484 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km) :: w3_ad
1485 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km+1),
INTENT(INOUT) :: gz
1486 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km+1),
INTENT(INOUT) :: &
1488 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km+1) :: pef
1489 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km+1) :: pef_ad
1490 REAL,
DIMENSION(is-1:ie+1, km) :: dm, dz2, w2, pm2, gm2, cp2
1491 REAL,
DIMENSION(is-1:ie+1, km) :: dm_ad, dz2_ad, w2_ad, pm2_ad
1492 REAL,
DIMENSION(is-1:ie+1, km+1) :: pem, pe2, peg
1493 REAL,
DIMENSION(is-1:ie+1, km+1) :: pem_ad, pe2_ad
1539 gz_ad(i, j, k+1) = gz_ad(i, j, k+1) + gz_ad(i, j, k)
1540 dz2_ad(i, k) = dz2_ad(i, k) -
grav*gz_ad(i, j, k)
1541 gz_ad(i, j, k) = 0.0
1546 gz_ad(i, j, km+1) = 0.0
1551 pe2_ad(i, k) = pe2_ad(i, k) + pef_ad(i, j, k)
1552 pem_ad(i, k) = pem_ad(i, k) + pef_ad(i, j, k)
1553 pef_ad(i, j, k) = 0.0
1557 IF (branch .EQ. 0)
THEN 1559 & akap, pe2, pe2_ad, dm, dm_ad, pm2, pm2_ad, pem, &
1560 & pem_ad, w2, w2_ad, dz2, dz2_ad, pt(is1:ie1, j, 1:&
1561 & km), pt_ad(is1:ie1, j, 1:km), ws(is1:ie1, j), &
1562 & ws_ad(is1:ie1, j), p_fac)
1563 ELSE IF (branch .EQ. 1)
THEN 1565 & pe2_ad, dm, dm_ad, pm2, pm2_ad, w2, w2_ad, dz2, dz2_ad&
1566 & , pt(is1:ie1, j, 1:km), pt_ad(is1:ie1, j, 1:km), ws(&
1567 & is1:ie1, j), ws_ad(is1:ie1, j), .true.)
1570 & , pe2_ad, dm, dm_ad, pem, pem_ad, w2, w2_ad, &
1571 & dz2, dz2_ad, pt(is1:ie1, j, 1:km), pt_ad(is1:&
1572 & ie1, j, 1:km), ws(is1:ie1, j), ws_ad(is1:ie1, j&
1573 & ), p_fac, scale_m)
1578 temp = pem(i, k+1)/temp1
1581 w3_ad(i, j, k) = w3_ad(i, j, k) + w2_ad(i, k)
1584 dm_ad(i, k) = pm2_ad(i, k)/temp0 + rgrav*dm_ad(i, k)
1586 temp_ad = -(dm(i, k)*pm2_ad(i, k)/(temp*temp0**2*temp1))
1587 pem_ad(i, k+1) = pem_ad(i, k+1) + temp_ad
1588 pem_ad(i, k) = pem_ad(i, k) - temp*temp_ad
1591 gz_ad(i, j, k+1) = gz_ad(i, j, k+1) + dz2_ad(i, k)
1592 gz_ad(i, j, k) = gz_ad(i, j, k) - dz2_ad(i, k)
1599 pem_ad(i, k-1) = pem_ad(i, k-1) + pem_ad(i, k)
1600 dm_ad(i, k-1) = dm_ad(i, k-1) + pem_ad(i, k)
1608 pef_ad(i, j, 1) = 0.0
1613 delp_ad(i, j, k) = delp_ad(i, j, k) + dm_ad(i, k)
1619 SUBROUTINE riem_solver_c(ms, dt, is, ie, js, je, km, ng, akap, cappa, &
1620 & cp, ptop, hs, w3, pt, q_con, delp, gz, pef, ws, p_fac, a_imp, &
1623 INTEGER,
INTENT(IN) :: is, ie, js, je, ng, km
1624 INTEGER,
INTENT(IN) :: ms
1625 REAL,
INTENT(IN) :: dt, akap, cp, ptop, p_fac, a_imp, scale_m
1626 REAL,
INTENT(IN) :: ws(is-ng:ie+ng, js-ng:je+ng)
1627 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: pt, &
1629 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: q_con, &
1631 REAL,
INTENT(IN) :: hs(is-ng:ie+ng, js-ng:je+ng)
1632 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km),
INTENT(IN) :: w3
1634 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km+1),
INTENT(INOUT) :: gz
1635 REAL,
DIMENSION(is-ng:ie+ng, js-ng:je+ng, km+1),
INTENT(OUT) :: pef
1637 REAL,
DIMENSION(is-1:ie+1, km) :: dm, dz2, w2, pm2, gm2, cp2
1638 REAL,
DIMENSION(is-1:ie+1, km+1) :: pem, pe2, peg
1653 dm(i, k) = delp(i, j, k)
1663 pem(i, k) = pem(i, k-1) + dm(i, k-1)
1668 dz2(i, k) = gz(i, j, k+1) - gz(i, j, k)
1669 pm2(i, k) = dm(i, k)/log(pem(i, k+1)/pem(i, k))
1670 dm(i, k) = dm(i, k)*rgrav
1671 w2(i, k) = w3(i, j, k)
1674 IF (a_imp .LT. -0.01)
THEN 1676 & , pem, w2, dz2, pt(is1:ie1, j, 1:km), ws(is1:ie1, j&
1677 & ), p_fac, scale_m)
1678 ELSE IF (a_imp .LE. 0.5)
THEN 1679 CALL rim_2d(ms, dt, is1, ie1, km,
rdgas, gama, gm2, pe2, dm, pm2&
1680 & , w2, dz2, pt(is1:ie1, j, 1:km), ws(is1:ie1, j), .true.)
1683 & pe2, dm, pm2, pem, w2, dz2, pt(is1:ie1, j, 1:km), ws(&
1684 & is1:ie1, j), p_fac)
1689 pef(i, j, k) = pe2(i, k) + pem(i, k)
1694 gz(i, j, km+1) = hs(i, j)
1698 gz(i, j, k) = gz(i, j, k+1) - dz2(i, k)*
grav 1705 SUBROUTINE riem_solver3test(ms, dt, is, ie, js, je, km, ng, isd, ied, &
1706 & jsd, jed, akap, cappa, cp, ptop, zs, q_con, w, delz, pt, delp, zh, &
1707 & pe, ppe, pk3, pk, peln, ws, scale_m, p_fac, a_imp, use_logp, &
1708 & last_call, fp_out)
1716 INTEGER,
INTENT(IN) :: ms, is, ie, js, je, km, ng
1717 INTEGER,
INTENT(IN) :: isd, ied, jsd, jed
1719 REAL,
INTENT(IN) :: dt
1720 REAL,
INTENT(IN) :: akap, cp, ptop, p_fac, a_imp, scale_m
1721 REAL,
INTENT(IN) :: zs(isd:ied, jsd:jed)
1722 LOGICAL,
INTENT(IN) :: last_call, use_logp, fp_out
1723 REAL,
INTENT(IN) :: ws(is:ie, js:je)
1724 REAL,
DIMENSION(isd:ied, jsd:jed, km),
INTENT(IN) :: q_con, cappa
1725 REAL,
DIMENSION(isd:ied, jsd:jed, km),
INTENT(IN) :: delp, pt
1726 REAL,
DIMENSION(isd:ied, jsd:jed, km+1),
INTENT(INOUT) :: zh
1727 REAL,
DIMENSION(isd:ied, jsd:jed, km),
INTENT(INOUT) :: w
1728 REAL,
INTENT(INOUT) :: pe(is-1:ie+1, km+1, js-1:je+1)
1730 REAL,
INTENT(OUT) :: peln(is:ie, km+1, js:je)
1731 REAL,
DIMENSION(isd:ied, jsd:jed, km+1),
INTENT(OUT) :: ppe
1732 REAL,
INTENT(OUT) :: delz(is-ng:ie+ng, js-ng:je+ng, km)
1733 REAL,
INTENT(OUT) :: pk(is:ie, js:je, km+1)
1734 REAL,
INTENT(OUT) :: pk3(isd:ied, jsd:jed, km+1)
1736 REAL,
DIMENSION(is:ie, km) :: dm, dz2, pm2, w2, gm2, cp2
1737 REAL,
DIMENSION(is:ie, km+1) :: pem, pe2, peln2, peg, pelng
1738 REAL :: gama, rgrav, ptk, peln1
1747 ptk = exp(akap*peln1)
1755 dm(i, k) = delp(i, j, k)
1765 pem(i, k) = pem(i, k-1) + dm(i, k-1)
1766 peln2(i, k) = log(pem(i, k))
1767 pk3(i, j, k) = exp(akap*peln2(i, k))
1772 pm2(i, k) = dm(i, k)/(peln2(i, k+1)-peln2(i, k))
1773 dm(i, k) = dm(i, k)*rgrav
1774 dz2(i, k) = zh(i, j, k+1) - zh(i, j, k)
1775 w2(i, k) = w(i, j, k)
1778 IF (a_imp .LT. -0.999)
THEN 1780 & pem, w2, dz2, pt(is:ie, j, 1:km), ws(is:ie, j), &
1782 ELSE IF (a_imp .LT. -0.5)
THEN 1783 IF (a_imp .GE. 0.)
THEN 1789 & , w2, dz2, pt(is:ie, j, 1:km), ws(is:ie, j), abs0, &
1791 ELSE IF (a_imp .LE. 0.5)
THEN 1792 CALL rim_2d(ms, dt, is, ie, km,
rdgas, gama, gm2, pe2, dm, pm2, &
1793 & w2, dz2, pt(is:ie, j, 1:km), ws(is:ie, j), .false.)
1794 ELSE IF (a_imp .GT. 0.999)
THEN 1796 & pe2, dm, pm2, pem, w2, dz2, pt(is:ie, j, 1:km), ws(is&
1800 & , dm, pm2, pem, w2, dz2, pt(is:ie, j, 1:km), ws(is:ie&
1801 & , j), a_imp, p_fac, scale_m)
1805 w(i, j, k) = w2(i, k)
1806 delz(i, j, k) = dz2(i, k)
1812 peln(i, k, j) = peln2(i, k)
1813 pk(i, j, k) = pk3(i, j, k)
1814 pe(i, k, j) = pem(i, k)
1821 ppe(i, j, k) = pe2(i, k) + pem(i, k)
1827 ppe(i, j, k) = pe2(i, k)
1834 pk3(i, j, k) = peln2(i, k)
1839 zh(i, j, km+1) = zs(i, j)
1843 zh(i, j, k) = zh(i, j, k+1) - dz2(i, k)
1848 SUBROUTINE imp_diff_w(j, is, ie, js, je, ng, km, cd, delz, ws, w, w3)
1850 INTEGER,
INTENT(IN) :: j, is, ie, js, je, km, ng
1851 REAL,
INTENT(IN) :: cd
1853 REAL,
INTENT(IN) :: delz(is-ng:ie+ng, km)
1855 REAL,
INTENT(IN) :: w(is:ie, km)
1856 REAL,
INTENT(IN) :: ws(is:ie)
1857 REAL,
INTENT(OUT) :: w3(is-ng:ie+ng, js-ng:je+ng, km)
1859 REAL,
DIMENSION(is:ie, km) :: c, gam, dz, wt
1865 dz(i, k) = 0.5*(delz(i, k-1)+delz(i, k))
1870 c(i, k) = -(cd/(dz(i, k+1)*delz(i, k)))
1876 bet(i) = 1. - c(i, 1)
1877 wt(i, 1) = w(i, 1)/bet(i)
1882 gam(i, k) = c(i, k-1)/bet(i)
1883 a = cd/(dz(i, k)*delz(i, k))
1884 bet(i) = 1. + a - c(i, k) + a*gam(i, k)
1885 wt(i, k) = (w(i, k)+a*wt(i, k-1))/bet(i)
1890 gam(i, km) = c(i, km-1)/bet(i)
1891 a = cd/(dz(i, km)*delz(i, km))
1892 wt(i, km) = (w(i, km)+2.*ws(i)*cd/delz(i, km)**2+a*wt(i, km-1))/(&
1893 & 1.+a+(cd+cd)/delz(i, km)**2+a*gam(i, km))
1897 wt(i, k) = wt(i, k) - gam(i, k+1)*wt(i, k+1)
1902 w3(i, j, k) = wt(i, k)
1926 SUBROUTINE rim_2d_fwd(ms, bdt, is, ie, km, rgas, gama, gm2, pe2, dm2, &
1927 & pm2, w2, dz2, pt2, ws, c_core)
1929 INTEGER,
INTENT(IN) :: ms, is, ie, km
1930 REAL,
INTENT(IN) :: bdt, gama, rgas
1931 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pm2, gm2
1932 LOGICAL,
INTENT(IN) :: c_core
1933 REAL,
INTENT(IN) :: pt2(is:ie, km)
1934 REAL,
INTENT(IN) :: ws(is:ie)
1936 REAL,
INTENT(INOUT) :: dz2(is:ie, km)
1937 REAL,
INTENT(INOUT) :: w2(is:ie, km)
1938 REAL :: pe2(is:ie, km+1)
1941 REAL,
DIMENSION(km+1) :: m_bot, m_top, r_bot, r_top, pe1, pbar, wbar
1942 REAL,
DIMENSION(km) :: r_hi, r_lo, dz, wm, dm, dts
1943 REAL,
DIMENSION(km) :: pf1, wc, cm, pp, pt1
1944 REAL :: dt, rdt, grg, z_frac, ptmp1, rden, pf, time_left
1946 INTEGER :: i, k, n, ke, kt1, ktop
1955 INTEGER :: ad_count0
1958 INTEGER :: ad_count1
1960 INTEGER :: ad_count2
1961 INTEGER :: ad_count3
2035 wm(k) = w2(i, k)*dm(k)
2044 IF (ms .GT. 1 .AND. ms .LT. 8)
THEN 2049 rden = -(rgas*dm(k)/dz(k))
2051 pf1(k) = exp(gama*log(rden*pt1(k)))
2053 dts(k) = -(dz(k)/sqrt(grg*pf1(k)/rden))
2054 IF (bdt .GT. dts(k))
THEN 2058 ad_count = ad_count + 1
2068 222
IF (ks0 .EQ. 1)
THEN 2074 cm(k) = dm(k)/dts(k)
2076 wc(k) = wm(k)/dts(k)
2078 pp(k) = pf1(k) - pm2(i, k)
2082 wbar(1) = (wc(1)+pp(1))/cm(1)
2085 wbar(k) = (wc(k-1)+wc(k)+pp(k)-pp(k-1))/(cm(k-1)+cm(k))
2087 pbar(k) = bdt*(cm(k-1)*wbar(k)-wc(k-1)+pp(k-1))
2091 IF (ks0 .EQ. km)
THEN 2093 pbar(km+1) = bdt*(cm(km)*wbar(km+1)-wc(km)+pp(km))
2097 dz2(i, k) = dz(k) + bdt*(wbar(k+1)-wbar(k))
2103 dz2(i, k) = dz(k) + bdt*(wbar(k+1)-wbar(k))
2105 w2(i, k) = (wm(k)+pbar(k+1)-pbar(k))/dm(k)
2113 pe2(i, k) = pbar(k)*rdt
2122 dz2(i, k) = dz(k) + bdt*(wbar(k+1)-wbar(k))
2129 dz2(i, k) = dz(k) + bdt*(wbar(k+1)-wbar(k))
2131 w2(i, k) = (wm(k)+pbar(k+1)-pbar(k))/dm(k)
2137 pbar(ks0) = pbar(ks0)/
REAL(ms)
2149 rden = -(rgas*dm(k)/dz(k))
2151 pf = exp(gama*log(rden*pt1(k)))
2153 dts(k) = -(dz(k)/sqrt(grg*pf/rden))
2154 ptmp1 = dts(k)*(pf-pm2(i, k))
2156 r_lo(k) = wm(k) + ptmp1
2158 r_hi(k) = wm(k) - ptmp1
2163 IF (dt .GT. dts(k))
THEN 2166 ad_count0 = ad_count0 + 1
2176 333
IF (ktop .GE. ks1)
THEN 2181 r_bot(k) = z_frac*r_lo(k)
2183 r_top(k+1) = z_frac*r_hi(k)
2185 m_bot(k) = z_frac*dm(k)
2187 m_top(k+1) = m_bot(k)
2191 IF (ktop .EQ. km)
THEN 2208 IF (1 .LT. ktop)
THEN 2213 DO 130 ke=km+1,ktop+2,-1
2219 IF (time_left .GT. dts(k))
THEN 2220 time_left = time_left - dts(k)
2222 m_top(ke) = m_top(ke) + dm(k)
2224 r_top(ke) = r_top(ke) + r_hi(k)
2226 ad_count1 = ad_count1 + 1
2237 z_frac = time_left/dts(k)
2239 m_top(ke) = m_top(ke) + z_frac*dm(k)
2241 r_top(ke) = r_top(ke) + z_frac*r_hi(k)
2256 DO 160 ke=ad_from3,km
2262 IF (time_left .GT. dts(k))
THEN 2263 time_left = time_left - dts(k)
2265 m_bot(ke) = m_bot(ke) + dm(k)
2267 r_bot(ke) = r_bot(ke) + r_lo(k)
2269 ad_count2 = ad_count2 + 1
2282 IF (time_left .GT. dts(k))
THEN 2283 time_left = time_left - dts(k)
2285 m_bot(ke) = m_bot(ke) + dm(k)
2287 r_bot(ke) = r_bot(ke) - r_hi(k)
2289 ad_count3 = ad_count3 + 1
2300 z_frac = time_left/dts(k)
2302 m_bot(ke) = m_bot(ke) + z_frac*dm(k)
2304 r_bot(ke) = r_bot(ke) - z_frac*r_hi(k) + (m_bot(ke)-m_surf)*&
2310 z_frac = time_left/dts(k)
2312 m_bot(ke) = m_bot(ke) + z_frac*dm(k)
2314 r_bot(ke) = r_bot(ke) + z_frac*r_lo(k)
2320 666
IF (ks1 .EQ. 1)
THEN 2322 wbar(1) = r_bot(1)/m_bot(1)
2331 wbar(k) = (r_bot(k)+r_top(k))/(m_top(k)+m_bot(k))
2338 pbar(k) = m_top(k)*wbar(k) - r_top(k)
2339 pe1(k) = pe1(k) + pbar(k)
2347 dz2(i, k) = dz(k) + dt*(wbar(k+1)-wbar(k))
2355 dz2(i, k) = dz(k) + dt*(wbar(k+1)-wbar(k))
2357 w2(i, k) = (wm(k)+pbar(k+1)-pbar(k))/dm(k)
2366 dz(k) = dz(k) + dt*(wbar(k+1)-wbar(k))
2368 wm(k) = wm(k) + pbar(k+1) - pbar(k)
2379 pe2(i, k) = pe1(k)*rdt
2429 SUBROUTINE rim_2d_bwd(ms, bdt, is, ie, km, rgas, gama, gm2, pe2, &
2430 & pe2_ad, dm2, dm2_ad, pm2, pm2_ad, w2, w2_ad, dz2, dz2_ad, pt2, &
2431 & pt2_ad, ws, ws_ad, c_core)
2433 INTEGER,
INTENT(IN) :: ms, is, ie, km
2434 REAL,
INTENT(IN) :: bdt, gama, rgas
2435 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pm2, gm2
2436 REAL,
DIMENSION(is:ie, km) :: dm2_ad, pm2_ad
2437 LOGICAL,
INTENT(IN) :: c_core
2438 REAL,
INTENT(IN) :: pt2(is:ie, km)
2439 REAL :: pt2_ad(is:ie, km)
2440 REAL,
INTENT(IN) :: ws(is:ie)
2441 REAL :: ws_ad(is:ie)
2442 REAL,
INTENT(INOUT) :: dz2(is:ie, km)
2443 REAL,
INTENT(INOUT) :: dz2_ad(is:ie, km)
2444 REAL,
INTENT(INOUT) :: w2(is:ie, km)
2445 REAL,
INTENT(INOUT) :: w2_ad(is:ie, km)
2446 REAL :: pe2(is:ie, km+1)
2447 REAL :: pe2_ad(is:ie, km+1)
2449 REAL :: ws2_ad(is:ie)
2450 REAL,
DIMENSION(km+1) :: m_bot, m_top, r_bot, r_top, pe1, pbar, wbar
2451 REAL,
DIMENSION(km+1) :: m_bot_ad, m_top_ad, r_bot_ad, r_top_ad, &
2452 & pe1_ad, pbar_ad, wbar_ad
2453 REAL,
DIMENSION(km) :: r_hi, r_lo, dz, wm, dm, dts
2454 REAL,
DIMENSION(km) :: r_hi_ad, r_lo_ad, dz_ad, wm_ad, dm_ad, dts_ad
2455 REAL,
DIMENSION(km) :: pf1, wc, cm, pp, pt1
2456 REAL,
DIMENSION(km) :: pf1_ad, wc_ad, cm_ad, pp_ad, pt1_ad
2457 REAL :: dt, rdt, grg, z_frac, ptmp1, rden, pf, time_left
2458 REAL :: z_frac_ad, ptmp1_ad, rden_ad, pf_ad, time_left_ad
2461 INTEGER :: i, k, n, ke, kt1, ktop
2502 INTEGER :: ad_count0
2507 INTEGER :: ad_count1
2511 INTEGER :: ad_count2
2513 INTEGER :: ad_count3
2624 IF (branch .EQ. 0)
THEN 2627 pbar_ad(k) = pbar_ad(k) + rdt*pe2_ad(i, k)
2633 IF (branch .EQ. 0)
THEN 2636 dz_ad(k) = dz_ad(k) + dz2_ad(i, k)
2637 wbar_ad(k+1) = wbar_ad(k+1) + bdt*dz2_ad(i, k)
2638 wbar_ad(k) = wbar_ad(k) - bdt*dz2_ad(i, k)
2644 temp_ad9 = w2_ad(i, k)/dm(k)
2645 wm_ad(k) = wm_ad(k) + temp_ad9
2646 pbar_ad(k+1) = pbar_ad(k+1) + temp_ad9
2647 pbar_ad(k) = pbar_ad(k) - temp_ad9
2648 dm_ad(k) = dm_ad(k) - (wm(k)+pbar(k+1)-pbar(k))*temp_ad9/dm(&
2652 dz_ad(k) = dz_ad(k) + dz2_ad(i, k)
2653 wbar_ad(k+1) = wbar_ad(k+1) + bdt*dz2_ad(i, k)
2654 wbar_ad(k) = wbar_ad(k) - bdt*dz2_ad(i, k)
2659 temp_ad8 = bdt*pbar_ad(km+1)
2660 cm_ad(km) = cm_ad(km) + wbar(km+1)*temp_ad8
2661 wbar_ad(km+1) = wbar_ad(km+1) + cm(km)*temp_ad8
2662 wc_ad(km) = wc_ad(km) - temp_ad8
2663 pp_ad(km) = pp_ad(km) + temp_ad8
2670 pe1_ad(k) = pe1_ad(k) + rdt*pe2_ad(i, k)
2678 IF (branch .EQ. 0)
THEN 2682 dz_ad(k) = dz_ad(k) + dz2_ad(i, k)
2683 wbar_ad(k+1) = wbar_ad(k+1) + dt*dz2_ad(i, k)
2684 wbar_ad(k) = wbar_ad(k) - dt*dz2_ad(i, k)
2687 ELSE IF (branch .EQ. 1)
THEN 2691 temp_ad20 = w2_ad(i, k)/dm(k)
2692 wm_ad(k) = wm_ad(k) + temp_ad20
2693 pbar_ad(k+1) = pbar_ad(k+1) + temp_ad20
2694 pbar_ad(k) = pbar_ad(k) - temp_ad20
2695 dm_ad(k) = dm_ad(k) - (wm(k)+pbar(k+1)-pbar(k))*temp_ad20/&
2699 dz_ad(k) = dz_ad(k) + dz2_ad(i, k)
2700 wbar_ad(k+1) = wbar_ad(k+1) + dt*dz2_ad(i, k)
2701 wbar_ad(k) = wbar_ad(k) - dt*dz2_ad(i, k)
2708 pbar_ad(k+1) = pbar_ad(k+1) + wm_ad(k)
2709 pbar_ad(k) = pbar_ad(k) - wm_ad(k)
2711 wbar_ad(k+1) = wbar_ad(k+1) + dt*dz_ad(k)
2712 wbar_ad(k) = wbar_ad(k) - dt*dz_ad(k)
2716 DO k=km+1,ad_from5,-1
2717 pbar_ad(k) = pbar_ad(k) + pe1_ad(k)
2719 m_top_ad(k) = m_top_ad(k) + wbar(k)*pbar_ad(k)
2720 wbar_ad(k) = wbar_ad(k) + m_top(k)*pbar_ad(k)
2721 r_top_ad(k) = r_top_ad(k) - pbar_ad(k)
2727 temp_ad18 = wbar_ad(k)/(m_top(k)+m_bot(k))
2728 temp_ad19 = -((r_bot(k)+r_top(k))*temp_ad18/(m_top(k)+m_bot(&
2730 r_bot_ad(k) = r_bot_ad(k) + temp_ad18
2731 r_top_ad(k) = r_top_ad(k) + temp_ad18
2732 m_top_ad(k) = m_top_ad(k) + temp_ad19
2733 m_bot_ad(k) = m_bot_ad(k) + temp_ad19
2738 IF (branch .NE. 0)
THEN 2740 temp_ad17 = wbar_ad(1)/m_bot(1)
2741 r_bot_ad(1) = r_bot_ad(1) + temp_ad17
2742 m_bot_ad(1) = m_bot_ad(1) - r_bot(1)*temp_ad17/m_bot(1)
2746 IF (branch .NE. 0)
THEN 2748 DO ke=km,ad_from3,-1
2750 IF (branch .EQ. 0)
THEN 2751 z_frac = time_left/dts(k)
2753 z_frac_ad = dm(k)*m_bot_ad(ke) + r_lo(k)*r_bot_ad(ke)
2754 r_lo_ad(k) = r_lo_ad(k) + z_frac*r_bot_ad(ke)
2756 dm_ad(k) = dm_ad(k) + z_frac*m_bot_ad(ke)
2757 temp_ad16 = z_frac_ad/dts(k)
2758 time_left_ad = temp_ad16
2759 dts_ad(k) = dts_ad(k) - time_left*temp_ad16/dts(k)
2761 IF (branch .EQ. 1)
THEN 2762 m_bot_ad(ke) = m_bot_ad(ke) + ws2(i)*r_bot_ad(ke)
2763 z_frac = time_left/dts(k)
2765 z_frac_ad = dm(k)*m_bot_ad(ke) - r_hi(k)*r_bot_ad(ke)
2766 r_hi_ad(k) = r_hi_ad(k) - z_frac*r_bot_ad(ke)
2767 m_surf_ad = -(ws2(i)*r_bot_ad(ke))
2768 ws2_ad(i) = ws2_ad(i) + (m_bot(ke)-m_surf)*r_bot_ad(ke&
2771 dm_ad(k) = dm_ad(k) + z_frac*m_bot_ad(ke)
2772 temp_ad15 = z_frac_ad/dts(k)
2773 time_left_ad = temp_ad15
2774 dts_ad(k) = dts_ad(k) - time_left*temp_ad15/dts(k)
2783 IF (branch .EQ. 0)
THEN 2789 r_hi_ad(k) = r_hi_ad(k) - r_bot_ad(ke)
2791 dm_ad(k) = dm_ad(k) + m_bot_ad(ke)
2792 dts_ad(k) = dts_ad(k) - time_left_ad
2797 m_bot_ad(ke) = m_bot_ad(ke) + m_surf_ad
2805 r_lo_ad(k) = r_lo_ad(k) + r_bot_ad(ke)
2807 dm_ad(k) = dm_ad(k) + m_bot_ad(ke)
2808 dts_ad(k) = dts_ad(k) - time_left_ad
2825 IF (branch .EQ. 0)
THEN 2826 z_frac = time_left/dts(k)
2828 z_frac_ad = dm(k)*m_top_ad(ke) + r_hi(k)*r_top_ad(ke)
2829 r_hi_ad(k) = r_hi_ad(k) + z_frac*r_top_ad(ke)
2831 dm_ad(k) = dm_ad(k) + z_frac*m_top_ad(ke)
2832 temp_ad14 = z_frac_ad/dts(k)
2833 time_left_ad = temp_ad14
2834 dts_ad(k) = dts_ad(k) - time_left*temp_ad14/dts(k)
2842 IF (branch .EQ. 0) time_left_ad = 0.0
2845 r_hi_ad(k) = r_hi_ad(k) + r_top_ad(ke)
2847 dm_ad(k) = dm_ad(k) + m_top_ad(ke)
2848 dts_ad(k) = dts_ad(k) - time_left_ad
2855 DO k=km+1,ad_from1,-1
2862 IF (branch .EQ. 0)
GOTO 100
2866 DO k=ad_to3,ad_from0,-1
2868 m_bot_ad(k) = m_bot_ad(k) + m_top_ad(k+1)
2872 z_frac_ad = r_hi(k)*r_top_ad(k+1) + r_lo(k)*r_bot_ad(k) + dm&
2874 dm_ad(k) = dm_ad(k) + z_frac*m_bot_ad(k)
2877 r_hi_ad(k) = r_hi_ad(k) + z_frac*r_top_ad(k+1)
2880 r_lo_ad(k) = r_lo_ad(k) + z_frac*r_bot_ad(k)
2882 dts_ad(k) = dts_ad(k) - dt*z_frac_ad/dts(k)**2
2891 wm_ad(k) = wm_ad(k) + r_lo_ad(k) + r_hi_ad(k)
2892 ptmp1_ad = r_lo_ad(k) - r_hi_ad(k)
2896 dts_ad(k) = dts_ad(k) + (pf-pm2(i, k))*ptmp1_ad
2897 pm2_ad(i, k) = pm2_ad(i, k) - dts(k)*ptmp1_ad
2898 rden = -(rgas*dm(k)/dz(k))
2900 temp2 = sqrt(grg*temp1)
2901 IF (grg*temp1 .EQ. 0.0)
THEN 2904 temp_ad11 = grg*dz(k)*dts_ad(k)/(2.0*temp2**3*rden)
2906 pf_ad = temp_ad11 + dts(k)*ptmp1_ad
2909 temp_ad13 = gama*exp(gama*log(rden*pt1(k)))*pf_ad/(rden*pt1(&
2911 rden_ad = pt1(k)*temp_ad13 - temp1*temp_ad11
2912 pt1_ad(k) = pt1_ad(k) + rden*temp_ad13
2913 temp_ad12 = -(rgas*rden_ad/dz(k))
2914 dz_ad(k) = dz_ad(k) - dm(k)*temp_ad12/dz(k) - dts_ad(k)/&
2917 dm_ad(k) = dm_ad(k) + temp_ad12
2922 IF (branch .EQ. 0)
THEN 2924 ELSE IF (branch .EQ. 1)
THEN 2926 pbar_ad(ks0) = pbar_ad(ks0)/
REAL(ms)
2928 IF (branch .EQ. 0)
THEN 2932 dz_ad(k) = dz_ad(k) + dz2_ad(i, k)
2933 wbar_ad(k+1) = wbar_ad(k+1) + bdt*dz2_ad(i, k)
2934 wbar_ad(k) = wbar_ad(k) - bdt*dz2_ad(i, k)
2941 temp_ad10 = w2_ad(i, k)/dm(k)
2942 wm_ad(k) = wm_ad(k) + temp_ad10
2943 pbar_ad(k+1) = pbar_ad(k+1) + temp_ad10
2944 pbar_ad(k) = pbar_ad(k) - temp_ad10
2945 dm_ad(k) = dm_ad(k) - (wm(k)+pbar(k+1)-pbar(k))*temp_ad10/&
2949 dz_ad(k) = dz_ad(k) + dz2_ad(i, k)
2950 wbar_ad(k+1) = wbar_ad(k+1) + bdt*dz2_ad(i, k)
2951 wbar_ad(k) = wbar_ad(k) - bdt*dz2_ad(i, k)
2961 pbar_ad(k) = pbar_ad(k) + pe1_ad(k)
2964 temp_ad5 = bdt*pbar_ad(k)
2965 wbar_ad(k) = wbar_ad(k) + cm(k-1)*temp_ad5
2966 pp_ad(k-1) = pp_ad(k-1) + temp_ad5
2968 temp_ad7 = wbar_ad(k)/(cm(k-1)+cm(k))
2969 wc_ad(k-1) = wc_ad(k-1) + temp_ad7 - temp_ad5
2970 temp_ad6 = -((wc(k-1)+wc(k)+pp(k)-pp(k-1))*temp_ad7/(cm(k-1)+cm(&
2972 cm_ad(k-1) = cm_ad(k-1) + temp_ad6 + wbar(k)*temp_ad5
2974 wc_ad(k) = wc_ad(k) + temp_ad7
2975 pp_ad(k) = pp_ad(k) + temp_ad7
2976 pp_ad(k-1) = pp_ad(k-1) - temp_ad7
2977 cm_ad(k) = cm_ad(k) + temp_ad6
2981 temp_ad4 = wbar_ad(1)/cm(1)
2982 wc_ad(1) = wc_ad(1) + temp_ad4
2983 pp_ad(1) = pp_ad(1) + temp_ad4
2984 cm_ad(1) = cm_ad(1) - (wc(1)+pp(1))*temp_ad4/cm(1)
2988 temp_ad3 = cm_ad(k)/dts(k)
2990 pf1_ad(k) = pf1_ad(k) + pp_ad(k)
2991 pm2_ad(i, k) = pm2_ad(i, k) - pp_ad(k)
2994 temp_ad2 = wc_ad(k)/dts(k)
2995 wm_ad(k) = wm_ad(k) + temp_ad2
2996 dts_ad(k) = dts_ad(k) - dm(k)*temp_ad3/dts(k) - wm(k)*temp_ad2/&
3000 dm_ad(k) = dm_ad(k) + temp_ad3
3008 IF (branch .EQ. 0)
GOTO 120
3010 rden = -(rgas*dm(k)/dz(k))
3013 temp0 = sqrt(grg*temp)
3014 IF (grg*temp .EQ. 0.0)
THEN 3017 temp_ad = grg*dz(k)*dts_ad(k)/(2.0*temp0**3*rden)
3019 pf1_ad(k) = pf1_ad(k) + temp_ad
3021 temp_ad1 = gama*exp(gama*log(rden*pt1(k)))*pf1_ad(k)/(rden*pt1(k&
3023 rden_ad = pt1(k)*temp_ad1 - temp*temp_ad
3024 pt1_ad(k) = pt1_ad(k) + rden*temp_ad1
3026 temp_ad0 = -(rgas*rden_ad/dz(k))
3027 dz_ad(k) = dz_ad(k) - dm(k)*temp_ad0/dz(k) - dts_ad(k)/temp0
3029 dm_ad(k) = dm_ad(k) + temp_ad0
3034 ws_ad(i) = ws_ad(i) + wbar_ad(km+1)
3038 pt2_ad(i, k) = pt2_ad(i, k) + pt1_ad(k)
3041 w2_ad(i, k) = w2_ad(i, k) + dm(k)*wm_ad(k)
3042 dm_ad(k) = dm_ad(k) + w2(i, k)*wm_ad(k)
3045 dm2_ad(i, k) = dm2_ad(i, k) + dm_ad(k)
3048 dz2_ad(i, k) = dz2_ad(i, k) + dz_ad(k)
3053 ws_ad(i) = ws_ad(i) + 2.*ws2_ad(i)
3057 SUBROUTINE rim_2d(ms, bdt, is, ie, km, rgas, gama, gm2, pe2, dm2, pm2&
3058 & , w2, dz2, pt2, ws, c_core)
3060 INTEGER,
INTENT(IN) :: ms, is, ie, km
3061 REAL,
INTENT(IN) :: bdt, gama, rgas
3062 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pm2, gm2
3063 LOGICAL,
INTENT(IN) :: c_core
3064 REAL,
INTENT(IN) :: pt2(is:ie, km)
3065 REAL,
INTENT(IN) :: ws(is:ie)
3067 REAL,
INTENT(INOUT) :: dz2(is:ie, km)
3068 REAL,
INTENT(INOUT) :: w2(is:ie, km)
3069 REAL,
INTENT(OUT) :: pe2(is:ie, km+1)
3072 REAL,
DIMENSION(km+1) :: m_bot, m_top, r_bot, r_top, pe1, pbar, wbar
3073 REAL,
DIMENSION(km) :: r_hi, r_lo, dz, wm, dm, dts
3074 REAL,
DIMENSION(km) :: pf1, wc, cm, pp, pt1
3075 REAL :: dt, rdt, grg, z_frac, ptmp1, rden, pf, time_left
3077 INTEGER :: i, k, n, ke, kt1, ktop
3097 wm(k) = w2(i, k)*dm(k)
3103 IF (ms .GT. 1 .AND. ms .LT. 8)
THEN 3106 rden = -(rgas*dm(k)/dz(k))
3107 pf1(k) = exp(gama*log(rden*pt1(k)))
3108 dts(k) = -(dz(k)/sqrt(grg*pf1(k)/rden))
3109 IF (bdt .GT. dts(k))
THEN 3115 222
IF (ks0 .NE. 1)
THEN 3117 cm(k) = dm(k)/dts(k)
3118 wc(k) = wm(k)/dts(k)
3119 pp(k) = pf1(k) - pm2(i, k)
3121 wbar(1) = (wc(1)+pp(1))/cm(1)
3123 wbar(k) = (wc(k-1)+wc(k)+pp(k)-pp(k-1))/(cm(k-1)+cm(k))
3124 pbar(k) = bdt*(cm(k-1)*wbar(k)-wc(k-1)+pp(k-1))
3127 IF (ks0 .EQ. km)
THEN 3128 pbar(km+1) = bdt*(cm(km)*wbar(km+1)-wc(km)+pp(km))
3131 dz2(i, k) = dz(k) + bdt*(wbar(k+1)-wbar(k))
3135 dz2(i, k) = dz(k) + bdt*(wbar(k+1)-wbar(k))
3136 w2(i, k) = (wm(k)+pbar(k+1)-pbar(k))/dm(k)
3141 pe2(i, k) = pbar(k)*rdt
3148 dz2(i, k) = dz(k) + bdt*(wbar(k+1)-wbar(k))
3152 dz2(i, k) = dz(k) + bdt*(wbar(k+1)-wbar(k))
3153 w2(i, k) = (wm(k)+pbar(k+1)-pbar(k))/dm(k)
3156 pbar(ks0) = pbar(ks0)/
REAL(ms)
3163 rden = -(rgas*dm(k)/dz(k))
3164 pf = exp(gama*log(rden*pt1(k)))
3165 dts(k) = -(dz(k)/sqrt(grg*pf/rden))
3166 ptmp1 = dts(k)*(pf-pm2(i, k))
3167 r_lo(k) = wm(k) + ptmp1
3168 r_hi(k) = wm(k) - ptmp1
3172 IF (dt .GT. dts(k))
THEN 3178 333
IF (ktop .GE. ks1)
THEN 3181 r_bot(k) = z_frac*r_lo(k)
3182 r_top(k+1) = z_frac*r_hi(k)
3183 m_bot(k) = z_frac*dm(k)
3184 m_top(k+1) = m_bot(k)
3186 IF (ktop .EQ. km)
GOTO 666
3192 IF (1 .LT. ktop)
THEN 3197 DO ke=km+1,ktop+2,-1
3200 IF (time_left .GT. dts(k))
THEN 3201 time_left = time_left - dts(k)
3202 m_top(ke) = m_top(ke) + dm(k)
3203 r_top(ke) = r_top(ke) + r_hi(k)
3205 z_frac = time_left/dts(k)
3206 m_top(ke) = m_top(ke) + z_frac*dm(k)
3207 r_top(ke) = r_top(ke) + z_frac*r_hi(k)
3221 IF (time_left .GT. dts(k))
THEN 3222 time_left = time_left - dts(k)
3223 m_bot(ke) = m_bot(ke) + dm(k)
3224 r_bot(ke) = r_bot(ke) + r_lo(k)
3226 z_frac = time_left/dts(k)
3227 m_bot(ke) = m_bot(ke) + z_frac*dm(k)
3228 r_bot(ke) = r_bot(ke) + z_frac*r_lo(k)
3235 IF (time_left .GT. dts(k))
THEN 3236 time_left = time_left - dts(k)
3237 m_bot(ke) = m_bot(ke) + dm(k)
3238 r_bot(ke) = r_bot(ke) - r_hi(k)
3240 z_frac = time_left/dts(k)
3241 m_bot(ke) = m_bot(ke) + z_frac*dm(k)
3242 r_bot(ke) = r_bot(ke) - z_frac*r_hi(k) + (m_bot(ke)-m_surf&
3250 666
IF (ks1 .EQ. 1) wbar(1) = r_bot(1)/m_bot(1)
3252 wbar(k) = (r_bot(k)+r_top(k))/(m_top(k)+m_bot(k))
3256 pbar(k) = m_top(k)*wbar(k) - r_top(k)
3257 pe1(k) = pe1(k) + pbar(k)
3262 dz2(i, k) = dz(k) + dt*(wbar(k+1)-wbar(k))
3266 dz2(i, k) = dz(k) + dt*(wbar(k+1)-wbar(k))
3267 w2(i, k) = (wm(k)+pbar(k+1)-pbar(k))/dm(k)
3272 dz(k) = dz(k) + dt*(wbar(k+1)-wbar(k))
3273 wm(k) = wm(k) + pbar(k+1) - pbar(k)
3279 pe2(i, k) = pe1(k)*rdt
3304 SUBROUTINE sim3_solver_fwd(dt, is, ie, km, rgas, gama, kappa, pe2, dm&
3305 & , pem, w2, dz2, pt2, ws, alpha, p_fac, scale_m)
3307 INTEGER,
INTENT(IN) :: is, ie, km
3308 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, alpha, p_fac, scale_m
3309 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm, pt2
3310 REAL,
INTENT(IN) :: ws(is:ie)
3311 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
3312 REAL :: pe2(is:ie, km+1)
3313 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
3315 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, wk, g_rat, gam
3316 REAL,
DIMENSION(is:ie, km+1) :: pp
3317 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
3318 REAL :: beta, t2, t1g, rdt, ra, capa1, r2g, r6g
3347 t1g = gama*2.*(alpha*dt)**2
3356 aa(i, k) = exp(gama*log(-(dm(i, k)/dz2(i, k)*rgas*pt2(i, k))))
3362 g_rat(i, k) = dm(i, k)/dm(i, k+1)
3363 bb(i, k) = 2.*(1.+g_rat(i, k))
3364 dd(i, k) = 3.*(aa(i, k)+g_rat(i, k)*aa(i, k+1))
3372 pe2(i, 1) = pem(i, 1)
3374 pe2(i, 2) = (dd(i, 1)-pem(i, 1))/bet(i)
3378 dd(i, km) = 3.*aa(i, km) + r2g*dm(i, km)
3382 gam(i, k) = g_rat(i, k-1)/bet(i)
3384 bet(i) = bb(i, k) - gam(i, k)
3386 pe2(i, k+1) = (dd(i, k)-pe2(i, k))/bet(i)
3392 pe2(i, k) = pe2(i, k) - gam(i, k)*pe2(i, k+1)
3399 pp(i, k) = pe2(i, k) - pem(i, k)
3405 aa(i, k) = t1g/(dz2(i, k-1)+dz2(i, k))*pe2(i, k)
3406 wk(i, k) = t2*aa(i, k)*(w1(i, k-1)-w1(i, k))
3408 aa(i, k) = aa(i, k) - scale_m*dm(i, 1)
3413 bet(i) = dm(i, 1) - aa(i, 2)
3415 w2(i, 1) = (dm(i, 1)*w1(i, 1)+dt*pp(i, 2)+wk(i, 2))/bet(i)
3420 gam(i, k) = aa(i, k)/bet(i)
3422 bet(i) = dm(i, k) - (aa(i, k)+aa(i, k+1)+aa(i, k)*gam(i, k))
3424 w2(i, k) = (dm(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, k))+wk(i, k+1&
3425 & )-wk(i, k)-aa(i, k)*w2(i, k-1))/bet(i)
3429 wk1(i) = t1g/dz2(i, km)*pe2(i, km+1)
3431 gam(i, km) = aa(i, km)/bet(i)
3433 bet(i) = dm(i, km) - (aa(i, km)+wk1(i)+aa(i, km)*gam(i, km))
3435 w2(i, km) = (dm(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i, km))-wk(i, &
3436 & km)+wk1(i)*(t2*w1(i, km)-ra*ws(i))-aa(i, km)*w2(i, km-1))/bet(i)
3441 w2(i, k) = w2(i, k) - gam(i, k+1)*w2(i, k+1)
3452 pe2(i, k+1) = pe2(i, k) + (dm(i, k)*(w2(i, k)-w1(i, k))*rdt-beta&
3453 & *(pp(i, k+1)-pp(i, k)))*ra
3459 pe2(i, 1) = pem(i, 1)
3463 IF (p_fac*pem(i, k) .LT. pe2(i, k) + pem(i, k))
THEN 3465 pe2(i, k) = pe2(i, k) + pem(i, k)
3469 pe2(i, k) = p_fac*pem(i, k)
3476 p1(i) = (pe2(i, km)+2.*pe2(i, km+1))*
r3 - r6g*dm(i, km)
3478 dz2(i, km) = -(dm(i, km)*rgas*pt2(i, km)*exp(capa1*log(p1(i))))
3483 p1(i) = (pe2(i, k)+bb(i, k)*pe2(i, k+1)+g_rat(i, k)*pe2(i, k+2))&
3484 & *
r3 - g_rat(i, k)*p1(i)
3486 dz2(i, k) = -(dm(i, k)*rgas*pt2(i, k)*exp(capa1*log(p1(i))))
3492 pe2(i, k) = pe2(i, k) - pem(i, k)
3493 pe2(i, k) = pe2(i, k) + beta*(pp(i, k)-pe2(i, k))
3537 & pe2_ad, dm, dm_ad, pem, pem_ad, w2, w2_ad, dz2, dz2_ad, pt2, pt2_ad&
3538 & , ws, ws_ad, alpha, p_fac, scale_m)
3540 INTEGER,
INTENT(IN) :: is, ie, km
3541 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, alpha, p_fac, scale_m
3542 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm, pt2
3543 REAL,
DIMENSION(is:ie, km) :: dm_ad, pt2_ad
3544 REAL,
INTENT(IN) :: ws(is:ie)
3545 REAL :: ws_ad(is:ie)
3546 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
3547 REAL,
DIMENSION(is:ie, km+1) :: pem_ad
3548 REAL :: pe2(is:ie, km+1)
3549 REAL :: pe2_ad(is:ie, km+1)
3550 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
3551 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2_ad, w2_ad
3552 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, wk, g_rat, gam
3553 REAL,
DIMENSION(is:ie, km) :: aa_ad, bb_ad, dd_ad, w1_ad, wk_ad, &
3555 REAL,
DIMENSION(is:ie, km+1) :: pp
3556 REAL,
DIMENSION(is:ie, km+1) :: pp_ad
3557 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
3558 REAL,
DIMENSION(is:ie) :: p1_ad, wk1_ad, bet_ad
3559 REAL :: beta, t2, t1g, rdt, ra, capa1, r2g, r6g
3636 pp_ad(i, k) = pp_ad(i, k) + beta*pe2_ad(i, k)
3637 pe2_ad(i, k) = (1.0-beta)*pe2_ad(i, k)
3639 pem_ad(i, k) = pem_ad(i, k) - pe2_ad(i, k)
3648 temp5 = capa1*log(p1(i))
3649 temp_ad17 = -(rgas*exp(temp5)*dz2_ad(i, k))
3650 dm_ad(i, k) = dm_ad(i, k) + pt2(i, k)*temp_ad17
3651 pt2_ad(i, k) = pt2_ad(i, k) + dm(i, k)*temp_ad17
3652 p1_ad(i) = p1_ad(i) - capa1*exp(temp5)*dm(i, k)*pt2(i, k)*rgas*&
3653 & dz2_ad(i, k)/p1(i)
3656 temp_ad18 =
r3*p1_ad(i)
3657 pe2_ad(i, k) = pe2_ad(i, k) + temp_ad18
3658 bb_ad(i, k) = bb_ad(i, k) + pe2(i, k+1)*temp_ad18
3659 pe2_ad(i, k+1) = pe2_ad(i, k+1) + bb(i, k)*temp_ad18
3660 g_rat_ad(i, k) = g_rat_ad(i, k) + pe2(i, k+2)*temp_ad18 - p1(i)*&
3662 pe2_ad(i, k+2) = pe2_ad(i, k+2) + g_rat(i, k)*temp_ad18
3663 p1_ad(i) = -(g_rat(i, k)*p1_ad(i))
3668 temp4 = capa1*log(p1(i))
3669 temp_ad16 = -(rgas*exp(temp4)*dz2_ad(i, km))
3670 pt2_ad(i, km) = pt2_ad(i, km) + dm(i, km)*temp_ad16
3671 p1_ad(i) = p1_ad(i) - capa1*exp(temp4)*dm(i, km)*pt2(i, km)*rgas*&
3672 & dz2_ad(i, km)/p1(i)
3673 dm_ad(i, km) = dm_ad(i, km) + pt2(i, km)*temp_ad16 - r6g*p1_ad(i)
3675 pe2_ad(i, km) = pe2_ad(i, km) +
r3*p1_ad(i)
3676 pe2_ad(i, km+1) = pe2_ad(i, km+1) +
r3*2.*p1_ad(i)
3682 IF (branch .EQ. 0)
THEN 3684 pem_ad(i, k) = pem_ad(i, k) + pe2_ad(i, k)
3687 pem_ad(i, k) = pem_ad(i, k) + p_fac*pe2_ad(i, k)
3694 pem_ad(i, 1) = pem_ad(i, 1) + pe2_ad(i, 1)
3701 temp_ad14 = ra*pe2_ad(i, k+1)
3702 temp_ad15 = rdt*dm(i, k)*temp_ad14
3703 pe2_ad(i, k) = pe2_ad(i, k) + pe2_ad(i, k+1)
3704 dm_ad(i, k) = dm_ad(i, k) + rdt*(w2(i, k)-w1(i, k))*temp_ad14
3705 w2_ad(i, k) = w2_ad(i, k) + temp_ad15
3706 w1_ad(i, k) = w1_ad(i, k) - temp_ad15
3707 pp_ad(i, k+1) = pp_ad(i, k+1) - beta*temp_ad14
3708 pp_ad(i, k) = pp_ad(i, k) + beta*temp_ad14
3709 pe2_ad(i, k+1) = 0.0
3720 gam_ad(i, k+1) = gam_ad(i, k+1) - w2(i, k+1)*w2_ad(i, k)
3721 w2_ad(i, k+1) = w2_ad(i, k+1) - gam(i, k+1)*w2_ad(i, k)
3730 temp_ad11 = w2_ad(i, km)/bet(i)
3731 temp3 = t2*w1(i, km) - ra*ws(i)
3732 w1_ad(i, km) = w1_ad(i, km) + (wk1(i)*t2+dm(i, km))*temp_ad11
3733 pp_ad(i, km+1) = pp_ad(i, km+1) + dt*temp_ad11
3734 pp_ad(i, km) = pp_ad(i, km) - dt*temp_ad11
3735 wk_ad(i, km) = wk_ad(i, km) - temp_ad11
3736 ws_ad(i) = ws_ad(i) - wk1(i)*ra*temp_ad11
3737 w2_ad(i, km-1) = w2_ad(i, km-1) - aa(i, km)*temp_ad11
3738 bet_ad(i) = bet_ad(i) - (dm(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i&
3739 & , km))-wk(i, km)+wk1(i)*temp3-aa(i, km)*w2(i, km-1))*temp_ad11/&
3741 dm_ad(i, km) = dm_ad(i, km) + bet_ad(i) + w1(i, km)*temp_ad11
3742 wk1_ad(i) = wk1_ad(i) + temp3*temp_ad11 - bet_ad(i)
3745 gam_ad(i, km) = gam_ad(i, km) - aa(i, km)*bet_ad(i)
3746 temp_ad12 = gam_ad(i, km)/bet(i)
3747 aa_ad(i, km) = aa_ad(i, km) + ((-1.0)-gam(i, km))*bet_ad(i) + &
3748 & temp_ad12 - w2(i, km-1)*temp_ad11
3749 bet_ad(i) = -(aa(i, km)*temp_ad12/bet(i))
3752 temp_ad13 = t1g*wk1_ad(i)/dz2(i, km)
3753 pe2_ad(i, km+1) = pe2_ad(i, km+1) + temp_ad13
3754 dz2_ad(i, km) = dz2_ad(i, km) - pe2(i, km+1)*temp_ad13/dz2(i, km)
3760 temp_ad9 = w2_ad(i, k)/bet(i)
3761 w1_ad(i, k) = w1_ad(i, k) + dm(i, k)*temp_ad9
3762 pp_ad(i, k+1) = pp_ad(i, k+1) + dt*temp_ad9
3763 pp_ad(i, k) = pp_ad(i, k) - dt*temp_ad9
3764 wk_ad(i, k+1) = wk_ad(i, k+1) + temp_ad9
3765 wk_ad(i, k) = wk_ad(i, k) - temp_ad9
3766 w2_ad(i, k-1) = w2_ad(i, k-1) - aa(i, k)*temp_ad9
3767 bet_ad(i) = bet_ad(i) - (dm(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, &
3768 & k))+wk(i, k+1)-wk(i, k)-aa(i, k)*w2(i, k-1))*temp_ad9/bet(i)
3769 dm_ad(i, k) = dm_ad(i, k) + bet_ad(i) + w1(i, k)*temp_ad9
3770 aa_ad(i, k) = aa_ad(i, k) + ((-1.0)-gam(i, k))*bet_ad(i) - w2(i&
3774 aa_ad(i, k+1) = aa_ad(i, k+1) - bet_ad(i)
3775 gam_ad(i, k) = gam_ad(i, k) - aa(i, k)*bet_ad(i)
3777 temp_ad10 = gam_ad(i, k)/bet(i)
3778 bet_ad(i) = -(aa(i, k)*temp_ad10/bet(i))
3779 aa_ad(i, k) = aa_ad(i, k) + temp_ad10
3785 temp_ad8 = w2_ad(i, 1)/bet(i)
3786 w1_ad(i, 1) = w1_ad(i, 1) + dm(i, 1)*temp_ad8
3787 pp_ad(i, 2) = pp_ad(i, 2) + dt*temp_ad8
3788 wk_ad(i, 2) = wk_ad(i, 2) + temp_ad8
3789 bet_ad(i) = bet_ad(i) - (dm(i, 1)*w1(i, 1)+dt*pp(i, 2)+wk(i, 2))*&
3791 dm_ad(i, 1) = dm_ad(i, 1) + bet_ad(i) + w1(i, 1)*temp_ad8
3794 aa_ad(i, 2) = aa_ad(i, 2) - bet_ad(i)
3800 dm_ad(i, 1) = dm_ad(i, 1) - scale_m*aa_ad(i, k)
3801 temp_ad5 = t2*aa(i, k)*wk_ad(i, k)
3802 aa_ad(i, k) = aa_ad(i, k) + t2*(w1(i, k-1)-w1(i, k))*wk_ad(i, k)
3803 w1_ad(i, k-1) = w1_ad(i, k-1) + temp_ad5
3804 w1_ad(i, k) = w1_ad(i, k) - temp_ad5
3807 temp2 = dz2(i, k-1) + dz2(i, k)
3808 temp_ad6 = t1g*aa_ad(i, k)/temp2
3809 temp_ad7 = -(pe2(i, k)*temp_ad6/temp2)
3810 pe2_ad(i, k) = pe2_ad(i, k) + temp_ad6
3811 dz2_ad(i, k-1) = dz2_ad(i, k-1) + temp_ad7
3812 dz2_ad(i, k) = dz2_ad(i, k) + temp_ad7
3818 pe2_ad(i, k) = pe2_ad(i, k) + pp_ad(i, k)
3819 pem_ad(i, k) = pem_ad(i, k) - pp_ad(i, k)
3826 gam_ad(i, k) = gam_ad(i, k) - pe2(i, k+1)*pe2_ad(i, k)
3827 pe2_ad(i, k+1) = pe2_ad(i, k+1) - gam(i, k)*pe2_ad(i, k)
3834 temp_ad3 = pe2_ad(i, k+1)/bet(i)
3835 dd_ad(i, k) = dd_ad(i, k) + temp_ad3
3836 pe2_ad(i, k) = pe2_ad(i, k) - temp_ad3
3837 bet_ad(i) = bet_ad(i) - (dd(i, k)-pe2(i, k))*temp_ad3/bet(i)
3838 pe2_ad(i, k+1) = 0.0
3840 bb_ad(i, k) = bb_ad(i, k) + bet_ad(i)
3841 gam_ad(i, k) = gam_ad(i, k) - bet_ad(i)
3842 temp_ad4 = gam_ad(i, k)/bet(i)
3843 bet_ad(i) = -(g_rat(i, k-1)*temp_ad4/bet(i))
3844 g_rat_ad(i, k-1) = g_rat_ad(i, k-1) + temp_ad4
3850 aa_ad(i, km) = aa_ad(i, km) + 3.*dd_ad(i, km)
3851 dm_ad(i, km) = dm_ad(i, km) + r2g*dd_ad(i, km)
3855 temp_ad2 = pe2_ad(i, 2)/bet(i)
3856 dd_ad(i, 1) = dd_ad(i, 1) + temp_ad2
3857 bet_ad(i) = bet_ad(i) - (dd(i, 1)-pem(i, 1))*temp_ad2/bet(i)
3859 pem_ad(i, 1) = pem_ad(i, 1) + pe2_ad(i, 1) - temp_ad2
3862 bb_ad(i, 1) = bb_ad(i, 1) + bet_ad(i)
3867 temp_ad0 = 3.*dd_ad(i, k)
3868 aa_ad(i, k) = aa_ad(i, k) + temp_ad0
3869 g_rat_ad(i, k) = g_rat_ad(i, k) + 2.*bb_ad(i, k) + aa(i, k+1)*&
3871 aa_ad(i, k+1) = aa_ad(i, k+1) + g_rat(i, k)*temp_ad0
3874 temp_ad1 = g_rat_ad(i, k)/dm(i, k+1)
3875 dm_ad(i, k) = dm_ad(i, k) + temp_ad1
3876 dm_ad(i, k+1) = dm_ad(i, k+1) - dm(i, k)*temp_ad1/dm(i, k+1)
3877 g_rat_ad(i, k) = 0.0
3883 temp0 = dm(i, k)*pt2(i, k)
3885 temp_ad = gama*exp(gama*log(-(rgas*temp)))*aa_ad(i, k)/(temp*&
3887 dm_ad(i, k) = dm_ad(i, k) + pt2(i, k)*temp_ad
3888 pt2_ad(i, k) = pt2_ad(i, k) + dm(i, k)*temp_ad
3889 dz2_ad(i, k) = dz2_ad(i, k) - temp*temp_ad
3891 w2_ad(i, k) = w2_ad(i, k) + w1_ad(i, k)
3896 SUBROUTINE sim3_solver(dt, is, ie, km, rgas, gama, kappa, pe2, dm, pem&
3897 & , w2, dz2, pt2, ws, alpha, p_fac, scale_m)
3899 INTEGER,
INTENT(IN) :: is, ie, km
3900 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, alpha, p_fac, scale_m
3901 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm, pt2
3902 REAL,
INTENT(IN) :: ws(is:ie)
3903 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
3904 REAL,
INTENT(OUT) :: pe2(is:ie, km+1)
3905 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
3907 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, wk, g_rat, gam
3908 REAL,
DIMENSION(is:ie, km+1) :: pp
3909 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
3910 REAL :: beta, t2, t1g, rdt, ra, capa1, r2g, r6g
3918 t1g = gama*2.*(alpha*dt)**2
3927 aa(i, k) = exp(gama*log(-(dm(i, k)/dz2(i, k)*rgas*pt2(i, k))))
3933 g_rat(i, k) = dm(i, k)/dm(i, k+1)
3934 bb(i, k) = 2.*(1.+g_rat(i, k))
3935 dd(i, k) = 3.*(aa(i, k)+g_rat(i, k)*aa(i, k+1))
3942 pe2(i, 1) = pem(i, 1)
3943 pe2(i, 2) = (dd(i, 1)-pem(i, 1))/bet(i)
3946 dd(i, km) = 3.*aa(i, km) + r2g*dm(i, km)
3950 gam(i, k) = g_rat(i, k-1)/bet(i)
3951 bet(i) = bb(i, k) - gam(i, k)
3952 pe2(i, k+1) = (dd(i, k)-pe2(i, k))/bet(i)
3957 pe2(i, k) = pe2(i, k) - gam(i, k)*pe2(i, k+1)
3964 pp(i, k) = pe2(i, k) - pem(i, k)
3969 aa(i, k) = t1g/(dz2(i, k-1)+dz2(i, k))*pe2(i, k)
3970 wk(i, k) = t2*aa(i, k)*(w1(i, k-1)-w1(i, k))
3971 aa(i, k) = aa(i, k) - scale_m*dm(i, 1)
3975 bet(i) = dm(i, 1) - aa(i, 2)
3976 w2(i, 1) = (dm(i, 1)*w1(i, 1)+dt*pp(i, 2)+wk(i, 2))/bet(i)
3980 gam(i, k) = aa(i, k)/bet(i)
3981 bet(i) = dm(i, k) - (aa(i, k)+aa(i, k+1)+aa(i, k)*gam(i, k))
3982 w2(i, k) = (dm(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, k))+wk(i, k+1&
3983 & )-wk(i, k)-aa(i, k)*w2(i, k-1))/bet(i)
3987 wk1(i) = t1g/dz2(i, km)*pe2(i, km+1)
3988 gam(i, km) = aa(i, km)/bet(i)
3989 bet(i) = dm(i, km) - (aa(i, km)+wk1(i)+aa(i, km)*gam(i, km))
3990 w2(i, km) = (dm(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i, km))-wk(i, &
3991 & km)+wk1(i)*(t2*w1(i, km)-ra*ws(i))-aa(i, km)*w2(i, km-1))/bet(i)
3995 w2(i, k) = w2(i, k) - gam(i, k+1)*w2(i, k+1)
4004 pe2(i, k+1) = pe2(i, k) + (dm(i, k)*(w2(i, k)-w1(i, k))*rdt-beta&
4005 & *(pp(i, k+1)-pp(i, k)))*ra
4010 pe2(i, 1) = pem(i, 1)
4014 IF (p_fac*pem(i, k) .LT. pe2(i, k) + pem(i, k))
THEN 4015 pe2(i, k) = pe2(i, k) + pem(i, k)
4017 pe2(i, k) = p_fac*pem(i, k)
4023 p1(i) = (pe2(i, km)+2.*pe2(i, km+1))*
r3 - r6g*dm(i, km)
4024 dz2(i, km) = -(dm(i, km)*rgas*pt2(i, km)*exp(capa1*log(p1(i))))
4028 p1(i) = (pe2(i, k)+bb(i, k)*pe2(i, k+1)+g_rat(i, k)*pe2(i, k+2))&
4029 & *
r3 - g_rat(i, k)*p1(i)
4030 dz2(i, k) = -(dm(i, k)*rgas*pt2(i, k)*exp(capa1*log(p1(i))))
4035 pe2(i, k) = pe2(i, k) - pem(i, k)
4036 pe2(i, k) = pe2(i, k) + beta*(pp(i, k)-pe2(i, k))
4061 & dm, pem, w2, dz2, pt2, ws, p_fac, scale_m)
4064 INTEGER,
INTENT(IN) :: is, ie, km
4065 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac, scale_m
4066 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm, pt2
4067 REAL,
INTENT(IN) :: ws(is:ie)
4068 REAL,
INTENT(IN) :: pem(is:ie, km+1)
4069 REAL :: pe2(is:ie, km+1)
4070 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
4072 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, g_rat, gam
4073 REAL,
DIMENSION(is:ie, km+1) :: pp
4074 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
4075 REAL :: t1g, rdt, capa1, r2g, r6g
4106 aa(i, k) = exp(gama*log(-(dm(i, k)/dz2(i, k)*rgas*pt2(i, k))))
4112 g_rat(i, k) = dm(i, k)/dm(i, k+1)
4113 bb(i, k) = 2.*(1.+g_rat(i, k))
4114 dd(i, k) = 3.*(aa(i, k)+g_rat(i, k)*aa(i, k+1))
4122 pe2(i, 1) = pem(i, 1)
4124 pe2(i, 2) = (dd(i, 1)-pem(i, 1))/bet(i)
4128 dd(i, km) = 3.*aa(i, km) + r2g*dm(i, km)
4132 gam(i, k) = g_rat(i, k-1)/bet(i)
4134 bet(i) = bb(i, k) - gam(i, k)
4136 pe2(i, k+1) = (dd(i, k)-pe2(i, k))/bet(i)
4142 pe2(i, k) = pe2(i, k) - gam(i, k)*pe2(i, k+1)
4149 pp(i, k) = pe2(i, k) - pem(i, k)
4155 aa(i, k) = t1g/(dz2(i, k-1)+dz2(i, k))*pe2(i, k) - scale_m*dm(i&
4161 bet(i) = dm(i, 1) - aa(i, 2)
4163 w2(i, 1) = (dm(i, 1)*w1(i, 1)+dt*pp(i, 2))/bet(i)
4168 gam(i, k) = aa(i, k)/bet(i)
4170 bet(i) = dm(i, k) - (aa(i, k)+aa(i, k+1)+aa(i, k)*gam(i, k))
4172 w2(i, k) = (dm(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, k))-aa(i, k)*&
4173 & w2(i, k-1))/bet(i)
4177 wk1(i) = t1g/dz2(i, km)*pe2(i, km+1)
4179 gam(i, km) = aa(i, km)/bet(i)
4181 bet(i) = dm(i, km) - (aa(i, km)+wk1(i)+aa(i, km)*gam(i, km))
4183 w2(i, km) = (dm(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i, km))-wk1(i)&
4184 & *ws(i)-aa(i, km)*w2(i, km-1))/bet(i)
4189 w2(i, k) = w2(i, k) - gam(i, k+1)*w2(i, k+1)
4200 pe2(i, k+1) = pe2(i, k) + dm(i, k)*(w2(i, k)-w1(i, k))*rdt
4206 pe2(i, 1) = pem(i, 1)
4210 IF (p_fac*pem(i, k) .LT. pe2(i, k) + pem(i, k))
THEN 4212 pe2(i, k) = pe2(i, k) + pem(i, k)
4216 pe2(i, k) = p_fac*pem(i, k)
4223 p1(i) = (pe2(i, km)+2.*pe2(i, km+1))*
r3 - r6g*dm(i, km)
4225 dz2(i, km) = -(dm(i, km)*rgas*pt2(i, km)*exp(capa1*log(p1(i))))
4230 p1(i) = (pe2(i, k)+bb(i, k)*pe2(i, k+1)+g_rat(i, k)*pe2(i, k+2))&
4231 & *
r3 - g_rat(i, k)*p1(i)
4233 dz2(i, k) = -(dm(i, k)*rgas*pt2(i, k)*exp(capa1*log(p1(i))))
4239 pe2(i, k) = pe2(i, k) - pem(i, k)
4279 & pe2_ad, dm, dm_ad, pem, pem_ad, w2, w2_ad, dz2, dz2_ad, pt2, pt2_ad&
4280 & , ws, ws_ad, p_fac, scale_m)
4282 INTEGER,
INTENT(IN) :: is, ie, km
4283 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac, scale_m
4284 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm, pt2
4285 REAL,
DIMENSION(is:ie, km) :: dm_ad, pt2_ad
4286 REAL,
INTENT(IN) :: ws(is:ie)
4287 REAL :: ws_ad(is:ie)
4288 REAL,
INTENT(IN) :: pem(is:ie, km+1)
4289 REAL :: pem_ad(is:ie, km+1)
4290 REAL :: pe2(is:ie, km+1)
4291 REAL :: pe2_ad(is:ie, km+1)
4292 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
4293 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2_ad, w2_ad
4294 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, g_rat, gam
4295 REAL,
DIMENSION(is:ie, km) :: aa_ad, bb_ad, dd_ad, w1_ad, g_rat_ad, &
4297 REAL,
DIMENSION(is:ie, km+1) :: pp
4298 REAL,
DIMENSION(is:ie, km+1) :: pp_ad
4299 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
4300 REAL,
DIMENSION(is:ie) :: p1_ad, wk1_ad, bet_ad
4301 REAL :: t1g, rdt, capa1, r2g, r6g
4367 pem_ad(i, k) = pem_ad(i, k) - pe2_ad(i, k)
4376 temp4 = capa1*log(p1(i))
4377 temp_ad15 = -(rgas*exp(temp4)*dz2_ad(i, k))
4378 dm_ad(i, k) = dm_ad(i, k) + pt2(i, k)*temp_ad15
4379 pt2_ad(i, k) = pt2_ad(i, k) + dm(i, k)*temp_ad15
4380 p1_ad(i) = p1_ad(i) - capa1*exp(temp4)*dm(i, k)*pt2(i, k)*rgas*&
4381 & dz2_ad(i, k)/p1(i)
4384 temp_ad16 =
r3*p1_ad(i)
4385 pe2_ad(i, k) = pe2_ad(i, k) + temp_ad16
4386 bb_ad(i, k) = bb_ad(i, k) + pe2(i, k+1)*temp_ad16
4387 pe2_ad(i, k+1) = pe2_ad(i, k+1) + bb(i, k)*temp_ad16
4388 g_rat_ad(i, k) = g_rat_ad(i, k) + pe2(i, k+2)*temp_ad16 - p1(i)*&
4390 pe2_ad(i, k+2) = pe2_ad(i, k+2) + g_rat(i, k)*temp_ad16
4391 p1_ad(i) = -(g_rat(i, k)*p1_ad(i))
4396 temp3 = capa1*log(p1(i))
4397 temp_ad14 = -(rgas*exp(temp3)*dz2_ad(i, km))
4398 pt2_ad(i, km) = pt2_ad(i, km) + dm(i, km)*temp_ad14
4399 p1_ad(i) = p1_ad(i) - capa1*exp(temp3)*dm(i, km)*pt2(i, km)*rgas*&
4400 & dz2_ad(i, km)/p1(i)
4401 dm_ad(i, km) = dm_ad(i, km) + pt2(i, km)*temp_ad14 - r6g*p1_ad(i)
4403 pe2_ad(i, km) = pe2_ad(i, km) +
r3*p1_ad(i)
4404 pe2_ad(i, km+1) = pe2_ad(i, km+1) +
r3*2.*p1_ad(i)
4410 IF (branch .EQ. 0)
THEN 4412 pem_ad(i, k) = pem_ad(i, k) + pe2_ad(i, k)
4415 pem_ad(i, k) = pem_ad(i, k) + p_fac*pe2_ad(i, k)
4422 pem_ad(i, 1) = pem_ad(i, 1) + pe2_ad(i, 1)
4429 temp_ad13 = rdt*dm(i, k)*pe2_ad(i, k+1)
4430 pe2_ad(i, k) = pe2_ad(i, k) + pe2_ad(i, k+1)
4431 dm_ad(i, k) = dm_ad(i, k) + rdt*(w2(i, k)-w1(i, k))*pe2_ad(i, k+&
4433 w2_ad(i, k) = w2_ad(i, k) + temp_ad13
4434 w1_ad(i, k) = w1_ad(i, k) - temp_ad13
4435 pe2_ad(i, k+1) = 0.0
4446 gam_ad(i, k+1) = gam_ad(i, k+1) - w2(i, k+1)*w2_ad(i, k)
4447 w2_ad(i, k+1) = w2_ad(i, k+1) - gam(i, k+1)*w2_ad(i, k)
4456 temp_ad10 = w2_ad(i, km)/bet(i)
4457 w1_ad(i, km) = w1_ad(i, km) + dm(i, km)*temp_ad10
4458 pp_ad(i, km+1) = pp_ad(i, km+1) + dt*temp_ad10
4459 pp_ad(i, km) = pp_ad(i, km) - dt*temp_ad10
4460 ws_ad(i) = ws_ad(i) - wk1(i)*temp_ad10
4461 w2_ad(i, km-1) = w2_ad(i, km-1) - aa(i, km)*temp_ad10
4462 bet_ad(i) = bet_ad(i) - (dm(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i&
4463 & , km))-wk1(i)*ws(i)-aa(i, km)*w2(i, km-1))*temp_ad10/bet(i)
4464 dm_ad(i, km) = dm_ad(i, km) + bet_ad(i) + w1(i, km)*temp_ad10
4465 wk1_ad(i) = wk1_ad(i) - bet_ad(i) - ws(i)*temp_ad10
4468 gam_ad(i, km) = gam_ad(i, km) - aa(i, km)*bet_ad(i)
4469 temp_ad11 = gam_ad(i, km)/bet(i)
4470 aa_ad(i, km) = aa_ad(i, km) + ((-1.0)-gam(i, km))*bet_ad(i) + &
4471 & temp_ad11 - w2(i, km-1)*temp_ad10
4472 bet_ad(i) = -(aa(i, km)*temp_ad11/bet(i))
4475 temp_ad12 = t1g*wk1_ad(i)/dz2(i, km)
4476 pe2_ad(i, km+1) = pe2_ad(i, km+1) + temp_ad12
4477 dz2_ad(i, km) = dz2_ad(i, km) - pe2(i, km+1)*temp_ad12/dz2(i, km)
4483 temp_ad8 = w2_ad(i, k)/bet(i)
4484 w1_ad(i, k) = w1_ad(i, k) + dm(i, k)*temp_ad8
4485 pp_ad(i, k+1) = pp_ad(i, k+1) + dt*temp_ad8
4486 pp_ad(i, k) = pp_ad(i, k) - dt*temp_ad8
4487 w2_ad(i, k-1) = w2_ad(i, k-1) - aa(i, k)*temp_ad8
4488 bet_ad(i) = bet_ad(i) - (dm(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, &
4489 & k))-aa(i, k)*w2(i, k-1))*temp_ad8/bet(i)
4490 dm_ad(i, k) = dm_ad(i, k) + bet_ad(i) + w1(i, k)*temp_ad8
4491 aa_ad(i, k) = aa_ad(i, k) + ((-1.0)-gam(i, k))*bet_ad(i) - w2(i&
4495 aa_ad(i, k+1) = aa_ad(i, k+1) - bet_ad(i)
4496 gam_ad(i, k) = gam_ad(i, k) - aa(i, k)*bet_ad(i)
4498 temp_ad9 = gam_ad(i, k)/bet(i)
4499 bet_ad(i) = -(aa(i, k)*temp_ad9/bet(i))
4500 aa_ad(i, k) = aa_ad(i, k) + temp_ad9
4506 temp_ad7 = w2_ad(i, 1)/bet(i)
4507 w1_ad(i, 1) = w1_ad(i, 1) + dm(i, 1)*temp_ad7
4508 pp_ad(i, 2) = pp_ad(i, 2) + dt*temp_ad7
4509 bet_ad(i) = bet_ad(i) - (dm(i, 1)*w1(i, 1)+dt*pp(i, 2))*temp_ad7/&
4511 dm_ad(i, 1) = dm_ad(i, 1) + bet_ad(i) + w1(i, 1)*temp_ad7
4514 aa_ad(i, 2) = aa_ad(i, 2) - bet_ad(i)
4520 temp2 = dz2(i, k-1) + dz2(i, k)
4521 temp_ad5 = t1g*aa_ad(i, k)/temp2
4522 temp_ad6 = -(pe2(i, k)*temp_ad5/temp2)
4523 pe2_ad(i, k) = pe2_ad(i, k) + temp_ad5
4524 dz2_ad(i, k-1) = dz2_ad(i, k-1) + temp_ad6
4525 dz2_ad(i, k) = dz2_ad(i, k) + temp_ad6
4526 dm_ad(i, 1) = dm_ad(i, 1) - scale_m*aa_ad(i, k)
4532 pe2_ad(i, k) = pe2_ad(i, k) + pp_ad(i, k)
4533 pem_ad(i, k) = pem_ad(i, k) - pp_ad(i, k)
4540 gam_ad(i, k) = gam_ad(i, k) - pe2(i, k+1)*pe2_ad(i, k)
4541 pe2_ad(i, k+1) = pe2_ad(i, k+1) - gam(i, k)*pe2_ad(i, k)
4548 temp_ad3 = pe2_ad(i, k+1)/bet(i)
4549 dd_ad(i, k) = dd_ad(i, k) + temp_ad3
4550 pe2_ad(i, k) = pe2_ad(i, k) - temp_ad3
4551 bet_ad(i) = bet_ad(i) - (dd(i, k)-pe2(i, k))*temp_ad3/bet(i)
4552 pe2_ad(i, k+1) = 0.0
4554 bb_ad(i, k) = bb_ad(i, k) + bet_ad(i)
4555 gam_ad(i, k) = gam_ad(i, k) - bet_ad(i)
4556 temp_ad4 = gam_ad(i, k)/bet(i)
4557 bet_ad(i) = -(g_rat(i, k-1)*temp_ad4/bet(i))
4558 g_rat_ad(i, k-1) = g_rat_ad(i, k-1) + temp_ad4
4564 aa_ad(i, km) = aa_ad(i, km) + 3.*dd_ad(i, km)
4565 dm_ad(i, km) = dm_ad(i, km) + r2g*dd_ad(i, km)
4569 temp_ad2 = pe2_ad(i, 2)/bet(i)
4570 dd_ad(i, 1) = dd_ad(i, 1) + temp_ad2
4571 bet_ad(i) = bet_ad(i) - (dd(i, 1)-pem(i, 1))*temp_ad2/bet(i)
4573 pem_ad(i, 1) = pem_ad(i, 1) + pe2_ad(i, 1) - temp_ad2
4576 bb_ad(i, 1) = bb_ad(i, 1) + bet_ad(i)
4581 temp_ad0 = 3.*dd_ad(i, k)
4582 aa_ad(i, k) = aa_ad(i, k) + temp_ad0
4583 g_rat_ad(i, k) = g_rat_ad(i, k) + 2.*bb_ad(i, k) + aa(i, k+1)*&
4585 aa_ad(i, k+1) = aa_ad(i, k+1) + g_rat(i, k)*temp_ad0
4588 temp_ad1 = g_rat_ad(i, k)/dm(i, k+1)
4589 dm_ad(i, k) = dm_ad(i, k) + temp_ad1
4590 dm_ad(i, k+1) = dm_ad(i, k+1) - dm(i, k)*temp_ad1/dm(i, k+1)
4591 g_rat_ad(i, k) = 0.0
4597 temp0 = dm(i, k)*pt2(i, k)
4599 temp_ad = gama*exp(gama*log(-(rgas*temp)))*aa_ad(i, k)/(temp*&
4601 dm_ad(i, k) = dm_ad(i, k) + pt2(i, k)*temp_ad
4602 pt2_ad(i, k) = pt2_ad(i, k) + dm(i, k)*temp_ad
4603 dz2_ad(i, k) = dz2_ad(i, k) - temp*temp_ad
4605 w2_ad(i, k) = w2_ad(i, k) + w1_ad(i, k)
4610 SUBROUTINE sim3p0_solver(dt, is, ie, km, rgas, gama, kappa, pe2, dm, &
4611 & pem, w2, dz2, pt2, ws, p_fac, scale_m)
4614 INTEGER,
INTENT(IN) :: is, ie, km
4615 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac, scale_m
4616 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm, pt2
4617 REAL,
INTENT(IN) :: ws(is:ie)
4618 REAL,
INTENT(IN) :: pem(is:ie, km+1)
4619 REAL,
INTENT(OUT) :: pe2(is:ie, km+1)
4620 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
4622 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, g_rat, gam
4623 REAL,
DIMENSION(is:ie, km+1) :: pp
4624 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
4625 REAL :: t1g, rdt, capa1, r2g, r6g
4639 aa(i, k) = exp(gama*log(-(dm(i, k)/dz2(i, k)*rgas*pt2(i, k))))
4645 g_rat(i, k) = dm(i, k)/dm(i, k+1)
4646 bb(i, k) = 2.*(1.+g_rat(i, k))
4647 dd(i, k) = 3.*(aa(i, k)+g_rat(i, k)*aa(i, k+1))
4654 pe2(i, 1) = pem(i, 1)
4655 pe2(i, 2) = (dd(i, 1)-pem(i, 1))/bet(i)
4658 dd(i, km) = 3.*aa(i, km) + r2g*dm(i, km)
4662 gam(i, k) = g_rat(i, k-1)/bet(i)
4663 bet(i) = bb(i, k) - gam(i, k)
4664 pe2(i, k+1) = (dd(i, k)-pe2(i, k))/bet(i)
4669 pe2(i, k) = pe2(i, k) - gam(i, k)*pe2(i, k+1)
4676 pp(i, k) = pe2(i, k) - pem(i, k)
4681 aa(i, k) = t1g/(dz2(i, k-1)+dz2(i, k))*pe2(i, k) - scale_m*dm(i&
4686 bet(i) = dm(i, 1) - aa(i, 2)
4687 w2(i, 1) = (dm(i, 1)*w1(i, 1)+dt*pp(i, 2))/bet(i)
4691 gam(i, k) = aa(i, k)/bet(i)
4692 bet(i) = dm(i, k) - (aa(i, k)+aa(i, k+1)+aa(i, k)*gam(i, k))
4693 w2(i, k) = (dm(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, k))-aa(i, k)*&
4694 & w2(i, k-1))/bet(i)
4698 wk1(i) = t1g/dz2(i, km)*pe2(i, km+1)
4699 gam(i, km) = aa(i, km)/bet(i)
4700 bet(i) = dm(i, km) - (aa(i, km)+wk1(i)+aa(i, km)*gam(i, km))
4701 w2(i, km) = (dm(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i, km))-wk1(i)&
4702 & *ws(i)-aa(i, km)*w2(i, km-1))/bet(i)
4706 w2(i, k) = w2(i, k) - gam(i, k+1)*w2(i, k+1)
4715 pe2(i, k+1) = pe2(i, k) + dm(i, k)*(w2(i, k)-w1(i, k))*rdt
4720 pe2(i, 1) = pem(i, 1)
4724 IF (p_fac*pem(i, k) .LT. pe2(i, k) + pem(i, k))
THEN 4725 pe2(i, k) = pe2(i, k) + pem(i, k)
4727 pe2(i, k) = p_fac*pem(i, k)
4733 p1(i) = (pe2(i, km)+2.*pe2(i, km+1))*
r3 - r6g*dm(i, km)
4734 dz2(i, km) = -(dm(i, km)*rgas*pt2(i, km)*exp(capa1*log(p1(i))))
4738 p1(i) = (pe2(i, k)+bb(i, k)*pe2(i, k+1)+g_rat(i, k)*pe2(i, k+2))&
4739 & *
r3 - g_rat(i, k)*p1(i)
4740 dz2(i, k) = -(dm(i, k)*rgas*pt2(i, k)*exp(capa1*log(p1(i))))
4745 pe2(i, k) = pe2(i, k) - pem(i, k)
4769 SUBROUTINE sim1_solver_fwd(dt, is, ie, km, rgas, gama, gm2, cp2, kappa&
4770 & , pe, dm2, pm2, pem, w2, dz2, pt2, ws, p_fac)
4772 INTEGER,
INTENT(IN) :: is, ie, km
4773 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac
4774 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pt2, pm2, gm2, cp2
4775 REAL,
INTENT(IN) :: ws(is:ie)
4776 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
4777 REAL :: pe(is:ie, km+1)
4778 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
4780 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, g_rat, gam
4781 REAL,
DIMENSION(is:ie, km+1) :: pp
4782 REAL,
DIMENSION(is:ie) :: p1, bet
4783 REAL :: t1g, rdt, capa1
4813 pe(i, k) = exp(gama*log(-(dm2(i, k)/dz2(i, k)*rgas*pt2(i, k)))) &
4819 g_rat(i, k) = dm2(i, k)/dm2(i, k+1)
4820 bb(i, k) = 2.*(1.+g_rat(i, k))
4821 dd(i, k) = 3.*(pe(i, k)+g_rat(i, k)*pe(i, k+1))
4827 pp(i, 2) = dd(i, 1)/bet(i)
4830 dd(i, km) = 3.*pe(i, km)
4834 gam(i, k) = g_rat(i, k-1)/bet(i)
4836 bet(i) = bb(i, k) - gam(i, k)
4838 pp(i, k+1) = (dd(i, k)-pp(i, k))/bet(i)
4844 pp(i, k) = pp(i, k) - gam(i, k)*pp(i, k+1)
4850 aa(i, k) = t1g/(dz2(i, k-1)+dz2(i, k))*(pem(i, k)+pp(i, k))
4855 bet(i) = dm2(i, 1) - aa(i, 2)
4857 w2(i, 1) = (dm2(i, 1)*w1(i, 1)+dt*pp(i, 2))/bet(i)
4862 gam(i, k) = aa(i, k)/bet(i)
4864 bet(i) = dm2(i, k) - (aa(i, k)+aa(i, k+1)+aa(i, k)*gam(i, k))
4866 w2(i, k) = (dm2(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, k))-aa(i, k)&
4867 & *w2(i, k-1))/bet(i)
4871 p1(i) = t1g/dz2(i, km)*(pem(i, km+1)+pp(i, km+1))
4873 gam(i, km) = aa(i, km)/bet(i)
4875 bet(i) = dm2(i, km) - (aa(i, km)+p1(i)+aa(i, km)*gam(i, km))
4877 w2(i, km) = (dm2(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i, km))-p1(i)&
4878 & *ws(i)-aa(i, km)*w2(i, km-1))/bet(i)
4883 w2(i, k) = w2(i, k) - gam(i, k+1)*w2(i, k+1)
4893 pe(i, k+1) = pe(i, k) + dm2(i, k)*(w2(i, k)-w1(i, k))*rdt
4898 p1(i) = (pe(i, km)+2.*pe(i, km+1))*
r3 4899 IF (p_fac*pm2(i, km) .LT. p1(i) + pm2(i, km))
THEN 4901 max1 = p1(i) + pm2(i, km)
4905 max1 = p_fac*pm2(i, km)
4909 dz2(i, km) = -(dm2(i, km)*rgas*pt2(i, km)*exp(capa1*log(max1)))
4914 p1(i) = (pe(i, k)+bb(i, k)*pe(i, k+1)+g_rat(i, k)*pe(i, k+2))*
r3&
4915 & - g_rat(i, k)*p1(i)
4916 IF (p_fac*pm2(i, k) .LT. p1(i) + pm2(i, k))
THEN 4918 max2 = p1(i) + pm2(i, k)
4922 max2 = p_fac*pm2(i, k)
4926 dz2(i, k) = -(dm2(i, k)*rgas*pt2(i, k)*exp(capa1*log(max2)))
4964 SUBROUTINE sim1_solver_bwd(dt, is, ie, km, rgas, gama, gm2, cp2, kappa&
4965 & , pe, pe_ad, dm2, dm2_ad, pm2, pm2_ad, pem, pem_ad, w2, w2_ad, dz2, &
4966 & dz2_ad, pt2, pt2_ad, ws, ws_ad, p_fac)
4968 INTEGER,
INTENT(IN) :: is, ie, km
4969 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac
4970 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pt2, pm2, gm2, cp2
4971 REAL,
DIMENSION(is:ie, km) :: dm2_ad, pt2_ad, pm2_ad
4972 REAL,
INTENT(IN) :: ws(is:ie)
4973 REAL :: ws_ad(is:ie)
4974 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
4975 REAL,
DIMENSION(is:ie, km+1) :: pem_ad
4976 REAL :: pe(is:ie, km+1)
4977 REAL :: pe_ad(is:ie, km+1)
4978 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
4979 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2_ad, w2_ad
4980 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, g_rat, gam
4981 REAL,
DIMENSION(is:ie, km) :: aa_ad, bb_ad, dd_ad, w1_ad, g_rat_ad, &
4983 REAL,
DIMENSION(is:ie, km+1) :: pp
4984 REAL,
DIMENSION(is:ie, km+1) :: pp_ad
4985 REAL,
DIMENSION(is:ie) :: p1, bet
4986 REAL,
DIMENSION(is:ie) :: p1_ad, bet_ad
4987 REAL :: t1g, rdt, capa1
5056 temp4 = capa1*log(max2)
5057 temp_ad16 = -(rgas*exp(temp4)*dz2_ad(i, k))
5058 dm2_ad(i, k) = dm2_ad(i, k) + pt2(i, k)*temp_ad16
5059 pt2_ad(i, k) = pt2_ad(i, k) + dm2(i, k)*temp_ad16
5060 max2_ad = -(capa1*exp(temp4)*dm2(i, k)*pt2(i, k)*rgas*dz2_ad(i, &
5064 IF (branch .EQ. 0)
THEN 5066 p1_ad(i) = p1_ad(i) + max2_ad
5067 pm2_ad(i, k) = pm2_ad(i, k) + max2_ad
5070 pm2_ad(i, k) = pm2_ad(i, k) + p_fac*max2_ad
5073 temp_ad15 =
r3*p1_ad(i)
5074 pe_ad(i, k) = pe_ad(i, k) + temp_ad15
5075 bb_ad(i, k) = bb_ad(i, k) + pe(i, k+1)*temp_ad15
5076 pe_ad(i, k+1) = pe_ad(i, k+1) + bb(i, k)*temp_ad15
5077 g_rat_ad(i, k) = g_rat_ad(i, k) + pe(i, k+2)*temp_ad15 - p1(i)*&
5079 pe_ad(i, k+2) = pe_ad(i, k+2) + g_rat(i, k)*temp_ad15
5080 p1_ad(i) = -(g_rat(i, k)*p1_ad(i))
5085 temp3 = capa1*log(max1)
5086 temp_ad14 = -(rgas*exp(temp3)*dz2_ad(i, km))
5087 dm2_ad(i, km) = dm2_ad(i, km) + pt2(i, km)*temp_ad14
5088 pt2_ad(i, km) = pt2_ad(i, km) + dm2(i, km)*temp_ad14
5089 max1_ad = -(capa1*exp(temp3)*dm2(i, km)*pt2(i, km)*rgas*dz2_ad(i, &
5093 IF (branch .EQ. 0)
THEN 5095 p1_ad(i) = p1_ad(i) + max1_ad
5096 pm2_ad(i, km) = pm2_ad(i, km) + max1_ad
5099 pm2_ad(i, km) = pm2_ad(i, km) + p_fac*max1_ad
5102 pe_ad(i, km) = pe_ad(i, km) +
r3*p1_ad(i)
5103 pe_ad(i, km+1) = pe_ad(i, km+1) +
r3*2.*p1_ad(i)
5110 temp_ad13 = rdt*dm2(i, k)*pe_ad(i, k+1)
5111 pe_ad(i, k) = pe_ad(i, k) + pe_ad(i, k+1)
5112 dm2_ad(i, k) = dm2_ad(i, k) + rdt*(w2(i, k)-w1(i, k))*pe_ad(i, k&
5114 w2_ad(i, k) = w2_ad(i, k) + temp_ad13
5115 w1_ad(i, k) = w1_ad(i, k) - temp_ad13
5127 gam_ad(i, k+1) = gam_ad(i, k+1) - w2(i, k+1)*w2_ad(i, k)
5128 w2_ad(i, k+1) = w2_ad(i, k+1) - gam(i, k+1)*w2_ad(i, k)
5136 temp_ad10 = w2_ad(i, km)/bet(i)
5137 w1_ad(i, km) = w1_ad(i, km) + dm2(i, km)*temp_ad10
5138 pp_ad(i, km+1) = pp_ad(i, km+1) + dt*temp_ad10
5139 pp_ad(i, km) = pp_ad(i, km) - dt*temp_ad10
5140 ws_ad(i) = ws_ad(i) - p1(i)*temp_ad10
5141 w2_ad(i, km-1) = w2_ad(i, km-1) - aa(i, km)*temp_ad10
5142 bet_ad(i) = bet_ad(i) - (dm2(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i&
5143 & , km))-p1(i)*ws(i)-aa(i, km)*w2(i, km-1))*temp_ad10/bet(i)
5144 dm2_ad(i, km) = dm2_ad(i, km) + bet_ad(i) + w1(i, km)*temp_ad10
5145 p1_ad(i) = p1_ad(i) - bet_ad(i) - ws(i)*temp_ad10
5148 gam_ad(i, km) = gam_ad(i, km) - aa(i, km)*bet_ad(i)
5149 temp_ad11 = gam_ad(i, km)/bet(i)
5150 aa_ad(i, km) = aa_ad(i, km) + ((-1.0)-gam(i, km))*bet_ad(i) + &
5151 & temp_ad11 - w2(i, km-1)*temp_ad10
5152 bet_ad(i) = -(aa(i, km)*temp_ad11/bet(i))
5155 temp_ad12 = t1g*p1_ad(i)/dz2(i, km)
5156 pem_ad(i, km+1) = pem_ad(i, km+1) + temp_ad12
5157 pp_ad(i, km+1) = pp_ad(i, km+1) + temp_ad12
5158 dz2_ad(i, km) = dz2_ad(i, km) - (pem(i, km+1)+pp(i, km+1))*&
5159 & temp_ad12/dz2(i, km)
5165 temp_ad8 = w2_ad(i, k)/bet(i)
5166 w1_ad(i, k) = w1_ad(i, k) + dm2(i, k)*temp_ad8
5167 pp_ad(i, k+1) = pp_ad(i, k+1) + dt*temp_ad8
5168 pp_ad(i, k) = pp_ad(i, k) - dt*temp_ad8
5169 w2_ad(i, k-1) = w2_ad(i, k-1) - aa(i, k)*temp_ad8
5170 bet_ad(i) = bet_ad(i) - (dm2(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i&
5171 & , k))-aa(i, k)*w2(i, k-1))*temp_ad8/bet(i)
5172 dm2_ad(i, k) = dm2_ad(i, k) + bet_ad(i) + w1(i, k)*temp_ad8
5173 aa_ad(i, k) = aa_ad(i, k) + ((-1.0)-gam(i, k))*bet_ad(i) - w2(i&
5177 aa_ad(i, k+1) = aa_ad(i, k+1) - bet_ad(i)
5178 gam_ad(i, k) = gam_ad(i, k) - aa(i, k)*bet_ad(i)
5180 temp_ad9 = gam_ad(i, k)/bet(i)
5181 bet_ad(i) = -(aa(i, k)*temp_ad9/bet(i))
5182 aa_ad(i, k) = aa_ad(i, k) + temp_ad9
5188 temp_ad7 = w2_ad(i, 1)/bet(i)
5189 w1_ad(i, 1) = w1_ad(i, 1) + dm2(i, 1)*temp_ad7
5190 pp_ad(i, 2) = pp_ad(i, 2) + dt*temp_ad7
5191 bet_ad(i) = bet_ad(i) - (dm2(i, 1)*w1(i, 1)+dt*pp(i, 2))*temp_ad7/&
5193 dm2_ad(i, 1) = dm2_ad(i, 1) + bet_ad(i) + w1(i, 1)*temp_ad7
5196 aa_ad(i, 2) = aa_ad(i, 2) - bet_ad(i)
5201 temp2 = dz2(i, k-1) + dz2(i, k)
5202 temp_ad5 = t1g*aa_ad(i, k)/temp2
5203 temp_ad6 = -((pem(i, k)+pp(i, k))*temp_ad5/temp2)
5204 pem_ad(i, k) = pem_ad(i, k) + temp_ad5
5205 pp_ad(i, k) = pp_ad(i, k) + temp_ad5
5206 dz2_ad(i, k-1) = dz2_ad(i, k-1) + temp_ad6
5207 dz2_ad(i, k) = dz2_ad(i, k) + temp_ad6
5214 gam_ad(i, k) = gam_ad(i, k) - pp(i, k+1)*pp_ad(i, k)
5215 pp_ad(i, k+1) = pp_ad(i, k+1) - gam(i, k)*pp_ad(i, k)
5222 temp_ad3 = pp_ad(i, k+1)/bet(i)
5223 dd_ad(i, k) = dd_ad(i, k) + temp_ad3
5224 pp_ad(i, k) = pp_ad(i, k) - temp_ad3
5225 bet_ad(i) = bet_ad(i) - (dd(i, k)-pp(i, k))*temp_ad3/bet(i)
5228 bb_ad(i, k) = bb_ad(i, k) + bet_ad(i)
5229 gam_ad(i, k) = gam_ad(i, k) - bet_ad(i)
5230 temp_ad4 = gam_ad(i, k)/bet(i)
5231 bet_ad(i) = -(g_rat(i, k-1)*temp_ad4/bet(i))
5232 g_rat_ad(i, k-1) = g_rat_ad(i, k-1) + temp_ad4
5238 pe_ad(i, km) = pe_ad(i, km) + 3.*dd_ad(i, km)
5241 temp_ad2 = pp_ad(i, 2)/bet(i)
5242 dd_ad(i, 1) = dd_ad(i, 1) + temp_ad2
5243 bet_ad(i) = bet_ad(i) - dd(i, 1)*temp_ad2/bet(i)
5246 bb_ad(i, 1) = bb_ad(i, 1) + bet_ad(i)
5251 temp_ad0 = 3.*dd_ad(i, k)
5252 pe_ad(i, k) = pe_ad(i, k) + temp_ad0
5253 g_rat_ad(i, k) = g_rat_ad(i, k) + 2.*bb_ad(i, k) + pe(i, k+1)*&
5255 pe_ad(i, k+1) = pe_ad(i, k+1) + g_rat(i, k)*temp_ad0
5258 temp_ad1 = g_rat_ad(i, k)/dm2(i, k+1)
5259 dm2_ad(i, k) = dm2_ad(i, k) + temp_ad1
5260 dm2_ad(i, k+1) = dm2_ad(i, k+1) - dm2(i, k)*temp_ad1/dm2(i, k+1)
5261 g_rat_ad(i, k) = 0.0
5268 temp0 = dm2(i, k)*pt2(i, k)
5270 temp_ad = gama*exp(gama*log(-(rgas*temp)))*pe_ad(i, k)/(temp*&
5272 dm2_ad(i, k) = dm2_ad(i, k) + pt2(i, k)*temp_ad
5273 pt2_ad(i, k) = pt2_ad(i, k) + dm2(i, k)*temp_ad
5274 dz2_ad(i, k) = dz2_ad(i, k) - temp*temp_ad
5275 pm2_ad(i, k) = pm2_ad(i, k) - pe_ad(i, k)
5277 w2_ad(i, k) = w2_ad(i, k) + w1_ad(i, k)
5282 SUBROUTINE sim1_solver(dt, is, ie, km, rgas, gama, gm2, cp2, kappa, pe&
5283 & , dm2, pm2, pem, w2, dz2, pt2, ws, p_fac)
5285 INTEGER,
INTENT(IN) :: is, ie, km
5286 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac
5287 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pt2, pm2, gm2, cp2
5288 REAL,
INTENT(IN) :: ws(is:ie)
5289 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
5290 REAL,
INTENT(OUT) :: pe(is:ie, km+1)
5291 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
5293 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, g_rat, gam
5294 REAL,
DIMENSION(is:ie, km+1) :: pp
5295 REAL,
DIMENSION(is:ie) :: p1, bet
5296 REAL :: t1g, rdt, capa1
5309 pe(i, k) = exp(gama*log(-(dm2(i, k)/dz2(i, k)*rgas*pt2(i, k)))) &
5315 g_rat(i, k) = dm2(i, k)/dm2(i, k+1)
5316 bb(i, k) = 2.*(1.+g_rat(i, k))
5317 dd(i, k) = 3.*(pe(i, k)+g_rat(i, k)*pe(i, k+1))
5323 pp(i, 2) = dd(i, 1)/bet(i)
5325 dd(i, km) = 3.*pe(i, km)
5329 gam(i, k) = g_rat(i, k-1)/bet(i)
5330 bet(i) = bb(i, k) - gam(i, k)
5331 pp(i, k+1) = (dd(i, k)-pp(i, k))/bet(i)
5336 pp(i, k) = pp(i, k) - gam(i, k)*pp(i, k+1)
5342 aa(i, k) = t1g/(dz2(i, k-1)+dz2(i, k))*(pem(i, k)+pp(i, k))
5346 bet(i) = dm2(i, 1) - aa(i, 2)
5347 w2(i, 1) = (dm2(i, 1)*w1(i, 1)+dt*pp(i, 2))/bet(i)
5351 gam(i, k) = aa(i, k)/bet(i)
5352 bet(i) = dm2(i, k) - (aa(i, k)+aa(i, k+1)+aa(i, k)*gam(i, k))
5353 w2(i, k) = (dm2(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, k))-aa(i, k)&
5354 & *w2(i, k-1))/bet(i)
5358 p1(i) = t1g/dz2(i, km)*(pem(i, km+1)+pp(i, km+1))
5359 gam(i, km) = aa(i, km)/bet(i)
5360 bet(i) = dm2(i, km) - (aa(i, km)+p1(i)+aa(i, km)*gam(i, km))
5361 w2(i, km) = (dm2(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i, km))-p1(i)&
5362 & *ws(i)-aa(i, km)*w2(i, km-1))/bet(i)
5366 w2(i, k) = w2(i, k) - gam(i, k+1)*w2(i, k+1)
5374 pe(i, k+1) = pe(i, k) + dm2(i, k)*(w2(i, k)-w1(i, k))*rdt
5378 p1(i) = (pe(i, km)+2.*pe(i, km+1))*
r3 5379 IF (p_fac*pm2(i, km) .LT. p1(i) + pm2(i, km))
THEN 5380 max1 = p1(i) + pm2(i, km)
5382 max1 = p_fac*pm2(i, km)
5384 dz2(i, km) = -(dm2(i, km)*rgas*pt2(i, km)*exp(capa1*log(max1)))
5388 p1(i) = (pe(i, k)+bb(i, k)*pe(i, k+1)+g_rat(i, k)*pe(i, k+2))*
r3&
5389 & - g_rat(i, k)*p1(i)
5390 IF (p_fac*pm2(i, k) .LT. p1(i) + pm2(i, k))
THEN 5391 max2 = p1(i) + pm2(i, k)
5393 max2 = p_fac*pm2(i, k)
5395 dz2(i, k) = -(dm2(i, k)*rgas*pt2(i, k)*exp(capa1*log(max2)))
5419 SUBROUTINE sim_solver_fwd(dt, is, ie, km, rgas, gama, gm2, cp2, kappa&
5420 & , pe2, dm2, pm2, pem, w2, dz2, pt2, ws, alpha, p_fac, scale_m)
5422 INTEGER,
INTENT(IN) :: is, ie, km
5423 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac, alpha, scale_m
5424 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pt2, pm2, gm2, cp2
5425 REAL,
INTENT(IN) :: ws(is:ie)
5426 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
5427 REAL :: pe2(is:ie, km+1)
5428 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
5430 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, wk, g_rat, gam
5431 REAL,
DIMENSION(is:ie, km+1) :: pp
5432 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
5433 REAL :: beta, t2, t1g, rdt, ra, capa1
5461 t1g = 2.*gama*(alpha*dt)**2
5469 pe2(i, k) = exp(gama*log(-(dm2(i, k)/dz2(i, k)*rgas*pt2(i, k))))&
5475 g_rat(i, k) = dm2(i, k)/dm2(i, k+1)
5476 bb(i, k) = 2.*(1.+g_rat(i, k))
5477 dd(i, k) = 3.*(pe2(i, k)+g_rat(i, k)*pe2(i, k+1))
5483 pp(i, 2) = dd(i, 1)/bet(i)
5486 dd(i, km) = 3.*pe2(i, km)
5490 gam(i, k) = g_rat(i, k-1)/bet(i)
5492 bet(i) = bb(i, k) - gam(i, k)
5494 pp(i, k+1) = (dd(i, k)-pp(i, k))/bet(i)
5500 pp(i, k) = pp(i, k) - gam(i, k)*pp(i, k+1)
5507 pe2(i, k) = pem(i, k) + pp(i, k)
5512 aa(i, k) = t1g/(dz2(i, k-1)+dz2(i, k))*pe2(i, k)
5513 wk(i, k) = t2*aa(i, k)*(w1(i, k-1)-w1(i, k))
5515 aa(i, k) = aa(i, k) - scale_m*dm2(i, 1)
5521 bet(i) = dm2(i, 1) - aa(i, 2)
5523 w2(i, 1) = (dm2(i, 1)*w1(i, 1)+dt*pp(i, 2)+wk(i, 2))/bet(i)
5529 gam(i, k) = aa(i, k)/bet(i)
5531 bet(i) = dm2(i, k) - (aa(i, k)+aa(i, k+1)+aa(i, k)*gam(i, k))
5533 w2(i, k) = (dm2(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, k))+wk(i, k+&
5534 & 1)-wk(i, k)-aa(i, k)*w2(i, k-1))/bet(i)
5539 wk1(i) = t1g/dz2(i, km)*pe2(i, km+1)
5541 gam(i, km) = aa(i, km)/bet(i)
5543 bet(i) = dm2(i, km) - (aa(i, km)+wk1(i)+aa(i, km)*gam(i, km))
5545 w2(i, km) = (dm2(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i, km))-wk(i&
5546 & , km)+wk1(i)*(t2*w1(i, km)-ra*ws(i))-aa(i, km)*w2(i, km-1))/bet(&
5552 w2(i, k) = w2(i, k) - gam(i, k+1)*w2(i, k+1)
5562 pe2(i, k+1) = pe2(i, k) + (dm2(i, k)*(w2(i, k)-w1(i, k))*rdt-&
5563 & beta*(pp(i, k+1)-pp(i, k)))*ra
5567 p1(i) = (pe2(i, km)+2.*pe2(i, km+1))*
r3 5568 IF (p_fac*pm2(i, km) .LT. p1(i) + pm2(i, km))
THEN 5570 max1 = p1(i) + pm2(i, km)
5574 max1 = p_fac*pm2(i, km)
5578 dz2(i, km) = -(dm2(i, km)*rgas*pt2(i, km)*exp(capa1*log(max1)))
5583 p1(i) = (pe2(i, k)+bb(i, k)*pe2(i, k+1)+g_rat(i, k)*pe2(i, k+2))&
5584 & *
r3 - g_rat(i, k)*p1(i)
5585 IF (p_fac*pm2(i, k) .LT. p1(i) + pm2(i, k))
THEN 5587 max2 = p1(i) + pm2(i, k)
5591 max2 = p_fac*pm2(i, k)
5596 dz2(i, k) = -(dm2(i, k)*rgas*pt2(i, k)*exp(capa1*log(max2)))
5602 pe2(i, k) = pe2(i, k) + beta*(pp(i, k)-pe2(i, k))
5644 SUBROUTINE sim_solver_bwd(dt, is, ie, km, rgas, gama, gm2, cp2, kappa&
5645 & , pe2, pe2_ad, dm2, dm2_ad, pm2, pm2_ad, pem, pem_ad, w2, w2_ad, dz2&
5646 & , dz2_ad, pt2, pt2_ad, ws, ws_ad, alpha, p_fac, scale_m)
5648 INTEGER,
INTENT(IN) :: is, ie, km
5649 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac, alpha, scale_m
5650 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pt2, pm2, gm2, cp2
5651 REAL,
DIMENSION(is:ie, km) :: dm2_ad, pt2_ad, pm2_ad
5652 REAL,
INTENT(IN) :: ws(is:ie)
5653 REAL :: ws_ad(is:ie)
5654 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
5655 REAL,
DIMENSION(is:ie, km+1) :: pem_ad
5656 REAL :: pe2(is:ie, km+1)
5657 REAL :: pe2_ad(is:ie, km+1)
5658 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
5659 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2_ad, w2_ad
5660 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, wk, g_rat, gam
5661 REAL,
DIMENSION(is:ie, km) :: aa_ad, bb_ad, dd_ad, w1_ad, wk_ad, &
5663 REAL,
DIMENSION(is:ie, km+1) :: pp
5664 REAL,
DIMENSION(is:ie, km+1) :: pp_ad
5665 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
5666 REAL,
DIMENSION(is:ie) :: p1_ad, wk1_ad, bet_ad
5667 REAL :: beta, t2, t1g, rdt, ra, capa1
5745 pp_ad(i, k) = pp_ad(i, k) + beta*pe2_ad(i, k)
5746 pe2_ad(i, k) = (1.0-beta)*pe2_ad(i, k)
5755 temp5 = capa1*log(max2)
5756 temp_ad18 = -(rgas*exp(temp5)*dz2_ad(i, k))
5757 dm2_ad(i, k) = dm2_ad(i, k) + pt2(i, k)*temp_ad18
5758 pt2_ad(i, k) = pt2_ad(i, k) + dm2(i, k)*temp_ad18
5759 max2_ad = -(capa1*exp(temp5)*dm2(i, k)*pt2(i, k)*rgas*dz2_ad(i, &
5763 IF (branch .EQ. 0)
THEN 5765 p1_ad(i) = p1_ad(i) + max2_ad
5766 pm2_ad(i, k) = pm2_ad(i, k) + max2_ad
5769 pm2_ad(i, k) = pm2_ad(i, k) + p_fac*max2_ad
5772 temp_ad17 =
r3*p1_ad(i)
5773 pe2_ad(i, k) = pe2_ad(i, k) + temp_ad17
5774 bb_ad(i, k) = bb_ad(i, k) + pe2(i, k+1)*temp_ad17
5775 pe2_ad(i, k+1) = pe2_ad(i, k+1) + bb(i, k)*temp_ad17
5776 g_rat_ad(i, k) = g_rat_ad(i, k) + pe2(i, k+2)*temp_ad17 - p1(i)*&
5778 pe2_ad(i, k+2) = pe2_ad(i, k+2) + g_rat(i, k)*temp_ad17
5779 p1_ad(i) = -(g_rat(i, k)*p1_ad(i))
5784 temp4 = capa1*log(max1)
5785 temp_ad16 = -(rgas*exp(temp4)*dz2_ad(i, km))
5786 dm2_ad(i, km) = dm2_ad(i, km) + pt2(i, km)*temp_ad16
5787 pt2_ad(i, km) = pt2_ad(i, km) + dm2(i, km)*temp_ad16
5788 max1_ad = -(capa1*exp(temp4)*dm2(i, km)*pt2(i, km)*rgas*dz2_ad(i, &
5792 IF (branch .EQ. 0)
THEN 5794 p1_ad(i) = p1_ad(i) + max1_ad
5795 pm2_ad(i, km) = pm2_ad(i, km) + max1_ad
5798 pm2_ad(i, km) = pm2_ad(i, km) + p_fac*max1_ad
5800 pe2_ad(i, km) = pe2_ad(i, km) +
r3*p1_ad(i)
5801 pe2_ad(i, km+1) = pe2_ad(i, km+1) +
r3*2.*p1_ad(i)
5808 temp_ad14 = ra*pe2_ad(i, k+1)
5809 temp_ad15 = rdt*dm2(i, k)*temp_ad14
5810 pe2_ad(i, k) = pe2_ad(i, k) + pe2_ad(i, k+1)
5811 dm2_ad(i, k) = dm2_ad(i, k) + rdt*(w2(i, k)-w1(i, k))*temp_ad14
5812 w2_ad(i, k) = w2_ad(i, k) + temp_ad15
5813 w1_ad(i, k) = w1_ad(i, k) - temp_ad15
5814 pp_ad(i, k+1) = pp_ad(i, k+1) - beta*temp_ad14
5815 pp_ad(i, k) = pp_ad(i, k) + beta*temp_ad14
5816 pe2_ad(i, k+1) = 0.0
5827 gam_ad(i, k+1) = gam_ad(i, k+1) - w2(i, k+1)*w2_ad(i, k)
5828 w2_ad(i, k+1) = w2_ad(i, k+1) - gam(i, k+1)*w2_ad(i, k)
5837 temp_ad11 = w2_ad(i, km)/bet(i)
5838 temp3 = t2*w1(i, km) - ra*ws(i)
5839 w1_ad(i, km) = w1_ad(i, km) + (wk1(i)*t2+dm2(i, km))*temp_ad11
5840 pp_ad(i, km+1) = pp_ad(i, km+1) + dt*temp_ad11
5841 pp_ad(i, km) = pp_ad(i, km) - dt*temp_ad11
5842 wk_ad(i, km) = wk_ad(i, km) - temp_ad11
5843 ws_ad(i) = ws_ad(i) - wk1(i)*ra*temp_ad11
5844 w2_ad(i, km-1) = w2_ad(i, km-1) - aa(i, km)*temp_ad11
5845 bet_ad(i) = bet_ad(i) - (dm2(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i&
5846 & , km))-wk(i, km)+wk1(i)*temp3-aa(i, km)*w2(i, km-1))*temp_ad11/&
5848 dm2_ad(i, km) = dm2_ad(i, km) + bet_ad(i) + w1(i, km)*temp_ad11
5849 wk1_ad(i) = wk1_ad(i) + temp3*temp_ad11 - bet_ad(i)
5852 gam_ad(i, km) = gam_ad(i, km) - aa(i, km)*bet_ad(i)
5853 temp_ad12 = gam_ad(i, km)/bet(i)
5854 aa_ad(i, km) = aa_ad(i, km) + ((-1.0)-gam(i, km))*bet_ad(i) + &
5855 & temp_ad12 - w2(i, km-1)*temp_ad11
5856 bet_ad(i) = -(aa(i, km)*temp_ad12/bet(i))
5859 temp_ad13 = t1g*wk1_ad(i)/dz2(i, km)
5860 pe2_ad(i, km+1) = pe2_ad(i, km+1) + temp_ad13
5861 dz2_ad(i, km) = dz2_ad(i, km) - pe2(i, km+1)*temp_ad13/dz2(i, km)
5867 temp_ad9 = w2_ad(i, k)/bet(i)
5868 w1_ad(i, k) = w1_ad(i, k) + dm2(i, k)*temp_ad9
5869 pp_ad(i, k+1) = pp_ad(i, k+1) + dt*temp_ad9
5870 pp_ad(i, k) = pp_ad(i, k) - dt*temp_ad9
5871 wk_ad(i, k+1) = wk_ad(i, k+1) + temp_ad9
5872 wk_ad(i, k) = wk_ad(i, k) - temp_ad9
5873 w2_ad(i, k-1) = w2_ad(i, k-1) - aa(i, k)*temp_ad9
5874 bet_ad(i) = bet_ad(i) - (dm2(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i&
5875 & , k))+wk(i, k+1)-wk(i, k)-aa(i, k)*w2(i, k-1))*temp_ad9/bet(i)
5876 dm2_ad(i, k) = dm2_ad(i, k) + bet_ad(i) + w1(i, k)*temp_ad9
5877 aa_ad(i, k) = aa_ad(i, k) + ((-1.0)-gam(i, k))*bet_ad(i) - w2(i&
5881 aa_ad(i, k+1) = aa_ad(i, k+1) - bet_ad(i)
5882 gam_ad(i, k) = gam_ad(i, k) - aa(i, k)*bet_ad(i)
5884 temp_ad10 = gam_ad(i, k)/bet(i)
5885 bet_ad(i) = -(aa(i, k)*temp_ad10/bet(i))
5886 aa_ad(i, k) = aa_ad(i, k) + temp_ad10
5892 temp_ad8 = w2_ad(i, 1)/bet(i)
5893 w1_ad(i, 1) = w1_ad(i, 1) + dm2(i, 1)*temp_ad8
5894 pp_ad(i, 2) = pp_ad(i, 2) + dt*temp_ad8
5895 wk_ad(i, 2) = wk_ad(i, 2) + temp_ad8
5896 bet_ad(i) = bet_ad(i) - (dm2(i, 1)*w1(i, 1)+dt*pp(i, 2)+wk(i, 2))*&
5898 dm2_ad(i, 1) = dm2_ad(i, 1) + bet_ad(i) + w1(i, 1)*temp_ad8
5901 aa_ad(i, 2) = aa_ad(i, 2) - bet_ad(i)
5907 dm2_ad(i, 1) = dm2_ad(i, 1) - scale_m*aa_ad(i, k)
5908 temp_ad5 = t2*aa(i, k)*wk_ad(i, k)
5909 aa_ad(i, k) = aa_ad(i, k) + t2*(w1(i, k-1)-w1(i, k))*wk_ad(i, k)
5910 w1_ad(i, k-1) = w1_ad(i, k-1) + temp_ad5
5911 w1_ad(i, k) = w1_ad(i, k) - temp_ad5
5913 temp2 = dz2(i, k-1) + dz2(i, k)
5914 temp_ad6 = t1g*aa_ad(i, k)/temp2
5915 temp_ad7 = -(pe2(i, k)*temp_ad6/temp2)
5916 pe2_ad(i, k) = pe2_ad(i, k) + temp_ad6
5917 dz2_ad(i, k-1) = dz2_ad(i, k-1) + temp_ad7
5918 dz2_ad(i, k) = dz2_ad(i, k) + temp_ad7
5925 pem_ad(i, k) = pem_ad(i, k) + pe2_ad(i, k)
5926 pp_ad(i, k) = pp_ad(i, k) + pe2_ad(i, k)
5933 gam_ad(i, k) = gam_ad(i, k) - pp(i, k+1)*pp_ad(i, k)
5934 pp_ad(i, k+1) = pp_ad(i, k+1) - gam(i, k)*pp_ad(i, k)
5941 temp_ad3 = pp_ad(i, k+1)/bet(i)
5942 dd_ad(i, k) = dd_ad(i, k) + temp_ad3
5943 pp_ad(i, k) = pp_ad(i, k) - temp_ad3
5944 bet_ad(i) = bet_ad(i) - (dd(i, k)-pp(i, k))*temp_ad3/bet(i)
5947 bb_ad(i, k) = bb_ad(i, k) + bet_ad(i)
5948 gam_ad(i, k) = gam_ad(i, k) - bet_ad(i)
5949 temp_ad4 = gam_ad(i, k)/bet(i)
5950 bet_ad(i) = -(g_rat(i, k-1)*temp_ad4/bet(i))
5951 g_rat_ad(i, k-1) = g_rat_ad(i, k-1) + temp_ad4
5957 pe2_ad(i, km) = pe2_ad(i, km) + 3.*dd_ad(i, km)
5960 temp_ad2 = pp_ad(i, 2)/bet(i)
5961 dd_ad(i, 1) = dd_ad(i, 1) + temp_ad2
5962 bet_ad(i) = bet_ad(i) - dd(i, 1)*temp_ad2/bet(i)
5965 bb_ad(i, 1) = bb_ad(i, 1) + bet_ad(i)
5970 temp_ad0 = 3.*dd_ad(i, k)
5971 pe2_ad(i, k) = pe2_ad(i, k) + temp_ad0
5972 g_rat_ad(i, k) = g_rat_ad(i, k) + 2.*bb_ad(i, k) + pe2(i, k+1)*&
5974 pe2_ad(i, k+1) = pe2_ad(i, k+1) + g_rat(i, k)*temp_ad0
5977 temp_ad1 = g_rat_ad(i, k)/dm2(i, k+1)
5978 dm2_ad(i, k) = dm2_ad(i, k) + temp_ad1
5979 dm2_ad(i, k+1) = dm2_ad(i, k+1) - dm2(i, k)*temp_ad1/dm2(i, k+1)
5980 g_rat_ad(i, k) = 0.0
5987 temp0 = dm2(i, k)*pt2(i, k)
5989 temp_ad = gama*exp(gama*log(-(rgas*temp)))*pe2_ad(i, k)/(temp*&
5991 dm2_ad(i, k) = dm2_ad(i, k) + pt2(i, k)*temp_ad
5992 pt2_ad(i, k) = pt2_ad(i, k) + dm2(i, k)*temp_ad
5993 dz2_ad(i, k) = dz2_ad(i, k) - temp*temp_ad
5994 pm2_ad(i, k) = pm2_ad(i, k) - pe2_ad(i, k)
5996 w2_ad(i, k) = w2_ad(i, k) + w1_ad(i, k)
6001 SUBROUTINE sim_solver(dt, is, ie, km, rgas, gama, gm2, cp2, kappa, pe2&
6002 & , dm2, pm2, pem, w2, dz2, pt2, ws, alpha, p_fac, scale_m)
6004 INTEGER,
INTENT(IN) :: is, ie, km
6005 REAL,
INTENT(IN) :: dt, rgas, gama, kappa, p_fac, alpha, scale_m
6006 REAL,
DIMENSION(is:ie, km),
INTENT(IN) :: dm2, pt2, pm2, gm2, cp2
6007 REAL,
INTENT(IN) :: ws(is:ie)
6008 REAL,
DIMENSION(is:ie, km+1),
INTENT(IN) :: pem
6009 REAL,
INTENT(OUT) :: pe2(is:ie, km+1)
6010 REAL,
DIMENSION(is:ie, km),
INTENT(INOUT) :: dz2, w2
6012 REAL,
DIMENSION(is:ie, km) :: aa, bb, dd, w1, wk, g_rat, gam
6013 REAL,
DIMENSION(is:ie, km+1) :: pp
6014 REAL,
DIMENSION(is:ie) :: p1, wk1, bet
6015 REAL :: beta, t2, t1g, rdt, ra, capa1
6025 t1g = 2.*gama*(alpha*dt)**2
6032 pe2(i, k) = exp(gama*log(-(dm2(i, k)/dz2(i, k)*rgas*pt2(i, k))))&
6038 g_rat(i, k) = dm2(i, k)/dm2(i, k+1)
6039 bb(i, k) = 2.*(1.+g_rat(i, k))
6040 dd(i, k) = 3.*(pe2(i, k)+g_rat(i, k)*pe2(i, k+1))
6046 pp(i, 2) = dd(i, 1)/bet(i)
6048 dd(i, km) = 3.*pe2(i, km)
6052 gam(i, k) = g_rat(i, k-1)/bet(i)
6053 bet(i) = bb(i, k) - gam(i, k)
6054 pp(i, k+1) = (dd(i, k)-pp(i, k))/bet(i)
6059 pp(i, k) = pp(i, k) - gam(i, k)*pp(i, k+1)
6065 pe2(i, k) = pem(i, k) + pp(i, k)
6070 aa(i, k) = t1g/(dz2(i, k-1)+dz2(i, k))*pe2(i, k)
6071 wk(i, k) = t2*aa(i, k)*(w1(i, k-1)-w1(i, k))
6072 aa(i, k) = aa(i, k) - scale_m*dm2(i, 1)
6077 bet(i) = dm2(i, 1) - aa(i, 2)
6078 w2(i, 1) = (dm2(i, 1)*w1(i, 1)+dt*pp(i, 2)+wk(i, 2))/bet(i)
6083 gam(i, k) = aa(i, k)/bet(i)
6084 bet(i) = dm2(i, k) - (aa(i, k)+aa(i, k+1)+aa(i, k)*gam(i, k))
6085 w2(i, k) = (dm2(i, k)*w1(i, k)+dt*(pp(i, k+1)-pp(i, k))+wk(i, k+&
6086 & 1)-wk(i, k)-aa(i, k)*w2(i, k-1))/bet(i)
6091 wk1(i) = t1g/dz2(i, km)*pe2(i, km+1)
6092 gam(i, km) = aa(i, km)/bet(i)
6093 bet(i) = dm2(i, km) - (aa(i, km)+wk1(i)+aa(i, km)*gam(i, km))
6094 w2(i, km) = (dm2(i, km)*w1(i, km)+dt*(pp(i, km+1)-pp(i, km))-wk(i&
6095 & , km)+wk1(i)*(t2*w1(i, km)-ra*ws(i))-aa(i, km)*w2(i, km-1))/bet(&
6100 w2(i, k) = w2(i, k) - gam(i, k+1)*w2(i, k+1)
6108 pe2(i, k+1) = pe2(i, k) + (dm2(i, k)*(w2(i, k)-w1(i, k))*rdt-&
6109 & beta*(pp(i, k+1)-pp(i, k)))*ra
6113 p1(i) = (pe2(i, km)+2.*pe2(i, km+1))*
r3 6114 IF (p_fac*pm2(i, km) .LT. p1(i) + pm2(i, km))
THEN 6115 max1 = p1(i) + pm2(i, km)
6117 max1 = p_fac*pm2(i, km)
6119 dz2(i, km) = -(dm2(i, km)*rgas*pt2(i, km)*exp(capa1*log(max1)))
6123 p1(i) = (pe2(i, k)+bb(i, k)*pe2(i, k+1)+g_rat(i, k)*pe2(i, k+2))&
6124 & *
r3 - g_rat(i, k)*p1(i)
6125 IF (p_fac*pm2(i, k) .LT. p1(i) + pm2(i, k))
THEN 6126 max2 = p1(i) + pm2(i, k)
6128 max2 = p_fac*pm2(i, k)
6131 dz2(i, k) = -(dm2(i, k)*rgas*pt2(i, k)*exp(capa1*log(max2)))
6136 pe2(i, k) = pe2(i, k) + beta*(pp(i, k)-pe2(i, k))
6143 INTEGER,
INTENT(IN) :: i1, i2, km
6145 INTEGER,
INTENT(IN) :: id
6146 REAL,
DIMENSION(i1:i2, km),
INTENT(IN) :: q1
6147 REAL,
DIMENSION(i1:i2, km+1),
INTENT(OUT) :: qe
6149 REAL,
PARAMETER :: r2o3=2./3.
6150 REAL,
PARAMETER :: r4o3=4./3.
6159 qe(i, 1) = r4o3*q1(i, 1) + r2o3*q1(i, 2)
6168 gak(k) = 1./(4.-gak(k-1))
6170 qe(i, k) = (3.*(q1(i, k-1)+q1(i, k))-qe(i, k-1))*gak(k)
6173 bet = 1./(1.5-3.5*gak(km))
6175 qe(i, km+1) = (4.*q1(i, km)+q1(i, km-1)-3.5*qe(i, km))*bet
6179 qe(i, k) = qe(i, k) - gak(k)*qe(i, k+1)
6203 SUBROUTINE edge_profile_fwd(q1, q2, q1e, q2e, i1, i2, j1, j2, j, km, &
6204 & dp0, uniform_grid, limiter)
6207 INTEGER,
INTENT(IN) :: i1, i2, j1, j2
6208 INTEGER,
INTENT(IN) :: j, km
6209 INTEGER,
INTENT(IN) :: limiter
6210 LOGICAL,
INTENT(IN) :: uniform_grid
6211 REAL,
INTENT(IN) :: dp0(km)
6212 REAL,
DIMENSION(i1:i2, j1:j2, km),
INTENT(IN) :: q1, q2
6213 REAL,
DIMENSION(i1:i2, j1:j2, km+1) :: q1e, q2e
6216 REAL,
DIMENSION(i1:i2, km+1) :: qe1, qe2, gam
6218 REAL :: bet, r2o3, r4o3
6219 REAL :: g0, gk, xt1, xt2, a_bot
6235 IF (uniform_grid)
THEN 6242 qe1(i, 1) = r4o3*q1(i, j, 1) + r2o3*q1(i, j, 2)
6243 qe2(i, 1) = r4o3*q2(i, j, 1) + r2o3*q2(i, j, 2)
6248 gak(k) = 1./(4.-gak(k-1))
6250 qe1(i, k) = (3.*(q1(i, j, k-1)+q1(i, j, k))-qe1(i, k-1))*gak(k&
6252 qe2(i, k) = (3.*(q2(i, j, k-1)+q2(i, j, k))-qe2(i, k-1))*gak(k&
6256 bet = 1./(1.5-3.5*gak(km))
6258 qe1(i, km+1) = (4.*q1(i, j, km)+q1(i, j, km-1)-3.5*qe1(i, km))*&
6260 qe2(i, km+1) = (4.*q2(i, j, km)+q2(i, j, km-1)-3.5*qe2(i, km))*&
6265 qe1(i, k) = qe1(i, k) - gak(k)*qe1(i, k+1)
6266 qe2(i, k) = qe2(i, k) - gak(k)*qe2(i, k+1)
6276 qe1(i, 1) = (xt1*q1(i, j, 1)+q1(i, j, 2))/bet
6277 qe2(i, 1) = (xt1*q2(i, j, 1)+q2(i, j, 2))/bet
6278 gam(i, 1) = (1.+g0*(g0+1.5))/bet
6282 gk = dp0(k-1)/dp0(k)
6285 bet = 2. + 2.*gk - gam(i, k-1)
6286 qe1(i, k) = (3.*(q1(i, j, k-1)+gk*q1(i, j, k))-qe1(i, k-1))/&
6288 qe2(i, k) = (3.*(q2(i, j, k-1)+gk*q2(i, j, k))-qe2(i, k-1))/&
6293 a_bot = 1. + gk*(gk+1.5)
6298 xt2 = gk*(gk+0.5) - a_bot*gam(i, km)
6299 qe1(i, km+1) = (xt1*q1(i, j, km)+q1(i, j, km-1)-a_bot*qe1(i, km)&
6301 qe2(i, km+1) = (xt1*q2(i, j, km)+q2(i, j, km-1)-a_bot*qe2(i, km)&
6306 qe1(i, k) = qe1(i, k) - gam(i, k)*qe1(i, k+1)
6307 qe2(i, k) = qe2(i, k) - gam(i, k)*qe2(i, k+1)
6315 IF (limiter .NE. 0)
THEN 6319 IF (q1(i, j, 1)*qe1(i, 1) .LT. 0.)
THEN 6325 IF (q2(i, j, 1)*qe2(i, 1) .LT. 0.)
THEN 6332 IF (q1(i, j, km)*qe1(i, km+1) .LT. 0.)
THEN 6338 IF (q2(i, j, km)*qe2(i, km+1) .LT. 0.)
THEN 6351 q1e(i, j, k) = qe1(i, k)
6352 q2e(i, j, k) = qe2(i, k)
6384 & q2e_ad, i1, i2, j1, j2, j, km, dp0, uniform_grid, limiter)
6386 INTEGER,
INTENT(IN) :: i1, i2, j1, j2
6387 INTEGER,
INTENT(IN) :: j, km
6388 INTEGER,
INTENT(IN) :: limiter
6389 LOGICAL,
INTENT(IN) :: uniform_grid
6390 REAL,
INTENT(IN) :: dp0(km)
6391 REAL,
DIMENSION(i1:i2, j1:j2, km),
INTENT(IN) :: q1, q2
6392 REAL,
DIMENSION(i1:i2, j1:j2, km) :: q1_ad, q2_ad
6393 REAL,
DIMENSION(i1:i2, j1:j2, km+1) :: q1e, q2e
6394 REAL,
DIMENSION(i1:i2, j1:j2, km+1) :: q1e_ad, q2e_ad
6395 REAL,
DIMENSION(i1:i2, km+1) :: qe1, qe2, gam
6396 REAL,
DIMENSION(i1:i2, km+1) :: qe1_ad, qe2_ad
6398 REAL :: bet, r2o3, r4o3
6399 REAL :: g0, gk, xt1, xt2, a_bot
6436 qe2_ad(i, k) = qe2_ad(i, k) + q2e_ad(i, j, k)
6437 q2e_ad(i, j, k) = 0.0
6438 qe1_ad(i, k) = qe1_ad(i, k) + q1e_ad(i, j, k)
6439 q1e_ad(i, j, k) = 0.0
6443 IF (branch .NE. 0)
THEN 6446 IF (branch .NE. 0) qe2_ad(i, km+1) = 0.0
6448 IF (branch .EQ. 0) qe1_ad(i, km+1) = 0.0
6450 IF (branch .EQ. 0) qe2_ad(i, 1) = 0.0
6452 IF (branch .EQ. 0) qe1_ad(i, 1) = 0.0
6456 IF (branch .EQ. 0)
THEN 6459 qe2_ad(i, k+1) = qe2_ad(i, k+1) - gak(k)*qe2_ad(i, k)
6460 qe1_ad(i, k+1) = qe1_ad(i, k+1) - gak(k)*qe1_ad(i, k)
6464 temp_ad1 = bet*qe2_ad(i, km+1)
6465 q2_ad(i, j, km) = q2_ad(i, j, km) + 4.*temp_ad1
6466 q2_ad(i, j, km-1) = q2_ad(i, j, km-1) + temp_ad1
6467 qe2_ad(i, km) = qe2_ad(i, km) - 3.5*temp_ad1
6468 qe2_ad(i, km+1) = 0.0
6469 temp_ad2 = bet*qe1_ad(i, km+1)
6470 q1_ad(i, j, km) = q1_ad(i, j, km) + 4.*temp_ad2
6471 q1_ad(i, j, km-1) = q1_ad(i, j, km-1) + temp_ad2
6472 qe1_ad(i, km) = qe1_ad(i, km) - 3.5*temp_ad2
6473 qe1_ad(i, km+1) = 0.0
6477 temp_ad = gak(k)*qe2_ad(i, k)
6478 q2_ad(i, j, k-1) = q2_ad(i, j, k-1) + 3.*temp_ad
6479 q2_ad(i, j, k) = q2_ad(i, j, k) + 3.*temp_ad
6480 qe2_ad(i, k-1) = qe2_ad(i, k-1) - temp_ad
6482 temp_ad0 = gak(k)*qe1_ad(i, k)
6483 q1_ad(i, j, k-1) = q1_ad(i, j, k-1) + 3.*temp_ad0
6484 q1_ad(i, j, k) = q1_ad(i, j, k) + 3.*temp_ad0
6485 qe1_ad(i, k-1) = qe1_ad(i, k-1) - temp_ad0
6493 q2_ad(i, j, 1) = q2_ad(i, j, 1) + r4o3*qe2_ad(i, 1)
6494 q2_ad(i, j, 2) = q2_ad(i, j, 2) + r2o3*qe2_ad(i, 1)
6496 q1_ad(i, j, 1) = q1_ad(i, j, 1) + r4o3*qe1_ad(i, 1)
6497 q1_ad(i, j, 2) = q1_ad(i, j, 2) + r2o3*qe1_ad(i, 1)
6503 qe2_ad(i, k+1) = qe2_ad(i, k+1) - gam(i, k)*qe2_ad(i, k)
6504 qe1_ad(i, k+1) = qe1_ad(i, k+1) - gam(i, k)*qe1_ad(i, k)
6508 temp_ad5 = qe2_ad(i, km+1)/xt2
6509 q2_ad(i, j, km) = q2_ad(i, j, km) + xt1*temp_ad5
6510 q2_ad(i, j, km-1) = q2_ad(i, j, km-1) + temp_ad5
6511 qe2_ad(i, km) = qe2_ad(i, km) - a_bot*temp_ad5
6512 qe2_ad(i, km+1) = 0.0
6513 temp_ad6 = qe1_ad(i, km+1)/xt2
6514 q1_ad(i, j, km) = q1_ad(i, j, km) + xt1*temp_ad6
6515 q1_ad(i, j, km-1) = q1_ad(i, j, km-1) + temp_ad6
6516 qe1_ad(i, km) = qe1_ad(i, km) - a_bot*temp_ad6
6517 qe1_ad(i, km+1) = 0.0
6523 temp_ad3 = qe2_ad(i, k)/bet
6524 q2_ad(i, j, k-1) = q2_ad(i, j, k-1) + 3.*temp_ad3
6525 q2_ad(i, j, k) = q2_ad(i, j, k) + 3.*gk*temp_ad3
6526 qe2_ad(i, k-1) = qe2_ad(i, k-1) - temp_ad3
6528 temp_ad4 = qe1_ad(i, k)/bet
6529 q1_ad(i, j, k-1) = q1_ad(i, j, k-1) + 3.*temp_ad4
6530 q1_ad(i, j, k) = q1_ad(i, j, k) + 3.*gk*temp_ad4
6531 qe1_ad(i, k-1) = qe1_ad(i, k-1) - temp_ad4
6538 q2_ad(i, j, 1) = q2_ad(i, j, 1) + xt1*qe2_ad(i, 1)/bet
6539 q2_ad(i, j, 2) = q2_ad(i, j, 2) + qe2_ad(i, 1)/bet
6541 q1_ad(i, j, 1) = q1_ad(i, j, 1) + xt1*qe1_ad(i, 1)/bet
6542 q1_ad(i, j, 2) = q1_ad(i, j, 2) + qe1_ad(i, 1)/bet
6547 SUBROUTINE edge_profile(q1, q2, q1e, q2e, i1, i2, j1, j2, j, km, dp0, &
6548 & uniform_grid, limiter)
6551 INTEGER,
INTENT(IN) :: i1, i2, j1, j2
6552 INTEGER,
INTENT(IN) :: j, km
6553 INTEGER,
INTENT(IN) :: limiter
6554 LOGICAL,
INTENT(IN) :: uniform_grid
6555 REAL,
INTENT(IN) :: dp0(km)
6556 REAL,
DIMENSION(i1:i2, j1:j2, km),
INTENT(IN) :: q1, q2
6557 REAL,
DIMENSION(i1:i2, j1:j2, km+1),
INTENT(OUT) :: q1e, q2e
6560 REAL,
DIMENSION(i1:i2, km+1) :: qe1, qe2, gam
6562 REAL :: bet, r2o3, r4o3
6563 REAL :: g0, gk, xt1, xt2, a_bot
6565 IF (uniform_grid)
THEN 6572 qe1(i, 1) = r4o3*q1(i, j, 1) + r2o3*q1(i, j, 2)
6573 qe2(i, 1) = r4o3*q2(i, j, 1) + r2o3*q2(i, j, 2)
6577 gak(k) = 1./(4.-gak(k-1))
6579 qe1(i, k) = (3.*(q1(i, j, k-1)+q1(i, j, k))-qe1(i, k-1))*gak(k&
6581 qe2(i, k) = (3.*(q2(i, j, k-1)+q2(i, j, k))-qe2(i, k-1))*gak(k&
6585 bet = 1./(1.5-3.5*gak(km))
6587 qe1(i, km+1) = (4.*q1(i, j, km)+q1(i, j, km-1)-3.5*qe1(i, km))*&
6589 qe2(i, km+1) = (4.*q2(i, j, km)+q2(i, j, km-1)-3.5*qe2(i, km))*&
6594 qe1(i, k) = qe1(i, k) - gak(k)*qe1(i, k+1)
6595 qe2(i, k) = qe2(i, k) - gak(k)*qe2(i, k+1)
6604 qe1(i, 1) = (xt1*q1(i, j, 1)+q1(i, j, 2))/bet
6605 qe2(i, 1) = (xt1*q2(i, j, 1)+q2(i, j, 2))/bet
6606 gam(i, 1) = (1.+g0*(g0+1.5))/bet
6609 gk = dp0(k-1)/dp0(k)
6611 bet = 2. + 2.*gk - gam(i, k-1)
6612 qe1(i, k) = (3.*(q1(i, j, k-1)+gk*q1(i, j, k))-qe1(i, k-1))/&
6614 qe2(i, k) = (3.*(q2(i, j, k-1)+gk*q2(i, j, k))-qe2(i, k-1))/&
6619 a_bot = 1. + gk*(gk+1.5)
6622 xt2 = gk*(gk+0.5) - a_bot*gam(i, km)
6623 qe1(i, km+1) = (xt1*q1(i, j, km)+q1(i, j, km-1)-a_bot*qe1(i, km)&
6625 qe2(i, km+1) = (xt1*q2(i, j, km)+q2(i, j, km-1)-a_bot*qe2(i, km)&
6630 qe1(i, k) = qe1(i, k) - gam(i, k)*qe1(i, k+1)
6631 qe2(i, k) = qe2(i, k) - gam(i, k)*qe2(i, k+1)
6638 IF (limiter .NE. 0)
THEN 6642 IF (q1(i, j, 1)*qe1(i, 1) .LT. 0.) qe1(i, 1) = 0.
6643 IF (q2(i, j, 1)*qe2(i, 1) .LT. 0.) qe2(i, 1) = 0.
6645 IF (q1(i, j, km)*qe1(i, km+1) .LT. 0.) qe1(i, km+1) = 0.
6646 IF (q2(i, j, km)*qe2(i, km+1) .LT. 0.) qe2(i, km+1) = 0.
6651 q1e(i, j, k) = qe1(i, k)
6652 q2e(i, j, k) = qe2(i, k)
6677 & phis, pkc, gz, pk3, npx, npy, npz, nested, pkc_pertn, computepk3, &
6682 INTEGER,
INTENT(IN) :: npx, npy, npz
6683 LOGICAL,
INTENT(IN) :: pkc_pertn, computepk3, fullhalo, nested
6684 REAL,
INTENT(IN) :: ptop, kappa, cp,
grav 6686 REAL,
INTENT(IN) :: phis(bd%isd:bd%ied, bd%jsd:bd%jed)
6687 REAL,
DIMENSION(bd%isd:bd%ied, bd%jsd:bd%jed, npz),
INTENT(IN) :: pt&
6689 REAL,
DIMENSION(bd%isd:bd%ied, bd%jsd:bd%jed, npz+1),
INTENT(INOUT) &
6694 REAL :: ptk, rgrav, rkap, peln1, rdg
6695 REAL,
DIMENSION(bd%isd:bd%ied, npz+1, bd%jsd:bd%jed) :: pe, peln
6696 REAL,
DIMENSION(bd%isd:bd%ied, npz) :: gam, bb, dd, pkz
6697 REAL,
DIMENSION(bd%isd:bd%ied, npz-1) :: g_rat
6698 REAL,
DIMENSION(bd%isd:bd%ied) :: bet
6700 INTEGER :: ifirst, ilast, jfirst, jlast
6701 INTEGER :: is, ie, js, je
6702 INTEGER :: isd, ied, jsd, jed
6741 IF (.NOT.nested)
THEN 6751 gama = 1./(1.-kappa)
6760 gz(i, j, npz+1) = phis(i, j)
6765 gz(i, j, k) = gz(i, j, k+1) - delz(i, j, k)*
grav 6771 peln(i, 1, j) = peln1
6776 pe(i, k, j) = pe(i, k-1, j) + delp(i, j, k-1)
6777 peln(i, k, j) = log(pe(i, k, j))
6785 pkz(i, k) = exp(gama*log(-(delp(i, j, k)*rgrav/delz(i, j, &
6786 & k)*
rdgas*pt(i, j, k))))
6788 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
6790 pkz(i, k) = pkz(i, k) - pm
6797 g_rat(i, k) = delp(i, j, k)/delp(i, j, k+1)
6798 bb(i, k) = 2.*(1.+g_rat(i, k))
6800 dd(i, k) = 3.*(pkz(i, k)+g_rat(i, k)*pkz(i, k+1))
6809 pkc(i, j, 2) = dd(i, 1)/bet(i)
6812 dd(i, npz) = 3.*pkz(i, npz)
6817 gam(i, k) = g_rat(i, k-1)/bet(i)
6819 bet(i) = bb(i, k) - gam(i, k)
6821 pkc(i, j, k+1) = (dd(i, k)-pkc(i, j, k))/bet(i)
6827 pkc(i, j, k) = pkc(i, j, k) - gam(i, k)*pkc(i, j, k+1)
6832 IF (.NOT.pkc_pertn)
THEN 6836 pkc(i, j, k) = pkc(i, j, k) + pe(i, k, j)
6844 IF (computepk3)
THEN 6852 pk3(i, j, k) = exp(kappa*log(pe(i, k, j)))
6864 IF (ie .EQ. npx - 1)
THEN 6869 gz(i, j, npz+1) = phis(i, j)
6874 gz(i, j, k) = gz(i, j, k+1) - delz(i, j, k)*
grav 6882 peln(i, 1, j) = peln1
6887 pe(i, k, j) = pe(i, k-1, j) + delp(i, j, k-1)
6889 peln(i, k, j) = log(pe(i, k, j))
6897 pkz(i, k) = exp(gama*log(-(delp(i, j, k)*rgrav/delz(i, j, &
6898 & k)*
rdgas*pt(i, j, k))))
6900 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
6902 pkz(i, k) = pkz(i, k) - pm
6909 g_rat(i, k) = delp(i, j, k)/delp(i, j, k+1)
6910 bb(i, k) = 2.*(1.+g_rat(i, k))
6912 dd(i, k) = 3.*(pkz(i, k)+g_rat(i, k)*pkz(i, k+1))
6921 pkc(i, j, 2) = dd(i, 1)/bet(i)
6924 dd(i, npz) = 3.*pkz(i, npz)
6929 gam(i, k) = g_rat(i, k-1)/bet(i)
6931 bet(i) = bb(i, k) - gam(i, k)
6933 pkc(i, j, k+1) = (dd(i, k)-pkc(i, j, k))/bet(i)
6939 pkc(i, j, k) = pkc(i, j, k) - gam(i, k)*pkc(i, j, k+1)
6944 IF (.NOT.pkc_pertn)
THEN 6948 pkc(i, j, k) = pkc(i, j, k) + pe(i, k, j)
6956 IF (computepk3)
THEN 6964 pk3(i, j, k) = exp(kappa*log(pe(i, k, j)))
6981 gz(i, j, npz+1) = phis(i, j)
6986 gz(i, j, k) = gz(i, j, k+1) - delz(i, j, k)*
grav 6994 peln(i, 1, j) = peln1
6999 pe(i, k, j) = pe(i, k-1, j) + delp(i, j, k-1)
7001 peln(i, k, j) = log(pe(i, k, j))
7009 pkz(i, k) = exp(gama*log(-(delp(i, j, k)*rgrav/delz(i, j, &
7010 & k)*
rdgas*pt(i, j, k))))
7013 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
7015 pkz(i, k) = pkz(i, k) - pm
7022 g_rat(i, k) = delp(i, j, k)/delp(i, j, k+1)
7023 bb(i, k) = 2.*(1.+g_rat(i, k))
7025 dd(i, k) = 3.*(pkz(i, k)+g_rat(i, k)*pkz(i, k+1))
7034 pkc(i, j, 2) = dd(i, 1)/bet(i)
7037 dd(i, npz) = 3.*pkz(i, npz)
7042 gam(i, k) = g_rat(i, k-1)/bet(i)
7044 bet(i) = bb(i, k) - gam(i, k)
7046 pkc(i, j, k+1) = (dd(i, k)-pkc(i, j, k))/bet(i)
7052 pkc(i, j, k) = pkc(i, j, k) - gam(i, k)*pkc(i, j, k+1)
7057 IF (.NOT.pkc_pertn)
THEN 7061 pkc(i, j, k) = pkc(i, j, k) + pe(i, k, j)
7069 IF (computepk3)
THEN 7077 pk3(i, j, k) = exp(kappa*log(pe(i, k, j)))
7089 IF (je .EQ. npy - 1)
THEN 7094 gz(i, j, npz+1) = phis(i, j)
7099 gz(i, j, k) = gz(i, j, k+1) - delz(i, j, k)*
grav 7107 peln(i, 1, j) = peln1
7112 pe(i, k, j) = pe(i, k-1, j) + delp(i, j, k-1)
7114 peln(i, k, j) = log(pe(i, k, j))
7122 pkz(i, k) = exp(gama*log(-(delp(i, j, k)*rgrav/delz(i, j, &
7123 & k)*
rdgas*pt(i, j, k))))
7126 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
7128 pkz(i, k) = pkz(i, k) - pm
7136 g_rat(i, k) = delp(i, j, k)/delp(i, j, k+1)
7137 bb(i, k) = 2.*(1.+g_rat(i, k))
7139 dd(i, k) = 3.*(pkz(i, k)+g_rat(i, k)*pkz(i, k+1))
7148 pkc(i, j, 2) = dd(i, 1)/bet(i)
7151 dd(i, npz) = 3.*pkz(i, npz)
7156 gam(i, k) = g_rat(i, k-1)/bet(i)
7158 bet(i) = bb(i, k) - gam(i, k)
7160 pkc(i, j, k+1) = (dd(i, k)-pkc(i, j, k))/bet(i)
7166 pkc(i, j, k) = pkc(i, j, k) - gam(i, k)*pkc(i, j, k+1)
7171 IF (.NOT.pkc_pertn)
THEN 7175 pkc(i, j, k) = pkc(i, j, k) + pe(i, k, j)
7183 IF (computepk3)
THEN 7191 pk3(i, j, k) = exp(kappa*log(pe(i, k, j)))
7202 CALL pushrealarray(pe, (bd%ied-bd%isd+1)*(npz+1)*(bd%jed-bd%jsd&
7210 CALL pushrealarray(peln, (bd%ied-bd%isd+1)*(npz+1)*(bd%jed-bd%&
7218 CALL pushrealarray(pe, (bd%ied-bd%isd+1)*(npz+1)*(bd%jed-bd%jsd&
7226 CALL pushrealarray(peln, (bd%ied-bd%isd+1)*(npz+1)*(bd%jed-bd%&
7254 & , delz_ad, pt, pt_ad, phis, pkc, pkc_ad, gz, gz_ad, pk3, pk3_ad, npx&
7255 & , npy, npz, nested, pkc_pertn, computepk3, fullhalo, bd)
7257 INTEGER,
INTENT(IN) :: npx, npy, npz
7258 LOGICAL,
INTENT(IN) :: pkc_pertn, computepk3, fullhalo, nested
7259 REAL,
INTENT(IN) :: ptop, kappa, cp,
grav 7261 REAL,
INTENT(IN) :: phis(bd%isd:bd%ied, bd%jsd:bd%jed)
7262 REAL,
DIMENSION(bd%isd:bd%ied, bd%jsd:bd%jed, npz),
INTENT(IN) :: pt&
7264 REAL,
DIMENSION(bd%isd:bd%ied, bd%jsd:bd%jed, npz) :: pt_ad, delp_ad&
7266 REAL,
DIMENSION(bd%isd:bd%ied, bd%jsd:bd%jed, npz+1),
INTENT(INOUT) &
7268 REAL,
DIMENSION(bd%isd:bd%ied, bd%jsd:bd%jed, npz+1),
INTENT(INOUT) &
7269 & :: gz_ad, pkc_ad, pk3_ad
7272 REAL :: ptk, rgrav, rkap, peln1, rdg
7273 REAL,
DIMENSION(bd%isd:bd%ied, npz+1, bd%jsd:bd%jed) :: pe, peln
7274 REAL,
DIMENSION(bd%isd:bd%ied, npz+1, bd%jsd:bd%jed) :: pe_ad, &
7276 REAL,
DIMENSION(bd%isd:bd%ied, npz) :: gam, bb, dd, pkz
7277 REAL,
DIMENSION(bd%isd:bd%ied, npz) :: gam_ad, bb_ad, dd_ad, pkz_ad
7278 REAL,
DIMENSION(bd%isd:bd%ied, npz-1) :: g_rat
7279 REAL,
DIMENSION(bd%isd:bd%ied, npz-1) :: g_rat_ad
7280 REAL,
DIMENSION(bd%isd:bd%ied) :: bet
7281 REAL,
DIMENSION(bd%isd:bd%ied) :: bet_ad
7284 INTEGER :: ifirst, ilast, jfirst, jlast
7285 INTEGER :: is, ie, js, je
7286 INTEGER :: isd, ied, jsd, jed
7368 IF (branch .NE. 0)
THEN 7369 IF (branch .EQ. 1)
THEN 7371 CALL poprealarray(peln, (bd%ied-bd%isd+1)*(npz+1)*(bd%jed-bd%&
7379 CALL poprealarray(pe, (bd%ied-bd%isd+1)*(npz+1)*(bd%jed-bd%jsd+&
7394 CALL poprealarray(peln, (bd%ied-bd%isd+1)*(npz+1)*(bd%jed-bd%&
7402 CALL poprealarray(pe, (bd%ied-bd%isd+1)*(npz+1)*(bd%jed-bd%jsd+&
7410 IF (branch .NE. 0)
THEN 7412 DO i=ilast,ifirst,-1
7414 pe_ad(i, k, j) = pe_ad(i, k, j) + kappa*exp(kappa*log(pe&
7415 & (i, k, j)))*pk3_ad(i, j, k)/pe(i, k, j)
7416 pk3_ad(i, j, k) = 0.0
7419 DO i=ilast,ifirst,-1
7421 pk3_ad(i, j, 1) = 0.0
7425 IF (branch .EQ. 0)
THEN 7427 DO i=ilast,ifirst,-1
7429 pe_ad(i, k, j) = pe_ad(i, k, j) + pkc_ad(i, j, k)
7443 DO i=ilast,ifirst,-1
7445 gam_ad(i, k) = gam_ad(i, k) - pkc(i, j, k+1)*pkc_ad(i, j, &
7447 pkc_ad(i, j, k+1) = pkc_ad(i, j, k+1) - gam(i, k)*pkc_ad(i&
7452 DO i=ilast,ifirst,-1
7454 temp_ad25 = pkc_ad(i, j, k+1)/bet(i)
7455 dd_ad(i, k) = dd_ad(i, k) + temp_ad25
7456 pkc_ad(i, j, k) = pkc_ad(i, j, k) - temp_ad25
7457 bet_ad(i) = bet_ad(i) - (dd(i, k)-pkc(i, j, k))*temp_ad25/&
7459 pkc_ad(i, j, k+1) = 0.0
7461 bb_ad(i, k) = bb_ad(i, k) + bet_ad(i)
7462 gam_ad(i, k) = gam_ad(i, k) - bet_ad(i)
7464 temp_ad26 = gam_ad(i, k)/bet(i)
7465 bet_ad(i) = -(g_rat(i, k-1)*temp_ad26/bet(i))
7466 g_rat_ad(i, k-1) = g_rat_ad(i, k-1) + temp_ad26
7470 DO i=ilast,ifirst,-1
7472 pkz_ad(i, npz) = pkz_ad(i, npz) + 3.*dd_ad(i, npz)
7476 temp_ad24 = pkc_ad(i, j, 2)/bet(i)
7477 dd_ad(i, 1) = dd_ad(i, 1) + temp_ad24
7478 bet_ad(i) = bet_ad(i) - dd(i, 1)*temp_ad24/bet(i)
7479 pkc_ad(i, j, 2) = 0.0
7481 pkc_ad(i, j, 1) = 0.0
7483 bb_ad(i, 1) = bb_ad(i, 1) + bet_ad(i)
7487 DO i=ilast,ifirst,-1
7489 temp_ad22 = 3.*dd_ad(i, k)
7490 pkz_ad(i, k) = pkz_ad(i, k) + temp_ad22
7491 g_rat_ad(i, k) = g_rat_ad(i, k) + 2.*bb_ad(i, k) + pkz(i, &
7493 pkz_ad(i, k+1) = pkz_ad(i, k+1) + g_rat(i, k)*temp_ad22
7497 temp_ad23 = g_rat_ad(i, k)/delp(i, j, k+1)
7498 delp_ad(i, j, k) = delp_ad(i, j, k) + temp_ad23
7499 delp_ad(i, j, k+1) = delp_ad(i, j, k+1) - delp(i, j, k)*&
7500 & temp_ad23/delp(i, j, k+1)
7501 g_rat_ad(i, k) = 0.0
7505 DO i=ilast,ifirst,-1
7506 temp17 = delz(i, j, k)
7507 temp16 = delp(i, j, k)*pt(i, j, k)
7508 temp14 = temp16/temp17
7509 temp15 = -(rgrav*
rdgas*temp14)
7510 temp_ad21 = -(rgrav*
rdgas*gama*exp(gama*log(temp15))*&
7511 & pkz_ad(i, k)/(temp15*temp17))
7512 pm_ad = -pkz_ad(i, k)
7513 temp18 = peln(i, k+1, j) - peln(i, k, j)
7514 temp_ad20 = -(delp(i, j, k)*pm_ad/temp18**2)
7515 delp_ad(i, j, k) = delp_ad(i, j, k) + pt(i, j, k)*&
7516 & temp_ad21 + pm_ad/temp18
7517 peln_ad(i, k+1, j) = peln_ad(i, k+1, j) + temp_ad20
7518 peln_ad(i, k, j) = peln_ad(i, k, j) - temp_ad20
7520 pt_ad(i, j, k) = pt_ad(i, j, k) + delp(i, j, k)*temp_ad21
7521 delz_ad(i, j, k) = delz_ad(i, j, k) - temp14*temp_ad21
7526 DO i=ilast,ifirst,-1
7528 pe_ad(i, k, j) = pe_ad(i, k, j) + peln_ad(i, k, j)/pe(i, k&
7530 peln_ad(i, k, j) = 0.0
7532 pe_ad(i, k-1, j) = pe_ad(i, k-1, j) + pe_ad(i, k, j)
7533 delp_ad(i, j, k-1) = delp_ad(i, j, k-1) + pe_ad(i, k, j)
7534 pe_ad(i, k, j) = 0.0
7537 DO i=ilast,ifirst,-1
7539 peln_ad(i, 1, j) = 0.0
7541 pe_ad(i, 1, j) = 0.0
7544 DO i=ilast,ifirst,-1
7546 gz_ad(i, j, k+1) = gz_ad(i, j, k+1) + gz_ad(i, j, k)
7547 delz_ad(i, j, k) = delz_ad(i, j, k) -
grav*gz_ad(i, j, k)
7548 gz_ad(i, j, k) = 0.0
7551 DO i=ilast,ifirst,-1
7553 gz_ad(i, j, npz+1) = 0.0
7560 IF (branch .EQ. 0)
THEN 7563 IF (branch .NE. 0)
THEN 7565 DO i=ilast,ifirst,-1
7567 pe_ad(i, k, j) = pe_ad(i, k, j) + kappa*exp(kappa*log(pe&
7568 & (i, k, j)))*pk3_ad(i, j, k)/pe(i, k, j)
7569 pk3_ad(i, j, k) = 0.0
7572 DO i=ilast,ifirst,-1
7574 pk3_ad(i, j, 1) = 0.0
7578 IF (branch .EQ. 0)
THEN 7580 DO i=ilast,ifirst,-1
7582 pe_ad(i, k, j) = pe_ad(i, k, j) + pkc_ad(i, j, k)
7589 DO i=ilast,ifirst,-1
7591 gam_ad(i, k) = gam_ad(i, k) - pkc(i, j, k+1)*pkc_ad(i, j, &
7593 pkc_ad(i, j, k+1) = pkc_ad(i, j, k+1) - gam(i, k)*pkc_ad(i&
7598 DO i=ilast,ifirst,-1
7600 temp_ad18 = pkc_ad(i, j, k+1)/bet(i)
7601 dd_ad(i, k) = dd_ad(i, k) + temp_ad18
7602 pkc_ad(i, j, k) = pkc_ad(i, j, k) - temp_ad18
7603 bet_ad(i) = bet_ad(i) - (dd(i, k)-pkc(i, j, k))*temp_ad18/&
7605 pkc_ad(i, j, k+1) = 0.0
7607 bb_ad(i, k) = bb_ad(i, k) + bet_ad(i)
7608 gam_ad(i, k) = gam_ad(i, k) - bet_ad(i)
7610 temp_ad19 = gam_ad(i, k)/bet(i)
7611 bet_ad(i) = -(g_rat(i, k-1)*temp_ad19/bet(i))
7612 g_rat_ad(i, k-1) = g_rat_ad(i, k-1) + temp_ad19
7616 DO i=ilast,ifirst,-1
7618 pkz_ad(i, npz) = pkz_ad(i, npz) + 3.*dd_ad(i, npz)
7622 temp_ad17 = pkc_ad(i, j, 2)/bet(i)
7623 dd_ad(i, 1) = dd_ad(i, 1) + temp_ad17
7624 bet_ad(i) = bet_ad(i) - dd(i, 1)*temp_ad17/bet(i)
7625 pkc_ad(i, j, 2) = 0.0
7627 pkc_ad(i, j, 1) = 0.0
7629 bb_ad(i, 1) = bb_ad(i, 1) + bet_ad(i)
7633 DO i=ilast,ifirst,-1
7635 temp_ad15 = 3.*dd_ad(i, k)
7636 pkz_ad(i, k) = pkz_ad(i, k) + temp_ad15
7637 g_rat_ad(i, k) = g_rat_ad(i, k) + 2.*bb_ad(i, k) + pkz(i, &
7639 pkz_ad(i, k+1) = pkz_ad(i, k+1) + g_rat(i, k)*temp_ad15
7643 temp_ad16 = g_rat_ad(i, k)/delp(i, j, k+1)
7644 delp_ad(i, j, k) = delp_ad(i, j, k) + temp_ad16
7645 delp_ad(i, j, k+1) = delp_ad(i, j, k+1) - delp(i, j, k)*&
7646 & temp_ad16/delp(i, j, k+1)
7647 g_rat_ad(i, k) = 0.0
7651 DO i=ilast,ifirst,-1
7652 temp12 = delz(i, j, k)
7653 temp11 = delp(i, j, k)*pt(i, j, k)
7654 temp9 = temp11/temp12
7655 temp10 = -(rgrav*
rdgas*temp9)
7656 temp_ad14 = -(rgrav*
rdgas*gama*exp(gama*log(temp10))*&
7657 & pkz_ad(i, k)/(temp10*temp12))
7658 pm_ad = -pkz_ad(i, k)
7659 temp13 = peln(i, k+1, j) - peln(i, k, j)
7660 temp_ad13 = -(delp(i, j, k)*pm_ad/temp13**2)
7661 delp_ad(i, j, k) = delp_ad(i, j, k) + pt(i, j, k)*&
7662 & temp_ad14 + pm_ad/temp13
7663 peln_ad(i, k+1, j) = peln_ad(i, k+1, j) + temp_ad13
7664 peln_ad(i, k, j) = peln_ad(i, k, j) - temp_ad13
7666 pt_ad(i, j, k) = pt_ad(i, j, k) + delp(i, j, k)*temp_ad14
7667 delz_ad(i, j, k) = delz_ad(i, j, k) - temp9*temp_ad14
7672 DO i=ilast,ifirst,-1
7674 pe_ad(i, k, j) = pe_ad(i, k, j) + peln_ad(i, k, j)/pe(i, k&
7676 peln_ad(i, k, j) = 0.0
7678 pe_ad(i, k-1, j) = pe_ad(i, k-1, j) + pe_ad(i, k, j)
7679 delp_ad(i, j, k-1) = delp_ad(i, j, k-1) + pe_ad(i, k, j)
7680 pe_ad(i, k, j) = 0.0
7683 DO i=ilast,ifirst,-1
7685 peln_ad(i, 1, j) = 0.0
7687 pe_ad(i, 1, j) = 0.0
7690 DO i=ilast,ifirst,-1
7692 gz_ad(i, j, k+1) = gz_ad(i, j, k+1) + gz_ad(i, j, k)
7693 delz_ad(i, j, k) = delz_ad(i, j, k) -
grav*gz_ad(i, j, k)
7694 gz_ad(i, j, k) = 0.0
7697 DO i=ilast,ifirst,-1
7699 gz_ad(i, j, npz+1) = 0.0
7704 IF (branch .EQ. 0)
THEN 7705 DO j=jlast,jfirst,-1
7707 IF (branch .NE. 0)
THEN 7711 pe_ad(i, k, j) = pe_ad(i, k, j) + kappa*exp(kappa*log(pe&
7712 & (i, k, j)))*pk3_ad(i, j, k)/pe(i, k, j)
7713 pk3_ad(i, j, k) = 0.0
7718 pk3_ad(i, j, 1) = 0.0
7722 IF (branch .EQ. 0)
THEN 7726 pe_ad(i, k, j) = pe_ad(i, k, j) + pkc_ad(i, j, k)
7731 DO j=jlast,jfirst,-1
7735 gam_ad(i, k) = gam_ad(i, k) - pkc(i, j, k+1)*pkc_ad(i, j, &
7737 pkc_ad(i, j, k+1) = pkc_ad(i, j, k+1) - gam(i, k)*pkc_ad(i&
7744 temp_ad11 = pkc_ad(i, j, k+1)/bet(i)
7745 dd_ad(i, k) = dd_ad(i, k) + temp_ad11
7746 pkc_ad(i, j, k) = pkc_ad(i, j, k) - temp_ad11
7747 bet_ad(i) = bet_ad(i) - (dd(i, k)-pkc(i, j, k))*temp_ad11/&
7749 pkc_ad(i, j, k+1) = 0.0
7751 bb_ad(i, k) = bb_ad(i, k) + bet_ad(i)
7752 gam_ad(i, k) = gam_ad(i, k) - bet_ad(i)
7754 temp_ad12 = gam_ad(i, k)/bet(i)
7755 bet_ad(i) = -(g_rat(i, k-1)*temp_ad12/bet(i))
7756 g_rat_ad(i, k-1) = g_rat_ad(i, k-1) + temp_ad12
7762 pkz_ad(i, npz) = pkz_ad(i, npz) + 3.*dd_ad(i, npz)
7766 temp_ad10 = pkc_ad(i, j, 2)/bet(i)
7767 dd_ad(i, 1) = dd_ad(i, 1) + temp_ad10
7768 bet_ad(i) = bet_ad(i) - dd(i, 1)*temp_ad10/bet(i)
7769 pkc_ad(i, j, 2) = 0.0
7771 pkc_ad(i, j, 1) = 0.0
7773 bb_ad(i, 1) = bb_ad(i, 1) + bet_ad(i)
7779 temp_ad8 = 3.*dd_ad(i, k)
7780 pkz_ad(i, k) = pkz_ad(i, k) + temp_ad8
7781 g_rat_ad(i, k) = g_rat_ad(i, k) + 2.*bb_ad(i, k) + pkz(i, &
7783 pkz_ad(i, k+1) = pkz_ad(i, k+1) + g_rat(i, k)*temp_ad8
7787 temp_ad9 = g_rat_ad(i, k)/delp(i, j, k+1)
7788 delp_ad(i, j, k) = delp_ad(i, j, k) + temp_ad9
7789 delp_ad(i, j, k+1) = delp_ad(i, j, k+1) - delp(i, j, k)*&
7790 & temp_ad9/delp(i, j, k+1)
7791 g_rat_ad(i, k) = 0.0
7796 temp7 = delz(i, j, k)
7797 temp6 = delp(i, j, k)*pt(i, j, k)
7799 temp5 = -(rgrav*
rdgas*temp4)
7800 temp_ad7 = -(rgrav*
rdgas*gama*exp(gama*log(temp5))*pkz_ad(&
7801 & i, k)/(temp5*temp7))
7802 pm_ad = -pkz_ad(i, k)
7803 temp8 = peln(i, k+1, j) - peln(i, k, j)
7804 temp_ad6 = -(delp(i, j, k)*pm_ad/temp8**2)
7805 delp_ad(i, j, k) = delp_ad(i, j, k) + pt(i, j, k)*temp_ad7&
7807 peln_ad(i, k+1, j) = peln_ad(i, k+1, j) + temp_ad6
7808 peln_ad(i, k, j) = peln_ad(i, k, j) - temp_ad6
7810 pt_ad(i, j, k) = pt_ad(i, j, k) + delp(i, j, k)*temp_ad7
7811 delz_ad(i, j, k) = delz_ad(i, j, k) - temp4*temp_ad7
7818 pe_ad(i, k, j) = pe_ad(i, k, j) + peln_ad(i, k, j)/pe(i, k&
7820 peln_ad(i, k, j) = 0.0
7822 pe_ad(i, k-1, j) = pe_ad(i, k-1, j) + pe_ad(i, k, j)
7823 delp_ad(i, j, k-1) = delp_ad(i, j, k-1) + pe_ad(i, k, j)
7824 pe_ad(i, k, j) = 0.0
7829 peln_ad(i, 1, j) = 0.0
7831 pe_ad(i, 1, j) = 0.0
7836 gz_ad(i, j, k+1) = gz_ad(i, j, k+1) + gz_ad(i, j, k)
7837 delz_ad(i, j, k) = delz_ad(i, j, k) -
grav*gz_ad(i, j, k)
7838 gz_ad(i, j, k) = 0.0
7843 gz_ad(i, j, npz+1) = 0.0
7848 IF (branch .EQ. 0)
THEN 7849 DO j=jlast,jfirst,-1
7851 IF (branch .NE. 0)
THEN 7855 pe_ad(i, k, j) = pe_ad(i, k, j) + kappa*exp(kappa*log(pe&
7856 & (i, k, j)))*pk3_ad(i, j, k)/pe(i, k, j)
7857 pk3_ad(i, j, k) = 0.0
7862 pk3_ad(i, j, 1) = 0.0
7866 IF (branch .EQ. 0)
THEN 7870 pe_ad(i, k, j) = pe_ad(i, k, j) + pkc_ad(i, j, k)
7875 DO j=jlast,jfirst,-1
7879 gam_ad(i, k) = gam_ad(i, k) - pkc(i, j, k+1)*pkc_ad(i, j, &
7881 pkc_ad(i, j, k+1) = pkc_ad(i, j, k+1) - gam(i, k)*pkc_ad(i&
7888 temp_ad4 = pkc_ad(i, j, k+1)/bet(i)
7889 dd_ad(i, k) = dd_ad(i, k) + temp_ad4
7890 pkc_ad(i, j, k) = pkc_ad(i, j, k) - temp_ad4
7891 bet_ad(i) = bet_ad(i) - (dd(i, k)-pkc(i, j, k))*temp_ad4/&
7893 pkc_ad(i, j, k+1) = 0.0
7895 bb_ad(i, k) = bb_ad(i, k) + bet_ad(i)
7896 gam_ad(i, k) = gam_ad(i, k) - bet_ad(i)
7898 temp_ad5 = gam_ad(i, k)/bet(i)
7899 bet_ad(i) = -(g_rat(i, k-1)*temp_ad5/bet(i))
7900 g_rat_ad(i, k-1) = g_rat_ad(i, k-1) + temp_ad5
7906 pkz_ad(i, npz) = pkz_ad(i, npz) + 3.*dd_ad(i, npz)
7910 temp_ad3 = pkc_ad(i, j, 2)/bet(i)
7911 dd_ad(i, 1) = dd_ad(i, 1) + temp_ad3
7912 bet_ad(i) = bet_ad(i) - dd(i, 1)*temp_ad3/bet(i)
7913 pkc_ad(i, j, 2) = 0.0
7915 pkc_ad(i, j, 1) = 0.0
7917 bb_ad(i, 1) = bb_ad(i, 1) + bet_ad(i)
7923 temp_ad1 = 3.*dd_ad(i, k)
7924 pkz_ad(i, k) = pkz_ad(i, k) + temp_ad1
7925 g_rat_ad(i, k) = g_rat_ad(i, k) + 2.*bb_ad(i, k) + pkz(i, &
7927 pkz_ad(i, k+1) = pkz_ad(i, k+1) + g_rat(i, k)*temp_ad1
7931 temp_ad2 = g_rat_ad(i, k)/delp(i, j, k+1)
7932 delp_ad(i, j, k) = delp_ad(i, j, k) + temp_ad2
7933 delp_ad(i, j, k+1) = delp_ad(i, j, k+1) - delp(i, j, k)*&
7934 & temp_ad2/delp(i, j, k+1)
7935 g_rat_ad(i, k) = 0.0
7940 temp2 = delz(i, j, k)
7941 temp1 = delp(i, j, k)*pt(i, j, k)
7943 temp0 = -(rgrav*
rdgas*temp)
7944 temp_ad0 = -(rgrav*
rdgas*gama*exp(gama*log(temp0))*pkz_ad(&
7945 & i, k)/(temp0*temp2))
7946 pm_ad = -pkz_ad(i, k)
7947 temp3 = peln(i, k+1, j) - peln(i, k, j)
7948 temp_ad = -(delp(i, j, k)*pm_ad/temp3**2)
7949 delp_ad(i, j, k) = delp_ad(i, j, k) + pt(i, j, k)*temp_ad0&
7951 peln_ad(i, k+1, j) = peln_ad(i, k+1, j) + temp_ad
7952 peln_ad(i, k, j) = peln_ad(i, k, j) - temp_ad
7954 pt_ad(i, j, k) = pt_ad(i, j, k) + delp(i, j, k)*temp_ad0
7955 delz_ad(i, j, k) = delz_ad(i, j, k) - temp*temp_ad0
7961 pe_ad(i, k, j) = pe_ad(i, k, j) + peln_ad(i, k, j)/pe(i, k&
7963 peln_ad(i, k, j) = 0.0
7965 pe_ad(i, k-1, j) = pe_ad(i, k-1, j) + pe_ad(i, k, j)
7966 delp_ad(i, j, k-1) = delp_ad(i, j, k-1) + pe_ad(i, k, j)
7967 pe_ad(i, k, j) = 0.0
7971 peln_ad(i, 1, j) = 0.0
7972 pe_ad(i, 1, j) = 0.0
7977 gz_ad(i, j, k+1) = gz_ad(i, j, k+1) + gz_ad(i, j, k)
7978 delz_ad(i, j, k) = delz_ad(i, j, k) -
grav*gz_ad(i, j, k)
7979 gz_ad(i, j, k) = 0.0
7984 gz_ad(i, j, npz+1) = 0.0
7990 SUBROUTINE nest_halo_nh(ptop, grav, kappa, cp, delp, delz, pt, phis, &
7991 & pkc, gz, pk3, npx, npy, npz, nested, pkc_pertn, computepk3, fullhalo&
7996 INTEGER,
INTENT(IN) :: npx, npy, npz
7997 LOGICAL,
INTENT(IN) :: pkc_pertn, computepk3, fullhalo, nested
7998 REAL,
INTENT(IN) :: ptop, kappa, cp,
grav 8000 REAL,
INTENT(IN) :: phis(bd%isd:bd%ied, bd%jsd:bd%jed)
8001 REAL,
DIMENSION(bd%isd:bd%ied, bd%jsd:bd%jed, npz),
INTENT(IN) :: pt&
8003 REAL,
DIMENSION(bd%isd:bd%ied, bd%jsd:bd%jed, npz+1),
INTENT(INOUT) &
8008 REAL :: ptk, rgrav, rkap, peln1, rdg
8009 REAL,
DIMENSION(bd%isd:bd%ied, npz+1, bd%jsd:bd%jed) :: pe, peln
8010 REAL,
DIMENSION(bd%isd:bd%ied, npz) :: gam, bb, dd, pkz
8011 REAL,
DIMENSION(bd%isd:bd%ied, npz-1) :: g_rat
8012 REAL,
DIMENSION(bd%isd:bd%ied) :: bet
8014 INTEGER :: ifirst, ilast, jfirst, jlast
8015 INTEGER :: is, ie, js, je
8016 INTEGER :: isd, ied, jsd, jed
8027 IF (.NOT.nested)
THEN 8037 gama = 1./(1.-kappa)
8041 rdg = -(
rdgas*rgrav)
8047 gz(i, j, npz+1) = phis(i, j)
8051 gz(i, j, k) = gz(i, j, k+1) - delz(i, j, k)*
grav 8057 peln(i, 1, j) = peln1
8061 pe(i, k, j) = pe(i, k-1, j) + delp(i, j, k-1)
8062 peln(i, k, j) = log(pe(i, k, j))
8069 pkz(i, k) = exp(gama*log(-(delp(i, j, k)*rgrav/delz(i, j, &
8070 & k)*
rdgas*pt(i, j, k))))
8072 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
8074 pkz(i, k) = pkz(i, k) - pm
8080 g_rat(i, k) = delp(i, j, k)/delp(i, j, k+1)
8081 bb(i, k) = 2.*(1.+g_rat(i, k))
8082 dd(i, k) = 3.*(pkz(i, k)+g_rat(i, k)*pkz(i, k+1))
8088 pkc(i, j, 2) = dd(i, 1)/bet(i)
8090 dd(i, npz) = 3.*pkz(i, npz)
8094 gam(i, k) = g_rat(i, k-1)/bet(i)
8095 bet(i) = bb(i, k) - gam(i, k)
8096 pkc(i, j, k+1) = (dd(i, k)-pkc(i, j, k))/bet(i)
8101 pkc(i, j, k) = pkc(i, j, k) - gam(i, k)*pkc(i, j, k+1)
8106 IF (.NOT.pkc_pertn)
THEN 8109 pkc(i, j, k) = pkc(i, j, k) + pe(i, k, j)
8114 IF (computepk3)
THEN 8120 pk3(i, j, k) = exp(kappa*log(pe(i, k, j)))
8126 IF (ie .EQ. npx - 1)
THEN 8130 gz(i, j, npz+1) = phis(i, j)
8134 gz(i, j, k) = gz(i, j, k+1) - delz(i, j, k)*
grav 8140 peln(i, 1, j) = peln1
8144 pe(i, k, j) = pe(i, k-1, j) + delp(i, j, k-1)
8145 peln(i, k, j) = log(pe(i, k, j))
8152 pkz(i, k) = exp(gama*log(-(delp(i, j, k)*rgrav/delz(i, j, &
8153 & k)*
rdgas*pt(i, j, k))))
8155 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
8157 pkz(i, k) = pkz(i, k) - pm
8163 g_rat(i, k) = delp(i, j, k)/delp(i, j, k+1)
8164 bb(i, k) = 2.*(1.+g_rat(i, k))
8165 dd(i, k) = 3.*(pkz(i, k)+g_rat(i, k)*pkz(i, k+1))
8171 pkc(i, j, 2) = dd(i, 1)/bet(i)
8173 dd(i, npz) = 3.*pkz(i, npz)
8177 gam(i, k) = g_rat(i, k-1)/bet(i)
8178 bet(i) = bb(i, k) - gam(i, k)
8179 pkc(i, j, k+1) = (dd(i, k)-pkc(i, j, k))/bet(i)
8184 pkc(i, j, k) = pkc(i, j, k) - gam(i, k)*pkc(i, j, k+1)
8189 IF (.NOT.pkc_pertn)
THEN 8192 pkc(i, j, k) = pkc(i, j, k) + pe(i, k, j)
8197 IF (computepk3)
THEN 8203 pk3(i, j, k) = exp(kappa*log(pe(i, k, j)))
8213 gz(i, j, npz+1) = phis(i, j)
8217 gz(i, j, k) = gz(i, j, k+1) - delz(i, j, k)*
grav 8223 peln(i, 1, j) = peln1
8227 pe(i, k, j) = pe(i, k-1, j) + delp(i, j, k-1)
8228 peln(i, k, j) = log(pe(i, k, j))
8235 pkz(i, k) = exp(gama*log(-(delp(i, j, k)*rgrav/delz(i, j, &
8236 & k)*
rdgas*pt(i, j, k))))
8238 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
8240 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
8242 pkz(i, k) = pkz(i, k) - pm
8248 g_rat(i, k) = delp(i, j, k)/delp(i, j, k+1)
8249 bb(i, k) = 2.*(1.+g_rat(i, k))
8250 dd(i, k) = 3.*(pkz(i, k)+g_rat(i, k)*pkz(i, k+1))
8256 pkc(i, j, 2) = dd(i, 1)/bet(i)
8258 dd(i, npz) = 3.*pkz(i, npz)
8262 gam(i, k) = g_rat(i, k-1)/bet(i)
8263 bet(i) = bb(i, k) - gam(i, k)
8264 pkc(i, j, k+1) = (dd(i, k)-pkc(i, j, k))/bet(i)
8269 pkc(i, j, k) = pkc(i, j, k) - gam(i, k)*pkc(i, j, k+1)
8274 IF (.NOT.pkc_pertn)
THEN 8277 pkc(i, j, k) = pkc(i, j, k) + pe(i, k, j)
8282 IF (computepk3)
THEN 8288 pk3(i, j, k) = exp(kappa*log(pe(i, k, j)))
8294 IF (je .EQ. npy - 1)
THEN 8298 gz(i, j, npz+1) = phis(i, j)
8302 gz(i, j, k) = gz(i, j, k+1) - delz(i, j, k)*
grav 8308 peln(i, 1, j) = peln1
8312 pe(i, k, j) = pe(i, k-1, j) + delp(i, j, k-1)
8313 peln(i, k, j) = log(pe(i, k, j))
8320 pkz(i, k) = exp(gama*log(-(delp(i, j, k)*rgrav/delz(i, j, &
8321 & k)*
rdgas*pt(i, j, k))))
8323 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
8325 pm = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
8327 pkz(i, k) = pkz(i, k) - pm
8334 g_rat(i, k) = delp(i, j, k)/delp(i, j, k+1)
8335 bb(i, k) = 2.*(1.+g_rat(i, k))
8336 dd(i, k) = 3.*(pkz(i, k)+g_rat(i, k)*pkz(i, k+1))
8342 pkc(i, j, 2) = dd(i, 1)/bet(i)
8344 dd(i, npz) = 3.*pkz(i, npz)
8348 gam(i, k) = g_rat(i, k-1)/bet(i)
8349 bet(i) = bb(i, k) - gam(i, k)
8350 pkc(i, j, k+1) = (dd(i, k)-pkc(i, j, k))/bet(i)
8355 pkc(i, j, k) = pkc(i, j, k) - gam(i, k)*pkc(i, j, k+1)
8360 IF (.NOT.pkc_pertn)
THEN 8363 pkc(i, j, k) = pkc(i, j, k) + pe(i, k, j)
8368 IF (computepk3)
THEN 8374 pk3(i, j, k) = exp(kappa*log(pe(i, k, j)))
subroutine, public sim_solver(dt, is, ie, km, rgas, gama, gm2, cp2, kappa, pe2, dm2, pm2, pem, w2, dz2, pt2, ws, alpha, p_fac, scale_m)
subroutine, public sim1_solver_fwd(dt, is, ie, km, rgas, gama, gm2, cp2, kappa, pe, dm2, pm2, pem, w2, dz2, pt2, ws, p_fac)
subroutine, public fill_4corners_fwd(q, dir, bd, npx, npy, sw_corner, se_corner, ne_corner, nw_corner)
subroutine, public update_dz_d_fwd(ndif, damp, hord, is, ie, js, je, km, ng, npx, npy, area, rarea, dp0, zs, zh, crx, cry, xfx, yfx, delz, ws, rdt, gridstruct, bd, hord_pert)
subroutine, public nest_halo_nh(ptop, grav, kappa, cp, delp, delz, pt, phis, pkc, gz, pk3, npx, npy, npz, nested, pkc_pertn, computepk3, fullhalo, bd)
subroutine, public rim_2d(ms, bdt, is, ie, km, rgas, gama, gm2, pe2, dm2, pm2, w2, dz2, pt2, ws, c_core)
subroutine edge_profile_bwd(q1, q1_ad, q2, q2_ad, q1e, q1e_ad, q2e, q2e_ad, i1, i2, j1, j2, j, km, dp0, uniform_grid, limiter)
subroutine, public pushcontrol(ctype, field)
subroutine, public fv_tp_2d_adm(q, q_ad, crx, crx_ad, cry, cry_ad, npx, npy, hord, fx, fx_ad, fy, fy_ad, xfx, xfx_ad, yfx, yfx_ad, gridstruct, bd, ra_x, ra_x_ad, ra_y, ra_y_ad, mfx, mfx_ad, mfy, mfy_ad, mass, mass_ad, nord, damp_c)
subroutine imp_diff_w(j, is, ie, js, je, ng, km, cd, delz, ws, w, w3)
subroutine edge_profile(q1, q2, q1e, q2e, i1, i2, j1, j2, j, km, dp0, uniform_grid, limiter)
subroutine, public rim_2d_fwd(ms, bdt, is, ie, km, rgas, gama, gm2, pe2, dm2, pm2, w2, dz2, pt2, ws, c_core)
real, parameter, public rdgas
Gas constant for dry air [J/kg/deg].
subroutine, public update_dz_d_bwd(ndif, damp, hord, is, ie, js, je, km, ng, npx, npy, area, rarea, dp0, zs, zh, zh_ad, crx, crx_ad, cry, cry_ad, xfx, xfx_ad, yfx, yfx_ad, delz, ws, ws_ad, rdt, gridstruct, bd, hord_pert)
subroutine, public fv_tp_2d_fwd(q, crx, cry, npx, npy, hord, fx, fy, xfx, yfx, gridstruct, bd, ra_x, ra_y, mfx, mfy, mass, nord, damp_c)
subroutine, public sim_solver_fwd(dt, is, ie, km, rgas, gama, gm2, cp2, kappa, pe2, dm2, pm2, pem, w2, dz2, pt2, ws, alpha, p_fac, scale_m)
subroutine, public sim1_solver_bwd(dt, is, ie, km, rgas, gama, gm2, cp2, kappa, pe, pe_ad, dm2, dm2_ad, pm2, pm2_ad, pem, pem_ad, w2, w2_ad, dz2, dz2_ad, pt2, pt2_ad, ws, ws_ad, p_fac)
subroutine edge_profile_fwd(q1, q2, q1e, q2e, i1, i2, j1, j2, j, km, dp0, uniform_grid, limiter)
subroutine, public sim3p0_solver_bwd(dt, is, ie, km, rgas, gama, kappa, pe2, pe2_ad, dm, dm_ad, pem, pem_ad, w2, w2_ad, dz2, dz2_ad, pt2, pt2_ad, ws, ws_ad, p_fac, scale_m)
subroutine, public update_dz_d(ndif, damp, hord, is, ie, js, je, km, ng, npx, npy, area, rarea, dp0, zs, zh, crx, cry, xfx, yfx, delz, ws, rdt, gridstruct, bd, hord_pert)
subroutine, public sim3p0_solver(dt, is, ie, km, rgas, gama, kappa, pe2, dm, pem, w2, dz2, pt2, ws, p_fac, scale_m)
subroutine, public fv_tp_2d(q, crx, cry, npx, npy, hord, fx, fy, xfx, yfx, gridstruct, bd, ra_x, ra_y, mfx, mfy, mass, nord, damp_c)
subroutine, public update_dz_c_fwd(is, ie, js, je, km, ng, dt, dp0, zs, area, ut, vt, gz, ws, npx, npy, sw_corner, se_corner, ne_corner, nw_corner, bd, grid_type)
subroutine, public sim3p0_solver_fwd(dt, is, ie, km, rgas, gama, kappa, pe2, dm, pem, w2, dz2, pt2, ws, p_fac, scale_m)
subroutine, public update_dz_c(is, ie, js, je, km, ng, dt, dp0, zs, area, ut, vt, gz, ws, npx, npy, sw_corner, se_corner, ne_corner, nw_corner, bd, grid_type)
subroutine, public del6_vt_flux(nord, npx, npy, damp, q, d2, fx2, fy2, gridstruct, bd)
subroutine, public sim3_solver_bwd(dt, is, ie, km, rgas, gama, kappa, pe2, pe2_ad, dm, dm_ad, pem, pem_ad, w2, w2_ad, dz2, dz2_ad, pt2, pt2_ad, ws, ws_ad, alpha, p_fac, scale_m)
subroutine, public fill_4corners_bwd(q, q_ad, dir, bd, npx, npy, sw_corner, se_corner, ne_corner, nw_corner)
real, parameter, public cp_air
Specific heat capacity of dry air at constant pressure [J/kg/deg].
subroutine, public sim_solver_bwd(dt, is, ie, km, rgas, gama, gm2, cp2, kappa, pe2, pe2_ad, dm2, dm2_ad, pm2, pm2_ad, pem, pem_ad, w2, w2_ad, dz2, dz2_ad, pt2, pt2_ad, ws, ws_ad, alpha, p_fac, scale_m)
subroutine, public sim1_solver(dt, is, ie, km, rgas, gama, gm2, cp2, kappa, pe, dm2, pm2, pem, w2, dz2, pt2, ws, p_fac)
subroutine, public fv_tp_2d_bwd(q, q_ad, crx, crx_ad, cry, cry_ad, npx, npy, hord, fx, fx_ad, fy, fy_ad, xfx, xfx_ad, yfx, yfx_ad, gridstruct, bd, ra_x, ra_x_ad, ra_y, ra_y_ad, mfx, mfx_ad, mfy, mfy_ad, mass, mass_ad, nord, damp_c)
subroutine, public nest_halo_nh_bwd(ptop, grav, kappa, cp, delp, delp_ad, delz, delz_ad, pt, pt_ad, phis, pkc, pkc_ad, gz, gz_ad, pk3, pk3_ad, npx, npy, npz, nested, pkc_pertn, computepk3, fullhalo, bd)
subroutine, public riem_solver_c_fwd(ms, dt, is, ie, js, je, km, ng, akap, cappa, cp, ptop, hs, w3, pt, q_con, delp, gz, pef, ws, p_fac, a_imp, scale_m)
subroutine, public sim3_solver(dt, is, ie, km, rgas, gama, kappa, pe2, dm, pem, w2, dz2, pt2, ws, alpha, p_fac, scale_m)
subroutine, public riem_solver_c_bwd(ms, dt, is, ie, js, je, km, ng, akap, cappa, cp, ptop, hs, w3, w3_ad, pt, pt_ad, q_con, delp, delp_ad, gz, gz_ad, pef, pef_ad, ws, ws_ad, p_fac, a_imp, scale_m)
real, parameter, public grav
Acceleration due to gravity [m/s^2].
subroutine, public rim_2d_bwd(ms, bdt, is, ie, km, rgas, gama, gm2, pe2, pe2_ad, dm2, dm2_ad, pm2, pm2_ad, w2, w2_ad, dz2, dz2_ad, pt2, pt2_ad, ws, ws_ad, c_core)
subroutine edge_scalar(q1, qe, i1, i2, km, id)
subroutine, public riem_solver_c(ms, dt, is, ie, js, je, km, ng, akap, cappa, cp, ptop, hs, w3, pt, q_con, delp, gz, pef, ws, p_fac, a_imp, scale_m)
subroutine, public sim3_solver_fwd(dt, is, ie, km, rgas, gama, kappa, pe2, dm, pem, w2, dz2, pt2, ws, alpha, p_fac, scale_m)
subroutine, public nest_halo_nh_fwd(ptop, grav, kappa, cp, delp, delz, pt, phis, pkc, gz, pk3, npx, npy, npz, nested, pkc_pertn, computepk3, fullhalo, bd)
subroutine, public update_dz_c_bwd(is, ie, js, je, km, ng, dt, dp0, zs, area, ut, ut_ad, vt, vt_ad, gz, gz_ad, ws, ws_ad, npx, npy, sw_corner, se_corner, ne_corner, nw_corner, bd, grid_type)
subroutine riem_solver3test(ms, dt, is, ie, js, je, km, ng, isd, ied, jsd, jed, akap, cappa, cp, ptop, zs, q_con, w, delz, pt, delp, zh, pe, ppe, pk3, pk, peln, ws, scale_m, p_fac, a_imp, use_logp, last_call, fp_out)
subroutine, public fill_4corners(q, dir, bd, npx, npy, sw_corner, se_corner, ne_corner, nw_corner)
subroutine, public popcontrol(ctype, field)
Derived type containing the data.
subroutine, public del6_vt_flux_adm(nord, npx, npy, damp, q, q_ad, d2, d2_ad, fx2, fx2_ad, fy2, fy2_ad, gridstruct, bd)