!======================================= ! 2D cumulus model - kaminari ! - subroutine eddymx_sub ! ! Author : TAKAHASHI Koko ! Date : 2003/11/21 最終更新 ! 2003/11/05 新規作成 ! Note : 2ステップ目のサブグリッド ! 運動エネルギー ! !======================================= subroutine eddymx_sub(i,k,im,km,nuh,nuv,eps,lv,rd,cp,cm,ml, & & grv,dtb,dx,dz,ptemp_bs,u,w,omg,temp, & & ptemp,qv,qc,g_sqrt,g_sqrt_z,g13x,g13z,g13xz, & & e_sub_old2,e_sub) implicit none integer, intent(in) :: i, k integer, intent(in) :: im, km real(8), intent(in) :: nuh,nuv,eps, lv, rd, cp, cm, ml real(8), intent(in) :: grv, dtb real(8), intent(in) :: dx, dz real(8), intent(in) :: u(-2:im+1,-2:km+1) real(8), intent(in) :: w(-2:im+1,-2:km+1) real(8), intent(in) :: omg(-2:im+1,-2:km+1) real(8), intent(in) :: ptemp_bs(-2:im+1,-2:km+1) real(8), intent(in) :: temp(-2:im+1,-2:km+1) real(8), intent(in) :: ptemp(-2:im+1,-2:km+1) real(8), intent(in) :: qv(-2:im+1,-2:km+1) real(8), intent(in) :: qc(-2:im+1,-2:km+1) real(8), intent(in) :: e_sub_old2(-2:im+1,-2:km+1) real(8), intent(in) :: g_sqrt(-2:im+1,-2:km+1) real(8), intent(in) :: g_sqrt_z(-2:im+1,-2:km+1) real(8), intent(in) :: g13x(-2:im+1,-2:km+1) real(8), intent(in) :: g13z(-2:im+1,-2:km+1) real(8), intent(in) :: g13xz(-2:im+1,-2:km+1) real(8), intent(inout) :: e_sub(-2:im+1,-2:km+1) real(8) :: e_sub_dif(-2:im+1,-2:km+1) real(8) :: e_sub_adv(-2:im+1,-2:km+1) real(8) :: a(-2:im+1,-2:km+1) real(8) :: ptemp_e(-2:im+1,-2:km+1) real(8) :: ql(-2:im+1,-2:km+1) ! real(8) :: (-2:im+1,-2:km+1) !--- サブグリッド運動エネルギーの数値粘性項 e_sub_dif(i,k) = nuh & & *( & & e_sub(i+1,k) - 2*e_sub(i,k) & & + e_sub(i-1,k) & & ) & & /(4*dx**2) & & + nuv & & *( & & e_sub(i,k+1) - 2*e_sub(i,k) & & + e_sub(i,k-1) & & )/(4*dz**2) e_sub_dif(i,k) = 0.0D0 !--- サブグリッド運動エネルギーの移流項 e_sub_adv(i,k) = (u(i+1,k) + u(i,k))/2 & & *(e_sub(i+1,k) - e_sub(i,k))/(4*dx) & & + (omg(i,k+1) + omg(i,k))/2 & & *(e_sub(i,k+1) - e_sub(i,k))/(4*dz) & & + e_sub_dif(i,k) e_sub_adv(i,k) = 0.0D0 !--- 係数 A a(i,k) = 1/ptemp_bs(i,k) & & *(1 + 1.61*eps*lv*qv(i,k)/(rd*temp(i,k))) & & /(1 + eps*(lv**2)*qv(i,k)/(cp*rd*temp(i,k)**2)) !--- 相当温位 ptemp_e(i,k) = ptemp(i,k)*(1 + lv*qv(i,k)/(cp*temp(i,k))) !--- 水蒸気と雲水混合比の和 ql(i,k) = qv(i,k) + qc(i,k) !-- サブグリッド運動エネルギー ! e_sub(i,k) は t での値. 本当は t - dt の値を代入しなければならない e_sub(i,k) = e_sub_old2(i,k) & & + 2*dtb & & *(- e_sub_adv(i,k) & & + 3*grv*cm*ml*e_sub(i,k)**(1/2)/g_sqrt(i,k) & & *(- a(i,k)*ptemp_e(i,k+1) & & - ptemp_e(i,k-1)/(4*dz) & & + (ql(i,k+1) - ql(i,k-1))/(4*dz) & & ) & & + 2*cm*ml*e_sub(i,k)**(1/2) & & *( & & ((u(i+1,k) - u(i,k))/(2*dx) & & + g13xz(i,k) & & *(u(i+1,k+1) + u(i,k+1) & & - u(i+1,k-1)- u(i,k-1) & & )/(8*dz) & & )**2 & & + ( & & 1/g_sqrt(i,k) & & *(w(i,k+1) - w(i,k)) & & /(2*dz) & & )**2 & & + (w(i+1,k+1) + w(i+1,k-1) & & - w(i-1,k+1) - w(i-1,k-1) & & )/(8*dx) & & + g13xz(i,k) & & *(w(i,k+1) - w(i,k)) & & /(2*dz) & & + 1/g_sqrt(i,k) & & *( & & u(i+1,k+1) + u(i-1,k+1) & & - u(i+1,k-1) - u(i-1,k-1) & & )/(8*dz) & & ) & & + cm*ml & & /(2*e_sub(i,k)**(1/2)) & & *( & & ( & & e_sub(i+1,k)**2 & & - e_sub(i,k)**2 & & )/dx**2 & & + ( & & g13z(i+1,k) & & *( & & e_sub(i+2,k+1)**2 & & + e_sub(i,k+1)**2 & & + e_sub(i+2,k-1)**2 & & + e_sub(i,k-1)**2 & & )/(8*dz) & & - g13z(i,k) & & *( & & e_sub(i+1,k+1)**2 & & + e_sub(i-1,k+1)**2 & & + e_sub(i+1,k-1)**2 & & + e_sub(i-1,k-1)**2 & & )/(8*dz) & & )/(2*dx) & & + ( & & g13x(i,k+1) & & *( & & e_sub(i+1,k+2)**2 & & + e_sub(i+1,k)**2 & & + e_sub(i-1,k+2)**2 & & + e_sub(i-1,k)**2 & & )/(8*dx) & & - g13x(i,k) & & *( & & e_sub(i+1,k+1)**2 & & + e_sub(i+1,k-1)**2 & & + e_sub(i-1,k+1)**2 & & + e_sub(i-1,k-1)**2 & & )/(8*dx) & & )/(2*dz) & & + g13xz(i,k) & & *( & & g13x(i+1,k) & & *( & & e_sub(i,k+2)**2 & & - e_sub(i,k+1)**2 & & )/(2*dz) & & - g13x(i,k) & & *( & & e_sub(i,k+1)**2 & & - e_sub(i,k)**2 & & )/(2*dz) & & )/(2*dz) & & + 1/g_sqrt(i,k) & & *( & & g_sqrt_z(i,k+1) & & *( & & e_sub(i,k+2)**2 & & - e_sub(i,k+1)**2 & & )/(2*dz) & & - g_sqrt_z(i,k+1) & & *( & & e_sub(i,k+1)**2 & & - e_sub(i,k)**2 & & )/(2*dz) & & ) & & - ( & & ( & & e_sub(i+1,k)**2 & & - e_sub(i-1,k) & & )/(4*dx) & & + g13xz(i,k) & & *( & & e_sub(i,k+1)**2 & & - e_sub(i,k-1) & & )/(4*dz) & & )**2 & & - ( & & 1/g_sqrt_z(i,k) & & *( & & e_sub(i,k+1) & & - e_sub(i,k-1) & & )/(2*dz) & & )**2 & & ) & & + cm/ml*e_sub(i,k)**(3/2) & & ) end subroutine eddymx_sub