| Class | sosi_utils |
| In: |
sosi/sosi_utils.f90
|
Note that Japanese and English are described in parallel.
| !$ ! PhyImplSDHTendency : | 時間変化率の計算 |
| !$ ! PhyImplSDHSetMethodFromMatthews : | SurfType から計算法インデクスの作成 |
| !$ ! PhyImplSDHInit : | 初期化 |
| !$ ! ——————————- : | ———— |
| !$ ! PhyImplSDHTendency : | Calculate tendency |
| !$ ! PhyImplSDHSetMethodFromMatthews : | Set index for calculation method |
| !$ ! PhyImplSDHInit : | Initialization |
| Subroutine : | |||
| xy_SOSeaIceMass(0:imax-1, 1:jmax) : | real(DP), intent(inout)
| ||
| xyz_SOSeaIceTemp(0:imax-1, 1:jmax, 1:ksimax) : | real(DP), intent(inout) | ||
| xy_DSOSeaIceMassDtPhyTop(0:imax-1, 1:jmax) : | real(DP), intent(in ) | ||
| xy_DSOSeaIceMassDtPhyBot(0:imax-1, 1:jmax) : | real(DP), intent(in ) | ||
| xyz_DSOSeaIceTempDtPhy(0:imax-1, 1:jmax, 1:ksimax) : | real(DP), intent(in ) |
Calculates slab sea ice horizontal transports by diffusion
subroutine SOSIUtilsAddPhysics( xy_SOSeaIceMass, xyz_SOSeaIceTemp, xy_DSOSeaIceMassDtPhyTop, xy_DSOSeaIceMassDtPhyBot, xyz_DSOSeaIceTempDtPhy )
!
! Calculates slab sea ice horizontal transports by diffusion
! ヒストリデータ出力
! History data output
!
use gtool_historyauto, only: HistoryAutoPut
use timeset , only : TimeN, DelTime
! $\Delta t$
! 雪と海氷の定数の設定
! Setting constants of snow and sea ice
!
use constants_snowseaice, only: TempBelowSeaIce, SeaIceDen
! 宣言文 ; Declaration statements
!
real(DP), intent(inout) :: xy_SOSeaIceMass (0:imax-1, 1:jmax)
! $ M_si $ . 海氷質量 (kg m-2)
! Slab ocean sea ice mass (kg m-2)
real(DP), intent(inout) :: xyz_SOSeaIceTemp (0:imax-1, 1:jmax, 1:ksimax)
real(DP), intent(in ) :: xy_DSOSeaIceMassDtPhyTop(0:imax-1, 1:jmax)
real(DP), intent(in ) :: xy_DSOSeaIceMassDtPhyBot(0:imax-1, 1:jmax)
real(DP), intent(in ) :: xyz_DSOSeaIceTempDtPhy(0:imax-1, 1:jmax, 1:ksimax)
! 作業変数
! Work variables
!
real(DP) :: xy_SeaIceThickness(0:imax-1, 1:jmax)
integer :: xy_SOSILocalKMax (0:imax-1, 1:jmax)
real(DP) :: xyr_SOSILocalDepth(0:imax-1, 1:jmax, 0:ksimax)
real(DP) :: xyz_SOSILocalDepth(0:imax-1, 1:jmax, 1:ksimax)
real(DP) :: xy_SOSeaIceMassTent1 (0:imax-1, 1:jmax)
real(DP) :: xy_SeaIceThicknessTent1(0:imax-1, 1:jmax)
integer :: xy_SOSILocalKMaxTent1 (0:imax-1, 1:jmax)
real(DP) :: xyr_SOSILocalDepthTent1(0:imax-1, 1:jmax, 0:ksimax)
real(DP) :: xyz_SOSILocalDepthTent1(0:imax-1, 1:jmax, 1:ksimax)
real(DP) :: xy_SOSeaIceMassTent2 (0:imax-1, 1:jmax)
real(DP) :: xy_SeaIceThicknessTent2(0:imax-1, 1:jmax)
integer :: xy_SOSILocalKMaxTent2 (0:imax-1, 1:jmax)
real(DP) :: xyr_SOSILocalDepthTent2(0:imax-1, 1:jmax, 0:ksimax)
real(DP) :: xyz_SOSILocalDepthTent2(0:imax-1, 1:jmax, 1:ksimax)
real(DP) :: xyr_SOSILocalDepthT2MapToT1(0:imax-1, 1:jmax, 0:ksimax)
real(DP) :: xyz_SOSeaIceTempTent1(0:imax-1, 1:jmax, 1:ksimax)
real(DP) :: xyz_SOSeaIceTempTent2(0:imax-1, 1:jmax, 1:ksimax)
!!$
!!$ real(DP) :: xyz_DSOSeaIceTempDtPhyUpdt(0:imax-1, 1:jmax, 1:ksimax)
integer:: i ! 東西方向に回る DO ループ用作業変数
! Work variables for DO loop in zonal direction
integer:: j ! 南北方向に回る DO ループ用作業変数
! Work variables for DO loop in meridional direction
integer:: k
integer:: kk
integer:: kTop
integer:: kBot
if ( .not. sosi_utils_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
! Calculate sea ice mass at next time step
!
! Add sea ice mass change at the bottom of sea ice
xy_SOSeaIceMassTent1 = xy_SOSeaIceMass + xy_DSOSeaIceMassDtPhyBot * DelTime
! Add sea ice mass change at the top of sea ice
xy_SOSeaIceMassTent2 = xy_SOSeaIceMassTent1 + xy_DSOSeaIceMassDtPhyTop * DelTime
! 海氷温度時間積分
! Time integration of sea ice temperature
!
xyz_SOSeaIceTempTent1 = xyz_SOSeaIceTemp + xyz_DSOSeaIceTempDtPhy * DelTime
!
! Adjust temperature
!
!
! Calcuate sea ice thickness
!
xy_SeaIceThickness = xy_SOSeaIceMass / SeaIceDen
!
! Set slab ocean sea ice levels
!
call SOSIUtilsSetSOSeaIceLevels( xy_SeaIceThickness, xy_SOSILocalKMax, xyr_SOSILocalDepth, xyz_SOSILocalDepth )
!
! Calcuate sea ice thickness
!
xy_SeaIceThicknessTent1 = xy_SOSeaIceMassTent1 / SeaIceDen
!
! Set slab ocean sea ice levels
!
call SOSIUtilsSetSOSeaIceLevels( xy_SeaIceThicknessTent1, xy_SOSILocalKMaxTent1, xyr_SOSILocalDepthTent1, xyz_SOSILocalDepthTent1 )
! Adjust levels
!
do j = 1, jmax
do i = 0, imax-1
if ( xy_SOSILocalKMaxTent1(i,j) > xy_SOSILocalKMax(i,j) ) then
! sea ice thickness increases
if ( xy_SOSILocalKMax(i,j) == 0 ) then
do k = 1, xy_SOSILocalKMaxTent1(i,j)
xyz_SOSeaIceTempTent1(i,j,k) = TempBelowSeaIce
end do
else
!!$ do k = 1, xy_SOSILocalKMaxB(i,j)
!!$ ! Do nothing
!!$ end do
do k = xy_SOSILocalKMax(i,j)+1, xy_SOSILocalKMaxTent1(i,j)
kk = xy_SOSILocalKMax(i,j)
xyz_SOSeaIceTempTent1(i,j,k) = xyz_SOSeaIceTempTent1(i,j,kk)
end do
end if
else if ( xy_SOSILocalKMaxTent1(i,j) < xy_SOSILocalKMax(i,j) ) then
! sea ice thickness decreases
! Do nothing
! Melted sea ice had freezing temperature
end if
end do
end do
!
! Calcuate sea ice thickness
!
xy_SeaIceThicknessTent2 = xy_SOSeaIceMassTent2 / SeaIceDen
!
! Set slab ocean sea ice levels
!
call SOSIUtilsSetSOSeaIceLevels( xy_SeaIceThicknessTent2, xy_SOSILocalKMaxTent2, xyr_SOSILocalDepthTent2, xyz_SOSILocalDepthTent2 )
! Adjust levels
!
do j = 1, jmax
do i = 0, imax-1
if ( xy_SOSILocalKMaxTent2(i,j) == 0 ) then
do k = 0, ksimax
xyr_SOSILocalDepthT2MapToT1(i,j,k) = 0.0_DP
end do
else
do k = 0, xy_SOSILocalKMaxTent2(i,j)
xyr_SOSILocalDepthT2MapToT1(i,j,k) = xyr_SOSILocalDepthTent2(i,j,k) - ( xy_SeaIceThicknessTent1(i,j) - xy_SeaIceThicknessTent2(i,j) )
end do
do k = xy_SOSILocalKMaxTent2(i,j)+1, ksimax
xyr_SOSILocalDepthT2MapToT1(i,j,k) = 0.0_DP
end do
end if
end do
end do
!
do j = 1, jmax
do i = 0, imax-1
if ( xy_SOSILocalKMaxTent2(i,j) == 0 ) then
do k = 1, ksimax
xyz_SOSeaIceTempTent2(i,j,k) = SOSeaIceTempMissingValue
end do
else
do k = 1, xy_SOSILocalKMaxTent2(i,j)
search_kTop : do kTop = 1, xy_SOSILocalKMaxTent1(i,j)
if ( xyr_SOSILocalDepthT2MapToT1(i,j,k-1) >= xyr_SOSILocalDepthTent1(i,j,kTop) ) then
exit search_kTop
end if
end do search_kTop
search_kBot : do kBot = kTop, xy_SOSILocalKMaxTent1(i,j)
if ( xyr_SOSILocalDepthT2MapToT1(i,j,k ) >= xyr_SOSILocalDepthTent1(i,j,kBot) ) then
exit search_kBot
end if
end do search_kBot
if ( kTop == kBot ) then
kk = kTop
xyz_SOSeaIceTempTent2(i,j,k) = + xyz_SOSeaIceTempTent1(i,j,kk) * ( xyr_SOSILocalDepthT2MapToT1(i,j,k-1) - xyr_SOSILocalDepthT2MapToT1(i,j,k ) )
else
xyz_SOSeaIceTempTent2(i,j,k) = 0.0_DP
kk = kTop
xyz_SOSeaIceTempTent2(i,j,k) = xyz_SOSeaIceTempTent2(i,j,k) + xyz_SOSeaIceTempTent1(i,j,kk) * ( xyr_SOSILocalDepthT2MapToT1(i,j,k-1) - xyr_SOSILocalDepthTent1(i,j,kk) )
do kk = kTop+1, kBot-1
xyz_SOSeaIceTempTent2(i,j,k) = xyz_SOSeaIceTempTent2(i,j,k) + xyz_SOSeaIceTempTent1(i,j,kk) * ( xyr_SOSILocalDepthTent1(i,j,kk-1) - xyr_SOSILocalDepthTent1(i,j,kk ) )
end do
kk = kBot
xyz_SOSeaIceTempTent2(i,j,k) = xyz_SOSeaIceTempTent2(i,j,k) + xyz_SOSeaIceTempTent1(i,j,kk) * ( xyr_SOSILocalDepthTent1(i,j,kk-1) - xyr_SOSILocalDepthT2MapToT1(i,j,k) )
end if
!
xyz_SOSeaIceTempTent2(i,j,k) = xyz_SOSeaIceTempTent2(i,j,k) / ( xyr_SOSILocalDepthTent2(i,j,k-1) - xyr_SOSILocalDepthTent2(i,j,k) )
end do
do k = xy_SOSILocalKMaxTent2(i,j)+1, ksimax
xyz_SOSeaIceTempTent2(i,j,k) = SOSeaIceTempMissingValue
end do
end if
end do
end do
! Update sea ice temperature
!
xyz_SOSeaIceTemp = xyz_SOSeaIceTempTent2
! Update sea ice mass
!
xy_SOSeaIceMass = xy_SOSeaIceMassTent2
! Check sea ice mass
!
do j = 1, jmax
do i = 0, imax-1
if ( xy_SOSeaIceMass(i,j) < 0 ) then
if ( xy_SOSeaIceMass(i,j) < SOSeaIceMassNegativeThreshold ) then
call MessageNotify( 'M', module_name, ' Slab sea ice mass is negative after physics, %f, and this is set to zero.', d = (/ xy_SOSeaIceMass(i,j) /) )
end if
xy_SOSeaIceMass(i,j) = 0.0_DP
end if
end do
end do
! Check
!
xy_SeaIceThickness = xy_SOSeaIceMass / SeaIceDen
!
call SOSIUtilsChkSOSeaIce( xy_SeaIceThickness, xyz_SOSeaIceTemp, "SOSIUtilsAddPhysics" )
end subroutine SOSIUtilsAddPhysics
| Subroutine : | |
| xy_SeaIceThickness(0:imax-1, 1:jmax) : | real(DP), intent(in ) |
| xyz_SOSeaIceTemp(0:imax-1, 1:jmax, 1:ksimax) : | real(DP), intent(in ) |
| ParentRoutine : | character(*), intent(in ), optional |
Set index for calculation method from Matthews’ index
subroutine SOSIUtilsChkSOSeaIce( xy_SeaIceThickness, xyz_SOSeaIceTemp, ParentRoutine )
!
!
!
! Set index for calculation method from Matthews' index
!
! モジュール引用 ; USE statements
!
! 宣言文 ; Declaration statements
!
real(DP), intent(in ) :: xy_SeaIceThickness(0:imax-1, 1:jmax)
real(DP), intent(in ) :: xyz_SOSeaIceTemp (0:imax-1, 1:jmax, 1:ksimax)
character(*), intent(in ), optional :: ParentRoutine
! 作業変数
! Work variables
!
integer :: xy_SOSILocalKMax(0:imax-1, 1:jmax)
real(DP) :: xyr_SOSILocalDepth(0:imax-1, 1:jmax, 0:ksimax)
real(DP) :: xyz_SOSILocalDepth(0:imax-1, 1:jmax, 1:ksimax)
integer:: i ! 経度方向に回る DO ループ用作業変数
! Work variables for DO loop in longitude
integer:: j ! 緯度方向に回る DO ループ用作業変数
! Work variables for DO loop in latitude
integer:: k
! 実行文 ; Executable statement
!
! 初期化確認
! Initialization check
!
if ( .not. sosi_utils_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
call SOSIUtilsSetSOSeaIceLevels( xy_SeaIceThickness, xy_SOSILocalKMax, xyr_SOSILocalDepth, xyz_SOSILocalDepth )
do i = 0, imax-1
do j = 1, jmax
do k = 1, xy_SOSILocalKMax(i,j)
if ( xyz_SOSeaIceTemp(i,j,k) < 0.0_DP ) then
if ( present( ParentRoutine ) ) then
call MessageNotify( 'M', module_name, 'Called from %c:', c1 = trim( ParentRoutine ) )
end if
call MessageNotify( 'M', module_name, 'xyz_SOSeaIceTemp(%d,%d,%d) = %f.', i = (/i,j,k/), d = (/xyz_SOSeaIceTemp(i,j,k)/) )
end if
end do
end do
end do
end subroutine SOSIUtilsChkSOSeaIce
| Subroutine : |
This procedure input/output NAMELIST#sosi_utils_nml .
subroutine SOSIUtilsInit
!
! sosi_utils モジュールの初期化を行います.
! NAMELIST#sosi_utils_nml の読み込みはこの手続きで行われます.
!
! "sosi_utils" module is initialized.
! "NAMELIST#sosi_utils_nml" is loaded in this procedure.
!
! モジュール引用 ; USE statements
!
! NAMELIST ファイル入力に関するユーティリティ
! Utilities for NAMELIST file input
!
use namelist_util, only: namelist_filename, NmlutilMsg, NmlutilAryValid
! ファイル入出力補助
! File I/O support
!
use dc_iounit, only: FileOpen
! 種別型パラメタ
! Kind type parameter
!
use dc_types, only: STDOUT ! 標準出力の装置番号. Unit number of standard output
! 文字列操作
! Character handling
!
use dc_string, only: StoA
! 宣言文 ; Declaration statements
!
! 作業変数
! Work variables
!
integer:: unit_nml ! NAMELIST ファイルオープン用装置番号.
! Unit number for NAMELIST file open
integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT.
! IOSTAT of NAMELIST read
! NAMELIST 変数群
! NAMELIST group name
!
namelist /sosi_utils_nml/ SOSeaIceMassNegativeThreshold
!!$ & SOMass, & ! Slab ocean heat capacity (J m-2 K-1)
!!$ & NumMaxItr, & ! Number of interation
!!$ & TempItrCrit, &
!!$ & FlagSublimation
!
! デフォルト値については初期化手続 "sosi_utils#SOSOUtilsInit"
! のソースコードを参照のこと.
!
! Refer to source codes in the initialization procedure
! "sosi_utils#SOSIUtilsInit" for the default values.
!
! 実行文 ; Executable statement
!
if ( sosi_utils_inited ) return
! デフォルト値の設定
! Default values settings
!
SOSeaIceMassNegativeThreshold = -1.0e-10_DP
! NAMELIST の読み込み
! NAMELIST is input
!
if ( trim(namelist_filename) /= '' ) then
call FileOpen( unit_nml, namelist_filename, mode = 'r' ) ! (in)
rewind( unit_nml )
read( unit_nml, nml = sosi_utils_nml, iostat = iostat_nml ) ! (out)
close( unit_nml )
call NmlutilMsg( iostat_nml, module_name ) ! (in)
end if
! Initialization of modules used in this model
!
! 印字 ; Print
!
call MessageNotify( 'M', module_name, '----- Initialization Messages -----' )
call MessageNotify( 'M', module_name, ' SOSeaIceMassNegativeThreshold = %f', d = (/ SOSeaIceMassNegativeThreshold /) )
!!$ call MessageNotify( 'M', module_name, ' SOMass = %f', d = (/ SOMass /) )
!!$ call MessageNotify( 'M', module_name, ' SOHeatCapacity = %f', d = (/ SOHeatCapacity /) )
!!$ call MessageNotify( 'M', module_name, ' NumMaxItr = %d', i = (/ NumMaxItr /) )
!!$ call MessageNotify( 'M', module_name, ' TempItrCrit = %f', d = (/ TempItrCrit /) )
!!$ call MessageNotify( 'M', module_name, ' FlagSublimation = %b', l = (/ FlagSublimation /) )
call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )
sosi_utils_inited = .true.
end subroutine SOSIUtilsInit
| Subroutine : | |||
| xy_SOSeaIceMass(0:imax-1, 1:jmax) : | real(DP), intent(in )
| ||
| xyz_SOSeaIceTemp(0:imax-1, 1:jmax, 1:ksimax) : | real(DP), intent(inout) | ||
| SOSeaIceValue : | real(DP), intent(in ), optional |
Set missing values
subroutine SOSIUtilsSetMissingValue( xy_SOSeaIceMass, xyz_SOSeaIceTemp, SOSeaIceValue )
!
! Set missing values
! 雪と海氷の定数の設定
! Setting constants of snow and sea ice
!
use constants_snowseaice, only: TempBelowSeaIce, SeaIceDen
! 宣言文 ; Declaration statements
!
real(DP), intent(in ) :: xy_SOSeaIceMass (0:imax-1, 1:jmax)
! $ M_si $ . 海氷質量 (kg m-2)
! Slab ocean sea ice mass (kg m-2)
real(DP), intent(inout) :: xyz_SOSeaIceTemp (0:imax-1, 1:jmax, 1:ksimax)
real(DP), intent(in ), optional :: SOSeaIceValue
! 作業変数
! Work variables
!
real(DP) :: xy_SeaIceThickness(0:imax-1, 1:jmax)
!
! Sea ice thickness
integer :: xy_SOSILocalKMax (0:imax-1, 1:jmax)
real(DP) :: xyr_SOSILocalDepth(0:imax-1, 1:jmax, 0:ksimax)
real(DP) :: xyz_SOSILocalDepth(0:imax-1, 1:jmax, 1:ksimax)
real(DP) :: xy_SeaIceThicknessA(0:imax-1, 1:jmax)
!
! Sea ice thickness
integer :: xy_SOSILocalKMaxA (0:imax-1, 1:jmax)
real(DP) :: xyr_SOSILocalDepthA(0:imax-1, 1:jmax, 0:ksimax)
real(DP) :: xyz_SOSILocalDepthA(0:imax-1, 1:jmax, 1:ksimax)
real(DP) :: SetValue
integer:: i ! 東西方向に回る DO ループ用作業変数
! Work variables for DO loop in zonal direction
integer:: j ! 南北方向に回る DO ループ用作業変数
! Work variables for DO loop in meridional direction
integer:: k
if ( .not. sosi_utils_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
!
! Calcuate sea ice thickness
!
xy_SeaIceThickness = xy_SOSeaIceMass / SeaIceDen
!
! Set slab ocean sea ice levels
!
call SOSIUtilsSetSOSeaIceLevels( xy_SeaIceThickness, xy_SOSILocalKMax, xyr_SOSILocalDepth, xyz_SOSILocalDepth )
if ( present ( SOSeaIceValue ) ) then
SetValue = SOSeaIceValue
else
SetValue = SOSeaIceTempMissingValue
end if
!
! Set missing values
!
do j = 1, jmax
do i = 0, imax-1
do k = xy_SOSILocalKMax(i,j)+1, ksimax
xyz_SOSeaIceTemp(i,j,k) = SetValue
end do
end do
end do
end subroutine SOSIUtilsSetMissingValue
| Subroutine : | |||
| xy_SeaIceThickness(0:imax-1, 1:jmax) : | real(DP), intent(in )
| ||
| xy_SOSILocalKMax(0:imax-1, 1:jmax) : | integer , intent(out) | ||
| xyr_SOSILocalDepth(0:imax-1, 1:jmax, 0:ksimax) : | real(DP), intent(out) | ||
| xyz_SOSILocalDepth(0:imax-1, 1:jmax, 1:ksimax) : | real(DP), intent(out) |
Set index for calculation method from Matthews’ index
subroutine SOSIUtilsSetSOSeaIceLevels( xy_SeaIceThickness, xy_SOSILocalKMax, xyr_SOSILocalDepth, xyz_SOSILocalDepth )
!
!
!
! Set index for calculation method from Matthews' index
!
! モジュール引用 ; USE statements
!
! 座標データ設定
! Axes data settings
!
use axesset, only: r_SIDepth ! sea ice grid on interface of layer
! 宣言文 ; Declaration statements
!
real(DP), intent(in ) :: xy_SeaIceThickness(0:imax-1, 1:jmax)
!
! Sea ice thickness
integer , intent(out) :: xy_SOSILocalKMax (0:imax-1, 1:jmax)
real(DP), intent(out) :: xyr_SOSILocalDepth(0:imax-1, 1:jmax, 0:ksimax)
real(DP), intent(out) :: xyz_SOSILocalDepth(0:imax-1, 1:jmax, 1:ksimax)
! 作業変数
! Work variables
!
integer:: i ! 経度方向に回る DO ループ用作業変数
! Work variables for DO loop in longitude
integer:: j ! 緯度方向に回る DO ループ用作業変数
! Work variables for DO loop in latitude
integer:: k
! 実行文 ; Executable statement
!
! 初期化確認
! Initialization check
!
if ( .not. sosi_utils_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
do i = 0, imax-1
do j = 1, jmax
if ( xy_SeaIceThickness(i,j) == 0.0_DP ) then
xy_SOSILocalKMax(i,j) = 0
else if ( - xy_SeaIceThickness(i,j) < r_SIDepth(ksimax) ) then
xy_SOSILocalKMax(i,j) = ksimax
else
xy_SOSILocalKMax(i,j) = 0
search_ksimax : do k = 0+1, ksimax
!!$ if ( - xy_SeaIceThickness(i,j) >= r_SIDepth(k) ) then
! This SIDepthMergin avoids very thin lowest layer.
if ( - xy_SeaIceThickness(i,j) >= r_SIDepth(k)-SIDepthMergin ) then
xy_SOSILocalKMax(i,j) = k
exit search_ksimax
end if
end do search_ksimax
end if
end do
end do
do j = 1, jmax
do i = 0, imax-1
do k = 0, xy_SOSILocalKMax(i,j)-1
xyr_SOSILocalDepth(i,j,k) = r_SIDepth(k)
end do
k = xy_SOSILocalKMax(i,j)
xyr_SOSILocalDepth(i,j,k) = - xy_SeaIceThickness(i,j)
do k = xy_SOSILocalKMax(i,j)+1, ksimax
xyr_SOSILocalDepth(i,j,k) = -1.0e100_DP
end do
!
do k = 1, xy_SOSILocalKMax(i,j)
xyz_SOSILocalDepth(i,j,k) = ( xyr_SOSILocalDepth(i,j,k-1) + xyr_SOSILocalDepth(i,j,k) ) / 2.0_DP
end do
do k = xy_SOSILocalKMax(i,j)+1, ksimax
xyz_SOSILocalDepth(i,j,k) = -1.0e100_DP
end do
end do
end do
end subroutine SOSIUtilsSetSOSeaIceLevels
| Constant : | |
| SOSeaIceTempMissingValue = -99999.0_DP : | real(DP), parameter, public |
| Constant : | |||
| SIDepthMergin = 1.0e-3_DP : | real(DP), parameter
|
| Variable : | |||
| sosi_utils_inited = .false. : | logical, save
|