SCALE-RM
Functions/Subroutines
scale_ocean_phy_slab Module Reference

module OCEAN / Physics Slab model More...

Functions/Subroutines

subroutine, public ocean_phy_slab_setup (OCEAN_TYPE)
 Setup. More...
 
subroutine, public ocean_phy_slab (OCEAN_TEMP_t, OCEAN_TEMP, OCEAN_SFLX_WH, OCEAN_SFLX_prec, OCEAN_SFLX_evap, dt)
 Slab ocean model. More...
 

Detailed Description

module OCEAN / Physics Slab model

Description
ocean physics module, slab model
Author
Team SCALE
NAMELIST
  • PARAM_OCEAN_PHY_SLAB
    nametypedefault valuecomment
    OCEAN_PHY_SLAB_DEPTH real(RP) 10.0_RP water depth of slab ocean [m]

History Output
No history output

Function/Subroutine Documentation

◆ ocean_phy_slab_setup()

subroutine, public scale_ocean_phy_slab::ocean_phy_slab_setup ( character(len=*), intent(in)  OCEAN_TYPE)

Setup.

Definition at line 51 of file scale_ocean_phy_slab.F90.

References scale_const::const_cl, scale_const::const_dwatr, scale_stdio::io_fid_conf, scale_stdio::io_fid_log, scale_stdio::io_fid_nml, scale_stdio::io_l, scale_stdio::io_nml, and scale_process::prc_mpistop().

Referenced by scale_ocean_phy::ocean_phy_setup().

51  use scale_process, only: &
53  use scale_const, only: &
54  dwatr => const_dwatr, &
55  cl => const_cl
56  implicit none
57 
58  character(len=*), intent(in) :: OCEAN_TYPE
59 
60  namelist / param_ocean_phy_slab / &
61  ocean_phy_slab_depth
62 
63  integer :: ierr
64  !---------------------------------------------------------------------------
65 
66  if( io_l ) write(io_fid_log,*)
67  if( io_l ) write(io_fid_log,*) '++++++ Module[SLAB] / Categ[OCEAN PHY] / Origin[SCALElib]'
68 
69  !--- read namelist
70  rewind(io_fid_conf)
71  read(io_fid_conf,nml=param_ocean_phy_slab,iostat=ierr)
72  if( ierr < 0 ) then !--- missing
73  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
74  elseif( ierr > 0 ) then !--- fatal error
75  write(*,*) 'xxx Not appropriate names in namelist PARAM_OCEAN_PHY_SLAB. Check!'
76  call prc_mpistop
77  endif
78  if( io_nml ) write(io_fid_nml,nml=param_ocean_phy_slab)
79 
80  ocean_phy_slab_heatcapacity = dwatr * cl * ocean_phy_slab_depth
81 
82  if( io_l ) write(io_fid_log,*)
83  if( io_l ) write(io_fid_log,*) '*** Slab ocean depth [m] : ', ocean_phy_slab_depth
84  if( io_l ) write(io_fid_log,*) '*** Ocean heat capacity [J/K/m2] : ', ocean_phy_slab_heatcapacity
85 
86  return
subroutine, public prc_mpistop
Abort MPI.
module PROCESS
module CONSTANT
Definition: scale_const.F90:14
integer, public io_fid_conf
Config file ID.
Definition: scale_stdio.F90:55
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ocean_phy_slab()

subroutine, public scale_ocean_phy_slab::ocean_phy_slab ( real(rp), dimension (ia,ja), intent(out)  OCEAN_TEMP_t,
real(rp), dimension (ia,ja), intent(in)  OCEAN_TEMP,
real(rp), dimension (ia,ja), intent(in)  OCEAN_SFLX_WH,
real(rp), dimension(ia,ja), intent(in)  OCEAN_SFLX_prec,
real(rp), dimension(ia,ja), intent(in)  OCEAN_SFLX_evap,
real(dp), intent(in)  dt 
)

Slab ocean model.

Definition at line 98 of file scale_ocean_phy_slab.F90.

References scale_grid_index::ie, scale_stdio::io_fid_log, scale_stdio::io_l, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, and scale_landuse::landuse_fact_ocean.

Referenced by scale_ocean_phy::ocean_phy_setup().

98  use scale_landuse, only: &
99  landuse_fact_ocean
100  implicit none
101 
102  real(RP), intent(out) :: OCEAN_TEMP_t (IA,JA)
103  real(RP), intent(in) :: OCEAN_TEMP (IA,JA)
104  real(RP), intent(in) :: OCEAN_SFLX_WH (IA,JA)
105  real(RP), intent(in) :: OCEAN_SFLX_prec(IA,JA)
106  real(RP), intent(in) :: OCEAN_SFLX_evap(IA,JA)
107  real(DP), intent(in) :: dt
108 
109  integer :: i, j
110  !---------------------------------------------------------------------------
111 
112  if( io_l ) write(io_fid_log,*) '*** Ocean physics step: Slab'
113 
114  do j = js, je
115  do i = is, ie
116  if( landuse_fact_ocean(i,j) > 0.0_rp ) then
117  ocean_temp_t(i,j) = - ocean_sflx_wh(i,j) / ocean_phy_slab_heatcapacity
118  else
119  ocean_temp_t(i,j) = 0.0_rp
120  endif
121  enddo
122  enddo
123 
124  return
module LANDUSE
Here is the caller graph for this function: