!= 初期値データファイル生成主プログラム ! != dcpam-primitive main program for generation of initial data file ! ! Authors:: Yoshiyuki O. Takahashi, Yasuhiro MORIKAWA ! modified by Shin-ichi Takehiro for primitive model. ! Version:: $Id: init_data_primitive_PST2004.f90,v 1.1.1.1 2010-08-17 05:24:50 takepiro Exp $ ! Tag Name:: $Name: $ ! Copyright:: Copyright (C) GFD Dennou Club, 2008. All rights reserved. ! License:: See COPYRIGHT[link:../../../COPYRIGHT] ! program init_data_primitive ! ! Note that Japanese and English are described in parallel. ! ! 初期値データファイルを生成します. ! ! Initial data file is created. ! ! モジュール引用 ; USE statements ! !!$ ! 初期値データ (リスタートデータ) 提供 !!$ ! Prepare initial data (restart data) !!$ ! !!$ use initial_data, only: SetInitData ! リスタートデータ入出力 ! Restart data input/output ! use restart_file_io, only: RestartFileOutPut ! 種別型パラメタ ! Kind type parameter ! use dc_types, only: DP, & ! 倍精度実数型. Double precision. & STRING, & ! 文字列. Strings. & TOKEN ! キーワード. Keywords. ! 宣言文 ; Declaration statements ! implicit none ! 予報変数 (ステップ $ t-\Delta t $ , $ t $ , $ t+\Delta t $ ) ! Prediction variables (Step $ t-\Delta t $ , $ t $ , $ t+\Delta t $ ) ! real(DP), allocatable:: xyz_U (:,:,:) ! $ u $ . 東西風速. Eastward wind real(DP), allocatable:: xyz_V (:,:,:) ! $ v $ . 南北風速. Northward wind real(DP), allocatable:: xyz_Temp (:,:,:) ! $ T $ . 温度. Temperature real(DP), allocatable:: xy_Ps (:,:) ! $ p_s $ . 地表面気圧. Surface pressure ! 作業変数 ! Work variables ! ! 実行文 ; Executable statement ! ! 主プログラムの初期化 (内部サブルーチン) ! Initialization for the main program (Internal subroutine) ! call MainInit ! 初期値データの作成 ! Generate initial data ! call SetInitData( & & xyz_U, xyz_V, xyz_Temp, xy_Ps ) ! (out) ! 初期値データの出力 ! Output initial data ! call RestartFileOutput( & & xyz_U, xyz_V, xyz_Temp, xy_Ps ) ! (in) ! 主プログラムの終了処理 (内部サブルーチン) ! Termination for the main program (Internal subroutine) ! call MainTerminate contains !------------------------------------------------------------------- subroutine MainInit ! ! 主プログラムの初期化手続き. ! ! Initialization procedure for the main program. ! ! MPI ! use mpi_wrapper, only : MPIWrapperInit use dc_message, only: MessageNotify ! コマンドライン引数処理 ! Command line option parser ! use option_parser, only: OptParseInit ! NAMELIST ファイル入力に関するユーティリティ ! Utilities for NAMELIST file input ! use namelist_util, only: NmlutilInit ! 時刻管理 ! Time control ! use timeset, only: TimesetInit ! ステップ $ t $ の時刻. Time of step $ t $. ! 出力ファイルの基本情報管理 ! Management basic information for output files ! use fileset, only: FilesetInit ! 格子点設定 ! Grid points settings ! use gridset, only: GridsetInit, & & imax, & ! 経度格子点数. ! Number of grid points in longitude & jmax, & ! 緯度格子点数. ! Number of grid points in latitude & kmax ! 鉛直層数. ! Number of vertical level ! 物理定数設定 ! Physical constants settings ! use constants, only: ConstantsInit ! 座標データ設定 ! Axes data settings ! use axesset, only: AxessetInit ! リスタートデータ入出力 ! Restart data input/output ! use restart_file_io, only: RestartFileOpen, RestartFileGet ! 宣言文 ; Declaration statements ! implicit none character(*), parameter:: prog_name = 'init_data' ! 主プログラム名. ! Main program name character(*), parameter:: version = & & '$Name: $' // & & '$Id: init_data_primitive_PST2004.f90,v 1.1.1.1 2010-08-17 05:24:50 takepiro Exp $' ! 主プログラムのバージョン ! Main program version character(STRING):: brief ! 実行ファイルの簡潔な説明 ! Brief account of executable file ! 実行文 ; Executable statement ! ! Initialize MPI ! call MPIWrapperInit brief = 'Initial data generation' call MessageNotify( 'M', prog_name, 'Run: %c', c1 = trim(brief) ) call MessageNotify( 'M', prog_name, '-- version = %c', c1 = trim(version) ) ! コマンドライン引数処理 ! Command line option parser ! call OptParseInit(prog_name, brief) ! NAMELIST ファイル名入力 ! Input NAMELIST file name ! call NmlutilInit ! 時刻管理 ! Time control ! call TimesetInit ! 出力ファイルの基本情報管理 ! Management basic information for output files ! call FilesetInit ! 格子点設定 ! Grid points settings ! call GridsetInit ! 物理定数設定 ! Physical constants settings ! call ConstantsInit ! 座標データ設定 ! Axes data settings ! call AxessetInit ! 予報変数の割付 ! Allocation of prediction variables ! allocate( xyz_U (0:imax-1, 1:jmax, 1:kmax) ) allocate( xyz_V (0:imax-1, 1:jmax, 1:kmax) ) allocate( xyz_Temp (0:imax-1, 1:jmax, 1:kmax) ) allocate( xy_Ps (0:imax-1, 1:jmax) ) ! 初期値データ出力 ! Initial data output ! call RestartFileOpen( flag_init_data = .true. ) ! (in) optional end subroutine MainInit !------------------------------------------------------------------- subroutine MainTerminate ! ! 主プログラムの終了処理手続き. ! ! Termination procedure for the main program. ! ! MPI ! use mpi_wrapper, only : MPIWrapperFinalize ! 時刻管理 ! Time control ! use timeset, only: TimesetClose ! リスタートデータ入出力 ! Restart data input/output ! use restart_file_io, only: RestartFileClose ! 宣言文 ; Declaration statements ! implicit none ! 実行文 ; Executable statement ! ! リスタートデータファイルクローズ ! Close restart data input ! call RestartFileClose ! 時刻管理終了処理 ! Termination of time control ! call TimesetClose ! 予報変数の割付解除 ! Deallocation of prediction variables ! deallocate( xyz_U ) deallocate( xyz_V ) deallocate( xyz_Temp ) deallocate( xy_Ps ) ! Finalize MPI ! call MPIWrapperFinalize end subroutine MainTerminate !------------------------------------------------------------------- subroutine SetInitData( & & xyz_U, xyz_V, xyz_Temp, xy_Ps ) ! (out) ! ! 初期データの設定 ! ! Setting initial values ! ! 格子点設定 ! Grid points settings ! use gridset, only: & & imax, & ! 経度格子点数. ! Number of grid points in longitude & jmax, & ! 緯度格子点数. ! Number of grid points in latitude & kmax ! 鉛直層数. ! Number of vertical level ! 物理定数設定 ! Physical constants settings ! use constants, only: PI, Omega, RPlanet, GasRDry ! 座標データ設定 ! Axes data settings ! use axesset, only: x_Lon, y_Lat, z_Sigma, y_Lat_Weight use wa_module, only : a_AvrLat_ya real(DP), intent(OUT) :: xyz_U (0:,:,:) ! $ u $ . 東西風速. Eastward wind real(DP), intent(OUT) :: xyz_V (0:,:,:) ! $ v $ . 南北風速. Northward wind real(DP), intent(OUT) :: xyz_Temp (0:,:,:) ! $ T $ . 温度. Temperature real(DP), intent(OUT) :: xy_Ps (0:,:) ! $ p_s $ . 地表面気圧. Surface pressure integer, parameter :: nmax=100 integer :: i, j, k, n real(DP) :: lat, dlat real(DP) :: z_Z(1:kmax) real(DP) :: z_F(1:kmax) real(DP) :: z_dFdZ(1:kmax) real(DP) :: z_Temp(1:kmax) real(DP) :: yz_Temp(1:jmax,1:kmax) real(DP), parameter :: u0 = 50.0d0 real(DP), parameter :: H = 7.34d3 real(DP), parameter :: z0 = 22.0d3 real(DP), parameter :: z1 = 30.0d3 real(DP), parameter :: dz0 = 5.0d3 real(DP), parameter :: p0 = 1.0d5 real(DP), parameter :: That = 1.0d0 real(DP), parameter :: lambda0 = 0.0d0 real(DP), parameter :: phi0 = PI/4 real(DP), parameter :: alpha = 1.0d0/3.0d0 real(DP), parameter :: beta = 1.0d0/6.0d0 xyz_V = 0.0D0 xy_Ps = p0 z_Z = - H*log(z_Sigma) z_F = 0.5D0 * ( 1 - tanh((z_Z-z0)/dz0)**3 ) * sin(PI*z_Z/z1) ! ! zonal velocity profile ! do k=1,kmax do j=1,jmax do i=0,imax-1 if ( y_Lat(j) > 0.0d0 ) then xyz_U(i,j,k) = u0 * sin(PI*sin(y_Lat(j))**2)**3 * z_F(k) else xyz_U(i,j,k) = 0.0d0 endif end do end do end do ! ! temprature profile ! z_dFdZ =-1.5d0 * tanh((z_Z-z0)/dz0)**2 / cosh((z_Z-z0)/dz0)**2 / dz0 & & * sin(PI*z_Z/z1) & & + 0.5D0 * (1 - tanh((z_Z-z0)/dz0)**3) & & * cos(PI*z_Z/z1) * (PI/z1) do k=1,kmax yz_Temp(1,k) = 0.0D0 do j=1,jmax-1 yz_Temp(j+1,k) = yz_Temp(j,k) dlat =(y_Lat(j+1)-y_Lat(j))/(nmax+1) do n=0,nmax lat = y_Lat(j) + dlat * (n+0.5)/(nmax+1) if ( lat > 0.0d0 ) then yz_Temp(j+1,k) = yz_Temp(j+1,k) & -H/GasRDry & * ( RPlanet * 2*Omega*sin(lat) & + 2 * u0*sin(PI*sin(lat)**2)**3*z_F(k) * tan(lat) ) & * u0*sin(PI*sin(lat)**2)**3*z_dFdz(k) * dlat else yz_Temp(j+1,k) = yz_Temp(j+1,k) endif enddo enddo enddo z_Temp = a_AvrLat_ya(yz_Temp) do k=1,kmax do j=1,jmax do i=0,imax-1 xyz_Temp(i,j,k) = yz_Temp(j,k)-z_Temp(k) + TempUSStandard(z_Z(k)) enddo end do end do ! ! adding temperature disturbance ! do k=1,kmax do j=1,jmax do i=0,imax-1 if ( x_Lon(i) <= PI ) then xyz_Temp(i,j,k) = xyz_Temp(i,j,k) & + That * cosh((x_Lon(i)-lambda0)/alpha)**(-2) & * cosh((y_Lat(j)-phi0)/beta)**(-2) else xyz_Temp(i,j,k) = xyz_Temp(i,j,k) & + That * cosh((x_Lon(i)-lambda0-2*PI)/alpha)**(-2) & * cosh((y_Lat(j)-phi0)/beta)**(-2) endif enddo end do end do end subroutine SetInitData function TempUSStandard(z) ! ! アメリカ標準大気(COESA 1976)温度の計算 ! ! Compute U.S. Standard Atmosphere temperature ! use dc_message, only: MessageNotify real(DP), intent(IN) :: z ! $ z $ . 高度 Height real(DP) :: TempUSStandard ! 米標準大気温度 U.S. Standard Atmosphere temperature if ( z >= 0 .AND. z < 11.0d3 ) then TempUSStandard = 288.15 - 6.5d-3*z else if ( z >= 11.0d3 .AND. z < 20.0d3 ) then TempUSStandard = 288.15 - 6.5d-3*11.0d3 else if ( z >= 20.0d3 .AND. z < 32.0d3 ) then TempUSStandard = 288.15 - 6.5d-3*11.0d3 + 1.0d-3*(z-20.0d3) else if ( z >= 32.0d3 .AND. z < 47.0d3 ) then TempUSStandard = 288.15 - 6.5d-3*11.0d3 + 1.0d-3*(32.0d3-20.0d3) & & + 2.8d-3 * (z-32.0d3) else if ( z >= 47.0d3 .AND. z < 51.0d3 ) then TempUSStandard = 288.15 - 6.5d-3*11.0d3 + 1.0d-3*(32.0d3-20.0d3) & & + 2.8d-3 * (47.0d3-32.0d3) else if ( z >= 51.0d3 .AND. z < 71.0d3 ) then TempUSStandard = 288.15 - 6.5d-3*11.0d3 + 1.0d-3*(32.0d3-20.0d3) & & + 2.8d-3 * (47.0d3-32.0d3) - 2.8d3 * (z-51.0d3) else if ( z >= 71.0d3 .AND. z < 80.0d3 ) then TempUSStandard = 288.15 - 6.5d-3*11.0d3 + 1.0d-3*(32.0d3-20.0d3) & & + 2.8d-3 * (47.0d3-32.0d3) - 2.8d3 * (71.0d3-51.0d3) & & - 2.0d-3 * (z-71.0d3) else if ( z >= 80.0d3 ) then TempUSStandard = 288.15 - 6.5d-3*11.0d3 + 1.0d-3*(32.0d3-20.0d3) & & + 2.8d-3 * (47.0d3-32.0d3) - 2.8d3 * (71.0d3-51.0d3) & & - 2.0d-3 * (80.0d3-71.0d3) else call MessageNotify('E','TempUSStandard','z must be larger than zero') endif end function TempUSStandard end program init_data_primitive