!c Description: !c 大気科学の諸々の物理量の計算 !c !c Current Code Owner: !c sugiyama@gfd-dennou.org !c !c Histry: !c Version Date Comment !c ------- ---------- -------- !c 1.0 2003-06-02 新規作成 !c !c Copyright (C) SUGIYAMA Ko-ichiro, 2003, All rights reserved program OboroAtmos use netcdf use gt4_history use prm_common, only: obrinit, property, AdbtFile, AtmosFile, & & grav, smax, emax use prm_phys, only: gasr use interface_netcdf_open use interface_netcdf_def use interface_netcdf_redef use interface_netcdf_get use interface_netcdf_put use interface_netcdf_close use interface_molweight use interface_cp_noblegas use interface_cp_spl use interface_lapserate use interface_clouddensity implicit none real(8), allocatable :: pres(:) !圧力 real(8), allocatable :: entr(:) !エントロピー real(8), allocatable :: temp(:) !温度 real(8), allocatable :: vtemp(:) !仮温度 real(8), allocatable :: mol(:,:) !化学種のモル数 real(8), allocatable :: emol(:,:) !元素のモル数 real(8), allocatable :: molfr(:,:) !化学種のモル比 real(8), allocatable :: molwt(:) !分子量 real(8), allocatable :: molwt_dry(:) !乾燥大気の分子量 real(8), allocatable :: dtdz(:) !温度減率 real(8), allocatable :: dtdz_low(:) ! 凝結成分少ない近似 real(8), allocatable :: dtdz_high(:) ! 凝結成分多い近似 real(8), allocatable :: dtdz_dry(:) ! 乾燥大気 real(8), allocatable :: vdtdz(:) !仮温度減率 real(8), allocatable :: vdtdz_low(:) ! 凝結成分少ない近似 real(8), allocatable :: vdtdz_high(:) ! 凝結成分多い近似 real(8), allocatable :: vdtdz_dry(:) ! 乾燥大気 real(8), allocatable :: stab(:) !静的安定度 real(8), allocatable :: vstab(:) !仮温度での静的安定度 real(8), allocatable :: vstab_low(:) ! 凝結成分少ない近似 real(8), allocatable :: vstab_high(:) ! 凝結成分多い近似 real(8), allocatable :: cp(:) !比熱 real(8), allocatable :: cp_dry(:) !乾燥大気の比熱 real(8), allocatable :: cp_tmp(:,:) !比熱(作業配列) real(8), allocatable :: cldns(:,:) !雲密度 real(8), allocatable :: mol_tmp(:,:) real(8), allocatable :: molfr_tmp(:,:) real(8), allocatable :: emol_tmp(:,:) real(8), allocatable :: A(:) real(8), allocatable :: B(:) integer :: smax_gas !気相の化学種の数 integer :: num !番号 integer :: snum !化学種番号 integer :: step !計算ステップ数 integer :: ncid character(30) :: long_name !名前 character(10) :: name !名前 character(10) :: dim_name1 character(10) :: dim_name2 character(10) :: units !名前 real(8), parameter :: LatentHeat = 40.66d3 integer, parameter :: cond = 3 !!! !!!データ呼び出し !!! call obrinit !!! !!!NetCDF ファイルのオープン !!!次元は pres, spc, elm であることを知っている !!! !ファイルのオープン call netcdf_open(trim(AdbtFile), ncid) !変数の取得. 取り出す変数は pres, spc, elm, temp, mol, molfr name = 'pres' call netcdf_get(ncid, name, pres) name = 'temp' call netcdf_get(ncid, name, temp) name = 'entr' call netcdf_get(ncid, name, entr) name = 'mol' call netcdf_get(ncid, name, mol_tmp) name = 'molfr' call netcdf_get(ncid, name, molfr_tmp) name = 'emol' call netcdf_get(ncid, name, emol_tmp) !配列の割り当て step = size(pres) allocate(mol(step, smax), & & molfr(step, smax), & & emol(step, emax)) mol = transpose(mol_tmp) molfr = transpose(molfr_tmp) emol = transpose(emol_tmp) !ファイルを閉じる call netcdf_close(ncid) !!! !!!配列の割り当て !!! smax_gas = count(property%phase == 1) allocate(molwt(step), & & molwt_dry(step), & & vtemp(step), & & dtdz(step), & & dtdz_low(step), & & dtdz_high(step), & & dtdz_dry(step), & & vdtdz(step), & & vdtdz_low(step), & & vdtdz_high(step), & & vdtdz_dry(step), & & stab(step), & & vstab(step), & & vstab_low(step), & & vstab_high(step), & & cp(step), & & cp_dry(step), & & cp_tmp(step, smax_gas), & & A(step), B(step)) !!! !!!分子量の計算. 乾燥大気の分子量は大気最上層での分子量とする. !!! !気塊の分子量 do num = 1, step call molweight(molfr(num, :), property, molwt(num)) end do !乾燥大気の分子量 molwt_dry = molwt(step) !!! !!!比熱の計算. 乾燥大気の比熱は大気最上層での分子量とする. !!! !気塊の比熱 do num = 1, step do snum = 1, smax if (property(snum)%method == 1) then call cp_noblegas(temp(num), cp=cp_tmp(num, snum)) elseif (property(snum)%method == 2) then call cp_spl(temp(num), property(snum), cp=cp_tmp(num, snum)) end if end do cp(num) = dot_product(molfr(num, 1: smax_gas), cp_tmp(num,:)) end do !乾燥大気の比熱. 比熱は温度(圧力)依存だが, モル比は一定. do num = 1, step cp_dry(num) = & & dot_product(molfr(step, 1: smax_gas), cp_tmp(num, 1: smax_gas)) end do !!! !!!仮温度の計算 !!! vtemp = temp * (molwt_dry / molwt) !!! !!!温度減率の計算 !!! !温度の減率 call lapserate(temp, pres, grav, molwt, dtdz) !仮温度の減率 call lapserate(vtemp, pres, grav, molwt_dry, vdtdz) !係数. 温度 Temp で評価 A(:) = LatentHeat * Molfr(:, cond) / ( GasR * Temp(:) ) B(:) = ( LatentHeat ** 2.0d0 ) * MolFr(:, cond) & & / ( cp(:) * GasR * ( Temp(:) ** 2.0d0 ) ) write(*,*) Molfr(:, cond) write(*,*) "" write(*,*) A write(*,*) "" write(*,*) B !温度減率 ! dtdz = Grav * MolWt * ( 1 + A ) / ( Cp * ( 1 + B ) ) dtdz_low = Grav * MolWt * ( 1 + A - B ) / Cp dtdz_high = Grav * MolWt * A / ( Cp * B ) !仮温度減率 vdtdz_low = dtdz_low * MolWt_dry / MolWt vdtdz_high = dtdz_high * MolWt_dry / MolWt !!! !!!乾燥大気の断熱温度減率 !!! !温度減率 dtdz_dry = grav * molwt / cp !仮温度減率 vdtdz_dry = grav * molwt_dry / cp_dry !!! !!!静的安定度 !!! !静的安定度 stab = - (grav * (dtdz - dtdz_dry) * 1.0D-3) / temp !仮温度での静的安定度 vstab = - (grav * (vdtdz - vdtdz_dry) * 1.0D-3) / vtemp !仮温度での静的安定度 vstab_low = - (grav * (vdtdz_low - vdtdz_dry) * 1.0D-3) / vtemp !仮温度での静的安定度 vstab_high = - (grav * (vdtdz_high - vdtdz_dry) * 1.0D-3) / vtemp !!! !!!雲量 !!! call clouddensity(temp, pres, property, molfr, cldns) !!! !!!NetCDF ファイルの用意 !!! call output_gtool4_init !!! !!!物理量の書き出し !!! call HistoryPut('entr', entr) call HistoryPut('temp', temp) call HistoryPut('vtemp', vtemp) call HistoryPut('mol', mol) call HistoryPut('emol', emol) call HistoryPut('molfr', molfr) call HistoryPut('molwt', molwt) call HistoryPut('molwt_dry', molwt_dry) call HistoryPut('dtdz', dtdz) call HistoryPut('dtdz_low', dtdz_low) call HistoryPut('dtdz_high', dtdz_high) call HistoryPut('dtdz_dry', dtdz_dry) call HistoryPut('vdtdz', vdtdz) call HistoryPut('vdtdz_low', vdtdz_low) call HistoryPut('vdtdz_high', vdtdz_high) call HistoryPut('vdtdz_dry', vdtdz_dry) call HistoryPut('stab', stab) call HistoryPut('vstab', vstab) call HistoryPut('vstab_low', vstab_low) call HistoryPut('vstab_high', vstab_high) call HistoryPut('cp', cp) call HistoryPut('cldns', cldns) !!! !!!NetCDF ファイルを閉じる !!! call HistoryClose !!! !!!名前を入れる !!! call netcdf_open(AtmosFile, ncid) name = "sname" dim_name1 = "size" dim_name2 = "spc" long_name = "name of chemical species" units = "1" call netcdf_redef(ncid, name, dim_name1, dim_name2, & & long_name, units, nf90_char) call netcdf_put(ncid, name, property%name) call netcdf_close(ncid) contains subroutine output_gtool4_init use prm_common, only: planet, AtmosFile implicit none real(8) :: spc(smax) real(8) :: elm(emax) real(8) :: sizes(20) integer :: num integer :: i character(1000) :: dataset !ファイルに書き出すべき初期値を変数に入れる write(dataset, *) planet, "; pres_ini=", pres(1), & & "; pres_end=", pres(step), & & "; temp_ini=", temp(1), & & "; step=", step num = len_trim(dataset) !大域属性, 座標軸の設定 !"size" は化学種名を入れるためのもの call HistoryCreate(& & file=AtmosFile, & & title="Equilibrium Calculation", & & source=dataset(:num+1), & & institution="GFD_Dennou_Club oboro Project", & & dims=(/'pres', 'spc ', 'elm ', 'size'/), & & dimsizes=(/step, smax, emax, 20/), & & longnames=(/"Pressure ", "chemical species", & & "element ", "size of name "/), & & units=(/"Pa", "1 ", "1 ", "1 "/), & & origin=0.0, & & interval=0.0) !次元の定義 spc = (/(i, i = 1, smax)/) call HistoryPut('spc', real(spc, 4)) elm = (/(i, i = 1, emax)/) call HistoryPut('elm', real(elm, 4)) call HistoryPut('pres', real(pres, 4)) sizes = (/(i, i = 1, 20)/) call HistoryPut('size', real(sizes, 4)) !変数定義, エントロピー call HistoryAddVariable( & & varname="entr", & & dims=(/'pres'/), & & longname="Entropy", & & units='J', & & xtype="double" ) !変数定義, 温度 call HistoryAddVariable( & & varname="temp", & & dims=(/'pres'/), & & longname="Temperature", & & units='K', & & xtype="double" ) !変数定義, 仮温度 call HistoryAddVariable( & & varname="vtemp", & & dims=(/'pres'/), & & longname="Vertial Temperature", & & units='K', & & xtype="double" ) !変数定義, モル call HistoryAddVariable( & & varname="mol", & & dims=(/'pres', 'spc '/), & & longname="Chemical Species Abundance", & & units='mol', & & xtype="double" ) !変数定義, 元素量 call HistoryAddVariable( & & varname="emol", & & dims=(/'pres', 'elm '/), & & longname="Element Abundance", & & units='mol', & & xtype="double" ) !変数定義, モル比 call HistoryAddVariable( & & varname="molfr", & & dims=(/'pres', 'spc '/), & & longname="Mol Fraction", & & units='1', & & xtype="double" ) !変数定義, 大気の分子量 call HistoryAddVariable( & & varname="molwt", & & dims=(/'pres'/), & & longname="Atmospheric Mol Weight", & & units='g/mol', & & xtype="double" ) !変数定義, 乾燥大気の分子量 call HistoryAddVariable( & & varname="molwt_dry", & & dims=(/'pres'/), & & longname="Atmospheric Mol Weight", & & units='g/mol', & & xtype="double" ) !変数定義, 温度減率 call HistoryAddVariable( & & varname="dtdz", & & dims=(/'pres'/), & & longname="Lapse Rate", & & units='K/km', & & xtype="double" ) !変数定義, 温度減率 call HistoryAddVariable( & & varname="dtdz_low", & & dims=(/'pres'/), & & longname="Lapse Rate", & & units='K/km', & & xtype="double" ) !変数定義, 温度減率 call HistoryAddVariable( & & varname="dtdz_high", & & dims=(/'pres'/), & & longname="Lapse Rate", & & units='K/km', & & xtype="double" ) !変数定義, 乾燥大気の温度減率 call HistoryAddVariable( & & varname="dtdz_dry", & & dims=(/'pres'/), & & longname="Dry Lapse Rate", & & units='K/km', & & xtype="double" ) !変数定義, 仮温度減率 call HistoryAddVariable( & & varname="vdtdz", & & dims=(/'pres'/), & & longname="Lapse Rate", & & units='K/km', & & xtype="double" ) !変数定義, 仮温度減率 call HistoryAddVariable( & & varname="vdtdz_low", & & dims=(/'pres'/), & & longname="Lapse Rate", & & units='K/km', & & xtype="double" ) !変数定義, 仮温度減率 call HistoryAddVariable( & & varname="vdtdz_high", & & dims=(/'pres'/), & & longname="Lapse Rate", & & units='K/km', & & xtype="double" ) !変数定義, 乾燥大気の仮温度減率 call HistoryAddVariable( & & varname="vdtdz_dry", & & dims=(/'pres'/), & & longname="Dry Lapse Rate", & & units='K/km', & & xtype="double" ) !変数定義, 静的安定度 call HistoryAddVariable( & & varname="stab", & & dims=(/'pres'/), & & longname="Stability", & & units='s^-2', & & xtype="double" ) !変数定義, 静的安定度 call HistoryAddVariable( & & varname="vstab", & & dims=(/'pres'/), & & longname="Stability", & & units='s^-2', & & xtype="double" ) !変数定義, 静的安定度 call HistoryAddVariable( & & varname="vstab_low", & & dims=(/'pres'/), & & longname="Stability", & & units='s^-2', & & xtype="double" ) !変数定義, 静的安定度 call HistoryAddVariable( & & varname="vstab_high", & & dims=(/'pres'/), & & longname="Stability", & & units='s^-2', & & xtype="double" ) !変数定義, 比熱 call HistoryAddVariable( & & varname="cp", & & dims=(/'pres'/), & & longname="Specific Heat", & & units='J/ k mol', & & xtype="double" ) !変数定義, 雲密度 call HistoryAddVariable( & & varname="cldns", & & dims=(/'pres', 'spc '/), & & longname="Cloud Density", & & units='g/l', & & xtype="double" ) end subroutine output_gtool4_init end program OboroAtmos