!======================================= ! 2D cumulus model - kaminari ! - subroutine base ! ! Author : TAKAHASHI Koko ! Date : 2003/12/04 新規作成 ! Note : 基本場の計算 ! !======================================= subroutine base(im,km,dz,prss_sfc,rd,cp,cv,grv,g_sqrt, & & qv_bs,prss_bs,pi_bs,temp_bs,ptemp_bs, & & vptemp_bs,dens_bs) implicit none integer(8), intent(in) :: im, km real(8), intent(in) :: dz real(8), intent(in) :: prss_sfc, rd, cp, cv, grv real(8), intent(in) :: g_sqrt(-2:im+1,-2:km+1) real(8), intent(out) :: qv_bs(-2:im+1,-2:km+1) real(8), intent(out) :: prss_bs(-2:im+1,-2:km+1) real(8), intent(out) :: pi_bs(-2:im+1,-2:km+1) real(8), intent(out) :: temp_bs(-2:im+1,-2:km+1) real(8), intent(out) :: ptemp_bs(-2:im+1,-2:km+1) real(8), intent(out) :: vptemp_bs(-2:im+1,-2:km+1) real(8), intent(out) :: dens_bs(-2:im+1,-2:km+1) integer(8) :: i, k !--- 全体の温度と水蒸気の混合比を与える do k = 0, km-1 do i = 0, im-2 !--- 温度 temp_bs(i,k) = 3.0d2! - 6.5d-3*z(i,k) !--- 水蒸気の混合比 qv_bs(i,k) = 0.0d0 end do end do !--- 地表面圧力を与える do i =0,im-2 prss_bs(i,0) = prss_sfc end do !--- 地表面での各値の計算 do i = 0, im-2 !--- 無次元圧力 pi_bs(i,0) = (prss_bs(i,0)/prss_sfc)**(rd/cp) !--- 温位 ptemp_bs(i,0) = temp_bs(i,0)/pi_bs(i,0) !--- 仮温位 vptemp_bs(i,0) = ptemp_bs(i,0)*(1.0d0 + 6.1d-2*qv_bs(i,0)) !--- 密度 dens_bs(i,0) = prss_sfc/rd*pi_bs(i,0)**(cv/rd)/vptemp_bs(i,0) end do !--- 地表面より上での各値 ! k = 1 do i = 0, im-2 !--- 静水圧平衡の式 -> 無次元圧力 pi_bs(i,1) = pi_bs(i,0) & & - dz*g_sqrt(i,0)*grv/(cp*vptemp_bs(i,0)) !--- 圧力 prss_bs(i,1) = prss_sfc*pi_bs(i,1)**(cp/rd) !--- 温位 ptemp_bs(i,1) = temp_bs(i,1)/pi_bs(i,1) !--- 仮温位 vptemp_bs(i,1) = ptemp_bs(i,1)*(1.0d0 + 6.1d-2*qv_bs(i,1)) !--- 密度 dens_bs(i,1) = prss_sfc*pi_bs(i,1)**(cv/cp)/(rd*vptemp_bs(i,1)) end do !--- 地表面より上での各値 do k = 2, km-1 do i = 0, im-2 !--- 静水圧平衡の式 -> 無次元圧力 pi_bs(i,k) = pi_bs(i,k-2) & & - 2.0d0*dz*g_sqrt(i,k-1) & & *grv/(cp*vptemp_bs(i,k-1)) !--- 圧力 prss_bs(i,k) = prss_sfc*pi_bs(i,k)**(cp/rd) !--- 温位 ptemp_bs(i,k) = temp_bs(i,k)/pi_bs(i,k) !--- 仮温位 vptemp_bs(i,k) = ptemp_bs(i,k)*(1.0d0 + 6.1d-2*qv_bs(i,k)) !--- 密度 dens_bs(i,k) = prss_sfc*pi_bs(i,k)**(cv/rd)/(rd*vptemp_bs(i,k)) end do end do end subroutine base