Class | initialdata_basic |
In: |
env/initialdata_basic.f90
|
デフォルトの基本場を作るためのルーチン
Constant : | |
ID_FlaginitDataBasic_isothermal = 3 : | integer, parameter, public |
Subroutine : | |
z_TempBZ(kmin:kmax) : | real(8), intent(out) |
z_PressBZ(kmin:kmax) : | real(8), intent(out) |
za_MolFr(kmin:kmax, 1:ncmax) : | real(8), intent(out) |
subroutine initialdata_basic_dry( z_TempBZ, z_PressBZ, za_MolFr ) implicit none real(8), intent(out) :: z_TempBZ(kmin:kmax) real(8), intent(out) :: z_PressBZ(kmin:kmax) real(8), intent(out) :: za_MolFr(kmin:kmax, 1:ncmax) call eccm_dry( SpcWetMolFr(1:ncmax), Humidity, z_TempBZ, z_PressBZ, za_MolFr ) end subroutine initialdata_basic_dry
Subroutine : | |
z_TempBZ(kmin:kmax) : | real(8), intent(out) |
z_PressBZ(kmin:kmax) : | real(8), intent(out) |
za_MolFr(kmin:kmax, 1:ncmax) : | real(8), intent(out) |
subroutine initialdata_basic_isothermal( z_TempBZ, z_PressBZ, za_MolFr ) implicit none real(8), intent(out) :: z_TempBZ(kmin:kmax) real(8), intent(out) :: z_PressBZ(kmin:kmax) real(8), intent(out) :: za_MolFr(kmin:kmax, 1:ncmax) z_TempBZ = TempSfc z_PressBZ = PressSfc za_MolFr = 0.0d0 end subroutine initialdata_basic_isothermal
Subroutine : | |
z_TempBZ(kmin:kmax) : | real(8), intent(out) |
z_PressBZ(kmin:kmax) : | real(8), intent(out) |
za_MolFr(kmin:kmax, 1:ncmax) : | real(8), intent(out) |
subroutine initialdata_basic_moist( z_TempBZ, z_PressBZ, za_MolFr ) implicit none real(8), intent(out) :: z_TempBZ(kmin:kmax) real(8), intent(out) :: z_PressBZ(kmin:kmax) real(8), intent(out) :: za_MolFr(kmin:kmax, 1:ncmax) call eccm_wet( SpcWetMolFr(1:ncmax), Humidity, z_TempBZ, z_PressBZ, za_MolFr ) end subroutine initialdata_basic_moist
Subroutine : | |
cfgfile : | character(STRING), intent(in) |
nml の読み込み
This procedure input/output NAMELIST#initialdata_basic .
subroutine initialdata_basic_nml(cfgfile) ! nml の読み込み ! use dc_iounit, only : FileOpen implicit none character(STRING), intent(in) :: cfgfile integer :: unit !装置番号 character(STRING) :: FlaginitdataBasic = "" NAMELIST /initialdata_basic/ FlaginitDataBasic, Humidity, TempStrat, HeightStrat, Dhight call FileOpen(unit, file=cfgfile, mode='r') read(unit, NML=initialdata_basic) close(unit) if (FlaginitdataBasic == "moist") then ID_FlaginitdataBasic = ID_FlaginitdataBasic_moist elseif (FlaginitdataBasic == "dry") then ID_FlaginitdataBasic = ID_FlaginitdataBasic_dry elseif (FlaginitdataBasic == "isothermal") then ID_FlaginitdataBasic = ID_FlaginitdataBasic_isothermal end if end subroutine initialdata_basic_nml
Subroutine : | |
z_TempBZ(kmin:kmax) : | real(8), intent(inout) |
z_PressBZ(kmin:kmax) : | real(8), intent(inout) |
subroutine initialdata_basic_strat(z_TempBZ, z_PressBZ) ! implicit none real(8), intent(inout) :: z_TempBZ(kmin:kmax) real(8), intent(inout) :: z_PressBZ(kmin:kmax) real(8) :: z_DTempDZ(kmin:kmax) real(8) :: DTempDZ real(8) :: Weight integer :: k do k = 1, kmax-1 z_DTempDZ(k) = (z_TempBZ(k) - z_TempBZ(k-1)) / r_dz(k-1) end do ! 対流圏界面より上の扱い ! 圏界面より上は等温大気とする. 等温位大気から等温大気への遷移は ! tanh を用いてなめらかにつなぐ do k = 1, kmax !重みつけの関数を用意. tanh を用いる Weight = ( tanh( (z_Z(k) - HeightStrat ) / Dhight ) + 1.0d0 ) * 5.0d-1 !仮値として温度を計算する. 圏界面より上では TempStrat の等温大気に近づける z_TempBZ(k) = z_TempBZ(k) * ( 1.0d0 - Weight ) + TempStrat * Weight !温度減率が断熱温度減率より小さくならないように DTempDZ = max( z_DTempDZ(k), (z_TempBZ(k) - z_TempBZ(k-1)) / r_dz(k-1) ) !基本場の温度を決める z_TempBZ(k) = z_TempBZ(k-1) + DTempDZ * r_dz(k-1) !圧力を静水圧平衡から計算 z_PressBZ(k) = z_PressBZ(k-1) * ( ( z_TempBZ(k-1) / z_TempBZ(k) ) ** (Grav / ( DTempDZ * GasRDry ) ) ) end do end subroutine initialdata_basic_strat