Class | WarmRainPrm |
In: |
moist/warmrainprm.f90
|
Copyright (C) GFD Dennou Club, 2005, 2006. All rights reserved.
* Developer: SUGIYAMA Ko-ichiro * Version: $Id: warmrainprm.f90,v 1.1.1.1 2006/04/25 03:43:58 deepconv Exp $ * Tag Name: $Name: $ * Change History:
暖かい雨のバルク法を用いた, 水蒸気と雨, 雲と雨の混合比の変換係数を求める.
* 中島健介 (1994) で利用した定式をそのまま利用.
subroutine WarmRainPrm_Init () !暗黙の型宣言禁止 implicit none !変数定義 integer :: s integer :: n1, n2, LoopNum2 !----------------------------------------------------------- ! 雨粒と雲粒と気体の ID の組を作る !----------------------------------------------------------- !初期化 allocate( RainSW(SpcNum) ) LoopNum = 0 LoopNum2 = 0 RainSW = 0.0d0 !化学種の中から雨粒を作るものを選び, その配列添え字と分子量を保管. SelectCloud: do s = 1, SpcNum !'Cloud' という文字列が含まれるものの個数を数える n1 = index(SpcWetSymbol(s), '-Cloud' ) if (n1 /= 0) then LoopNum = LoopNum + 1 GasNum(LoopNum) = minloc(SpcWetID, 1, SpcWetID == ChemData_OneSpcID(SpcWetSymbol(s)(1:n1-3) // '-g')) CloudNum(LoopNum) = s end if !'Rain' という文字列が含まれるものの個数を数える n2 = index(SpcWetSymbol(s), '-Rain' ) if (n2 /= 0) then LoopNum2 = LoopNum2 + 1 RainNum(LoopNum2) = s RainSW(s) = 1.0d0 end if ! NH4SH が存在する場合は LoopNum を 1 つ減らす if ( trim(SpcWetSymbol(s)) == 'NH4SH-s-Cloud' ) then LoopNum = LoopNum - 1 end if end do SelectCloud !----------------------------------------------------------- ! 硫化アンモニウム, およびアンモニアと硫化水素の ID を取得 ! 'Cloud' の方が 'Rain' よりも先に存在すると仮定している !----------------------------------------------------------- NH3Num = minloc(SpcWetID, 1, SpcWetID == ChemData_OneSpcID('NH3-g')) H2SNum = minloc(SpcWetID, 1, SpcWetID == ChemData_OneSpcID('H2S-g')) NH4SHCloudNum = minloc(SpcWetID, 1, SpcWetID == ChemData_OneSpcID('NH4SH-s')) NH4SHRainNum = NH4SHCloudNum + 1 !----------------------------------------------------------- ! 確認 !----------------------------------------------------------- if ( LoopNum == 0 ) then write(*,*) "WarmRainPrm: CloudNum = 0, please comment out of WarmRainPrm" ! stop end if write(*,*) "WarmRainPrm_Init, LoopNum: ", LoopNum write(*,*) "WarmRainPrm_Init, GasNum: ", GasNum write(*,*) "WarmRainPrm_Init, CloudNum: ", CloudNum write(*,*) "WarmRainPrm_Init, RainNum: ", RainNum write(*,*) "WarmRainPrm_Init, RainSW: ", RainSW write(*,*) "WarmRainPrm_Init, NH3Num: ", NH3Num write(*,*) "WarmRainPrm_Init, H2SNum: ", H2SNum write(*,*) "WarmRainPrm_Init, NH4SHNum: ", NH4SHCloudNum write(*,*) "WarmRainPrm_Init, NH4SHNum: ", NH4SHRainNum end subroutine