!======================================= ! 2D cumulus model - kaminari ! - subroutine vel_u_tb ! ! Author : TAKAHASHI Koko ! Date : 2004/01/06 最終更新 ! 2003/11/11 新規作成 ! Note : x 方向運動方程式の長い時間 ! 間隔で計算する部分 ! !======================================= !--- 1,2,3行目入力変数, 4行目出力変数 subroutine vel_u_tb(im,km,dx,dz,nuh,nuv,cm,ml, & & g_sqrt_x,g_sqrt_xz,g13,g13z,g13xz, & & u_old,omg_old,u_old2,w_old2,km_sub_old2,& & u_adv,u_cori,u_trb) implicit none integer(8), intent(in) :: im, km real(8), intent(in) :: dx, dz real(8), intent(in) :: nuh, nuv ,cm, ml real(8), intent(in) :: g_sqrt_x(-2:im+2,-2:km+2) real(8), intent(in) :: g_sqrt_xz(-2:im+2,-2:km+2) real(8), intent(in) :: g13(-2:im+2,-2:km+2) real(8), intent(in) :: g13z(-2:im+2,-2:km+2) real(8), intent(in) :: g13xz(-2:im+2,-2:km+2) real(8), intent(in) :: u_old(-2:im+2,-2:km+2) 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) :: km_sub_old2(-2:im+2,-2:km+2) real(8), intent(out) :: u_adv(-2:im+2,-2:km+2) real(8), intent(out) :: u_cori(-2:im+2,-2:km+2) real(8), intent(out) :: u_trb(-2:im+2,-2:km+2) integer(8) :: i, k real(8) :: u_dif(-2:im+2,-2:km+2) do k = 0,km do i = 0,im !--- 数値粘性項 u_dif(i,k) = nuh*( & & u_old(i,k) - 2.0d0*u_old(i-1,k) & & + u_old(i-2,k) & & )/(dx**2.0d0) & & + nuv*( & & u_old(i,k) - 2.0d0*u_old(i,k-1) & & + u_old(i,k-2) & & )/(dz**2.0d0) !--- 移流項 u_adv(i,k) = u_old(i,k) & & *(u_old(i+1,k) - u_old(i-1,k)) & & /(2.0d0*dx) & & + ( & & omg_old(i,k+1) + omg_old(i,k) & & + omg_old(i-1,k+1) + omg_old(i-1,k) & & )/4.0d0 & & *(u_old(i,k+1) - u_old(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) = 2.0d0 & & *( & & km_sub_old2(i+1,k) & & *( & & (u_old2(i+1,k) - u_old2(i,k))/dx & & + g13xz(i,k) & & *( & & u_old2(i+1,k+1) & & + u_old2(i,k+1) & & - u_old2(i+1,k-1) & & - u_old2(i,k-1) & & )/(4.0d0*dz) & & ) & & - ( & & km_sub_old2(i,k) & & *( & & (u_old2(i,k) - u_old2(i-1,k))/dx& & + g13xz(i-1,k) & & *( & & u_old2(i,k+1) & & + u_old2(i-1,k+1) & & - u_old2(i,k-1) & & - u_old2(i-1,k-1) & & )/(4.0d0*dz) & & ) & & ) & & )/dx & & + 2.0d0*g13z(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 & & *( & & ( & & u_old2(i+1,k+1) & & + u_old2(i+1,k) & & - u_old2(i-1,k+1) & & - u_old2(i-1,k) & & )/(4.0d0*dx) & & + g13(i,k+1) & & *(u_old2(i,k+1) - u_old2(i,k)) & & /dz & & ) & & - ( & & km_sub_old2(i+1,k) & & + km_sub_old2(i,k) & & + km_sub_old2(i+1,k-1) & & + km_sub_old2(i,k-1) & & )/4.0d0 & & *( & & ( & & u_old2(i+1,k) & & + u_old2(i+1,k-1) & & - u_old2(i-1,k) & & - u_old2(i-1,k-1) & & )/(4.0d0*dx) & & + g13(i,k) & & *(u_old2(i,k) - u_old2(i,k-1))& & /dz & & ) & & )/dz & & + 1.0d0/g_sqrt_x(i,k) & & *( & & ( & & 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+1) -w_old2(i-1,k+1)) & & /dx & & + 1.0d0/g_sqrt_xz(i,k+1) & & *(u_old2(i,k+1) - u_old2(i,k)) & & /dz & & + g13(i,k+1) & & *( & & w_old2(i,k+2) & & + w_old2(i-1,k+2) & & - w_old2(i,k) & & - w_old2(i-1,k) & & )/(4.0d0*dz) & & ) & & - ( & & km_sub_old2(i-1,k) & & + km_sub_old2(i-2,k) & & + km_sub_old2(i-1,k-1) & & + km_sub_old2(i-2,k-1) & & )/4.0d0 & & *( & & (w_old2(i,k) -w_old2(i-1,k))/dx& & + 1.0d0/g_sqrt_xz(i,k) & & *(u_old2(i,k) - u_old2(i,k-1))& & /dz & & + g13(i,k) & & *( & & w_old2(i,k+1) & & + w_old2(i-1,k+1) & & - w_old2(i,k-2) & & - w_old2(i-1,k-2) & & )/(4.0d0*dz) & & ) & & )/dz & & - 2.0d0/(3.0d0*(cm*ml)**2.0d0) & & *( & & ( & & km_sub_old2(i,k)**2.0d0 & & - km_sub_old2(i-1,k)**2.0d0 & & )/dx & & + g13z(i,k) & & *( & & km_sub_old2(i,k+1)**2.0d0 & & + km_sub_old2(i-1,k+1)**2.0d0 & & - km_sub_old2(i,k-1)**2.0d0 & & - km_sub_old2(i-1,k-1)**2.0d0 & & )/(4.0d0*dz) & & + 1.0d0/g_sqrt_x(i,k) & & *( & & km_sub_old2(i,k+1)**2.0d0 & & + km_sub_old2(i-1,k+1)**2.0d0 & & - km_sub_old2(i,k-1)**2.0d0 & & - km_sub_old2(i-1,k-1)**2.0d0 & & )/(4.0d0*dz) & & ) end do end do end subroutine vel_u_tb