!======================================= ! 2D cumulus model - kaminari ! - subroutine vel_u_tb ! ! Author : TAKAHASHI Koko ! Date : 2003/11/20 最終更新 ! 2003/11/11 新規作成 ! Note : x 方向運動方程式の長い時間 ! 間隔で計算する部分 ! !======================================= subroutine vel_u_tb(i,k,im,km,nuh,nuv,cm,ml,dx,dz, & & g_sqrt_x,g_sqrt_xz,g13,g13z,g13xz,u,omg,& & u_old,w_old,km_sub,u_adv,u_cori,u_trb) implicit none integer(8), intent(in) :: i, k integer(8), intent(in) :: im, km real(8), intent(in) :: nuh,nuv,cm,ml real(8), intent(in) :: dx, dz real(8), intent(in) :: g_sqrt_x(0:im+1,0:km+1) real(8), intent(in) :: g_sqrt_xz(0:im+1,0:km+1) real(8), intent(in) :: g13(0:im+1,0:km+1) real(8), intent(in) :: g13z(0:im+1,0:km+1) real(8), intent(in) :: g13xz(0:im+1,0:km+1) real(8), intent(in) :: u(0:im+1,0:km+1) real(8), intent(in) :: omg(0:im+1,0:km+1) real(8), intent(in) :: u_old(0:im+1,0:km+1) real(8), intent(in) :: w_old(0:im+1,0:km+1) real(8), intent(in) :: km_sub(0:im+1,0:km+1) real(8), intent(out) :: u_adv(0:im+1,0:km+1) real(8), intent(out) :: u_cori(0:im+1,0:km+1) real(8), intent(out) :: u_trb(0:im+1,0:km+1) real(8) :: u_dif(0:im+1,0:km+1) !--- 音波減衰項 u_dif(i,k) = nuh*(u(i+1,k) - 2*u(i,k) + u(i-1,k)) & & /dx**2.0D0 & & + nuv*(u(i,k+1) - 2*u(i,k) + u(i,k-1)) & & /dz**2.0D0 !--- 移流項 u_adv(i,k) = u(i,k)*(u(i+1,k) - u(i-1,k))/(4.0D0*dx) & & + ( & & omg(i+1,k+1) + omg(i+1,k-1) & & + omg(i-1,k+1) + omg(i-1,k-1) & & )/4.0D0 & & *(u(i,k+1) - u(i,k-1))/(2.0D0*dz) & & - u_dif(i,k) !--- コリオリ項 u_cori(i,k) = 0.0D0 !--- 拡散項 ! 与える値は ts - dts のときの値 ***_old(i,k) ! u_trb(i,k) = 0.0D0 u_trb(i,k) = 2.0D0 & & *( & & km_sub(i,k) & & *( & & (u_old(i+1,k) - u_old(i,k))/dx & & + g13xz(i+1,k) & & *( & & u_old(i+1,k+1) + u_old(i-1,k+1) & & + u_old(i+1,k-1) + u_old(i-1,k-1) & & )/(4.0d0*dz) & & ) & & - ( & & km_sub(i-1,k) & & *( & & (u_old(i,k) - u_old(i-1,k))/dx & & + g13xz(i,k) & & *( & & u_old(i,k+1) + u_old(i-2,k+1) & & + u_old(i,k-1) + u_old(i-2,k-1) & & )/(4.0d0*dz) & & ) & & ) & & )/dx & & + 2.0d0*g13z(i,k) & & *( & & ( & & km_sub(i+1,k+1) + km_sub(i,k+1) & & + km_sub(i+1,k-1) + km_sub(i,k-1) & & )/4.0d0 & & *( & & ( & & u_old(i+1,k+1) + u_old(i+1,k-1) & & - u_old(i-1,k+1) - u_old(i-1,k-1) & & )/(4.0d0*dx) & & + g13(i,k+1) & & *(u_old(i,k+1) - u_old(i,k))/dz & & ) & & - ( & & km_sub(i+1,k) + km_sub(i,k) & & + km_sub(i+1,k-2) + km_sub(i,k-2) & & )/4.0d0 & & *( & & ( & & u_old(i+1,k) + u_old(i+1,k-2) & & - u_old(i-1,k) - u_old(i-1,k-2) & & )/(4.0d0*dx) & & + g13(i,k) & & *(u_old(i,k) - u_old(i,k-1))/dz & & ) & & )/dz & & + 1.0d0/g_sqrt_x(i,k) & & *( & & ( & & km_sub(i+1,k+1) + km_sub(i-1,k+1) & & + km_sub(i+1,k) + km_sub(i-1,k) & & )/4.0d0 & & *( & & (w_old(i+1,k+1) -w_old(i,k+1))/dx & & + 1.0d0/g_sqrt_xz(i,k+1) & & *(u_old(i,k+1) - u_old(i,k))/dz & & + g13(i,k+1) & & *( & & w_old(i+1,k+1) + w_old(i-1,k+1) & & - w_old(i+1,k-1) - w_old(i-1,k-1) & & )/(4.0d0*dz) & & ) & & - ( & & km_sub(i+1,k) + km_sub(i-1,k) & & + km_sub(i+1,k-1) + km_sub(i-1,k-1) & & )/4.0d0 & & *( & & (w_old(i+1,k) -w_old(i,k))/dx & & + 1.0d0/g_sqrt_xz(i,k) & & *(u_old(i,k) - u_old(i,k-1))/dz & & + g13(i,k) & & *( & & w_old(i+1,k) + w_old(i-1,k) & & - w_old(i+1,k-2) + w_old(i-1,k-2)& & )/(4.0d0*dz) & & ) & & )/dz & & - 2.0d0/(3.0d0*(cm*ml)**2.0d0) & & *( & & (km_sub(i,k)**2.0d0 - km_sub(i,k-1)**2.0d0)& & /dx & & + g13z(i,k) & & *( & & km_sub(i+1,k+1) + km_sub(i-1,k+1) & & - km_sub(i+1,k-1) - km_sub(i-1,k-1) & & )/(4.0d0*dz) & & ) end subroutine vel_u_tb