!======================================= ! 2D cumulus model - kaminari ! - subroutine vel_w_tb ! ! Author : TAKAHASHI Koko ! Date : 2004/01/06 最終更新 ! 2003/11/11 新規作成 ! Note : z 方向運動方程式の長い時間 ! 間隔で計算する部分 ! !======================================= !--- 1,2,3,4行目入力変数, 5行目出力変数 subroutine vel_w_tb(im,km,dx,dz,nuh,nuv,grv,cm,ml,ptemp_bs, & & g_sqrt,g_sqrt_x,g_sqrt_z,g13,g13x,g13xz,& & u_old,w_old,omg_old,u_old2,w_old2, & & ptemp,km_sub_old2, & & w_adv,w_buoy,w_trb) implicit none integer(8), intent(in) :: im, km real(8), intent(in) :: dx, dz real(8), intent(in) :: nuh,nuv,grv,cm,ml real(8), intent(in) :: g_sqrt(-2:im+2,-2:km+2) real(8), intent(in) :: g_sqrt_x(-2:im+2,-2:km+2) real(8), intent(in) :: g_sqrt_z(-2:im+2,-2:km+2) real(8), intent(in) :: g13(-2:im+2,-2:km+2) real(8), intent(in) :: g13x(-2:im+2,-2:km+2) real(8), intent(in) :: g13xz(-2:im+2,-2:km+2) real(8), intent(in) :: ptemp_bs(-2:im+2,-2:km+2) ! real(8), intent(in) :: qv_bs(-2:im+2,-2:km+2) real(8), intent(in) :: u_old(-2:im+2,-2:km+2) real(8), intent(in) :: w_old(-2:im+2,-2:km+3) real(8), intent(in) :: omg_old(-2:im+2,-2:km+3) real(8), intent(in) :: u_old2(-2:im+2,-2:km+2) real(8), intent(in) :: w_old2(-2:im+2,-2:km+3) real(8), intent(in) :: ptemp(-2:im+2,-2:km+2) ! real(8), intent(in) :: qv(-2:im+2,-2:km+2) ! real(8), intent(in) :: qc(-2:im+2,-2:km+2) ! real(8), intent(in) :: qr(-2:im+2,-2:km+2) real(8), intent(in) :: km_sub_old2(-2:im+2,-2:km+2) real(8), intent(out) :: w_adv(-2:im+2,-2:km+3) real(8), intent(out) :: w_buoy(-2:im+2,-2:km+3) real(8), intent(out) :: w_trb(-2:im+2,-2:km+3) integer(8) :: i, k real(8) :: w_dif(-2:im+2,-2:km+3) do k = 1,km do i = 0,im !--- 数値粘性項 w_dif(i,k) = nuh*( & & w_old(i,k) - 2.0d0*w_old(i-1,k) & & + w_old(i-2,k) & & )/(dx**2.0d0) & & + nuv*( & & w_old(i,k) - 2.0d0*w_old(i,k-1) & & + w_old(i,k-2) & & )/(dz**2.0d0) !--- 移流項 w_adv(i,k) = ( & & u_old(i,k) + u_old(i,k-1) & & + u_old(i+1,k) + u_old(i+1,k-1) & & )/4.0d0 & & *(w_old(i+1,k) - w_old(i-1,k)) & & /(2.0d0*dx) & & + omg_old(i,k) & & *(w_old(i,k+1) - w_old(i,k-1)) & & /(2.0d0*dz) & & + w_dif(i,k) !--- 浮力項 w_buoy(i,k) = grv & & *( & & (ptemp(i,k) + ptemp(i,k-1))/2.0d0 & & /( & & (ptemp_bs(i,k) + ptemp_bs(i,k-1)) & & /2.0d0 & & ) & & - 1.0d0 & ! & + 0.61d0 & ! & *( & ! & (qv(i,k+1) + qv(i,k))/2.0d0 & ! & - (qv_bs(i,k+1) - qv_bs(i,k)) & ! & /2.0d0 & ! & ) & ! & - (qc(i,k+1) + qc(i,k))/2.0d0 & ! & - (qr(i,k+1) + qr(i,k))/2.0d0 & & ) !--- 拡散項 w_trb(i,k) = ( & & ( & & ( & & km_sub_old2(i+1,k+1) & & + km_sub_old2(i,k+1) & & + km_sub_old2(i+1,k) & & + km_sub_old2(i,k) & & )/4.0d0 & & *( & & (w_old2(i+1,k) - w_old2(i,k))/dx & & + 1.0d0/g_sqrt_x(i+1,k) & & *( & & u_old2(i+1,k) + u_old2(i,k) & & - u_old2(i+1,k-1) & & - u_old2(i,k-1) & & )/(4.0d0*dz) & & + g13(i+1,k) & & *( & & w_old2(i+1,k+1) + w_old2(i,k+1)& & - w_old2(i+1,k-1) & & - w_old2(i,k-1) & & )/(4.0d0*dz) & & ) & & ) & & - ( & & ( & & km_sub_old2(i,k+1) & & + km_sub_old2(i-1,k+1) & & + km_sub_old2(i,k) & & + km_sub_old2(i-1,k) & & )/4.0d0 & & *( & & (w_old2(i,k) - w_old2(i-1,k))/dx & & + 1.0d0/g_sqrt_x(i,k) & & *( & & u_old2(i,k) + u_old2(i-1,k) & & - u_old2(i,k-1) & & - u_old2(i-1,k-1) & & )/(4.0d0*dz) & & + g13(i,k) & & *( & & w_old2(i,k) + w_old2(i-1,k) & & - w_old2(i,k-1) & & - w_old2(i-1,k-1) & & )/(4.0d0*dz) & & ) & & ) & & )/dx & & + g13x(i,k) & & *( & & km_sub_old2(i,k) & & *( & & ( & & w_old2(i,k+1) & & + w_old2(i,k) & & - w_old2(i-1,k+1) & & - w_old2(i-1,k) & & )/(4.0d0*dx) & & + 1.0d0/g_sqrt(i,k) & & *( & & u_old2(i+1,k+1) & & + u_old2(i,k+1) & & - u_old2(i+1,k) & & - u_old2(i,k) & & )/(4.0d0*dz) & & + g13xz(i,k+1) & & *(w_old2(i,k+1) - w_old2(i,k)) & & /dz & & ) & & - km_sub_old2(i,k-1) & & *( & & ( & & w_old2(i,k) & & + w_old2(i,k-1) & & - w_old2(i-1,k) & & - w_old2(i-1,k-1) & & )/(4.0d0*dx) & & + 1.0d0/g_sqrt(i,k) & & *( & & u_old2(i+1,k) & & + u_old2(i,k) & & - u_old2(i+1,k-1) & & - u_old2(i,k-1) & & )/(4.0d0*dz) & & + g13xz(i,k) & & *(w_old2(i,k) - w_old2(i,k-1))& & /dz & & ) & & )/dz & & + 2.0d0/g_sqrt_z(i,k) & & *( & & km_sub_old2(i,k)/g_sqrt(i,k) & & *(w_old2(i,k+1) -w_old2(i,k))/dz & & - km_sub_old2(i,k-1)/g_sqrt(i,k-1) & & *(w_old2(i,k) -w_old2(i,k-1))/dz & & )/dz & & - 2.0d0/(3.0d0*(cm*ml)**2.0d0) & & *( & & ( & & km_sub_old2(i,k)**2.0d0 & & + km_sub_old2(i,k-1)**2.0d0 & & - km_sub_old2(i-1,k)**2.0d0 & & - km_sub_old2(i-1,k-1)**2.0d0 & & )/(4.0d0*dx) & & + g13x(i,k) & & *( & & km_sub_old2(i,k)**2.0d0 & & - km_sub_old2(i,k-1)**20.d0 & & )/dz & & + 1.0d0/g_sqrt_z(i,k) & & *( & & km_sub_old2(i,k)**2.0d0 & & - km_sub_old2(i,k-1)**2.0d0 & & )/dz & & ) end do end do end subroutine vel_w_tb