mod_ocean_component.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!> module Oceanic component
!!
!! @par Description
!!          Oceanic component module
!!
!! @author Yuta Kawai, Team SCALE
!!
!<
!-------------------------------------------------------------------------------
#include "scaleFElib.h"
module mod_ocean_component
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_io
  use scale_prof
  use scale_prc
 
  use scale_element_base, only: &
    elementbase2d, elementbase3d
  use scale_mesh_base, only: meshbase
  use scale_mesh_base2d, only: meshbase2d
  use scale_mesh_base3d, only: meshbase3d
  use scale_localmesh_2d, only: localmesh2d
  use scale_localmesh_3d, only: localmesh3d
  use scale_localmeshfield_base, only: &
    localmeshfieldbase, localmeshfieldbaselist
  use scale_model_component, only: modelcomponent
  use scale_model_mesh_manager, only: modelmesh3d
 
  use mod_ocean_mesh, only: oceanmesh
  use mod_ocean_mesh_rm, only: oceanmeshrm
  use mod_ocean_mesh_gm, only: oceanmeshgm
 
  use mod_cpl_component, only: couplercomponent
 
  use mod_ocean_vars, only: &
    oceanvars,                                       &
    prgvar_therm_id => prgvar_therm_id,              &
    sfc_temp_id => auxvar2d_sfc_temp_id,             &
    atm_sfc_dens_id => atmvar2d_sfc_dens_id,         &
    atm_sfc_pres_id => atmvar2d_sfc_pres_id,         &
    atm_temp_id => atmvar2d_atm_temp_id,             &
    atm_dens_id => atmvar2d_atm_dens_id,             &
    atm_pres_id => atmvar2d_atm_pres_id,             &
    atm_w_id    => atmvar2d_atm_w_id,                &
    atm_u_id    => atmvar2d_atm_u_id,                &
    atm_v_id    => atmvar2d_atm_v_id,                &
    atm_qv_id   => atmvar2d_atm_qv_id,               &
    sflx_rd_sw_dir_id => atmvar2d_sflx_rd_sw_dir_id, &
    sflx_rd_lw_dif_id => atmvar2d_sflx_rd_lw_dif_id, &
    dz_a_id => atmvar2d_lowest_layer_dz_id
 
 
  use mod_ocean_dyn, only: oceandyn
 
  !-----------------------------------------------------------------------------
  implicit none
  private
  !-----------------------------------------------------------------------------
  !
  !++ Public type & procedure
  !
 
  !> Derived type to manage oceanic component
  type, extends(modelcomponent), public :: oceancomponent
    type(oceanvars) :: vars             !< Object to mange variables with oceanic component
 
    character(len=H_SHORT) :: mesh_type !< Type name of mesh oceanic component
    class(oceanmesh), pointer :: mesh   !< Pointer of mesh oceanic component
    type(oceanmeshrm) :: mesh_rm        !< Object to manage mesh for the case of regional mode
    type(oceanmeshgm) :: mesh_gm        !< Object to manage mesh for the case of global mode
 
    type(oceandyn) :: dyn_proc          !< Object to manage dynamical process
 
    type(couplercomponent), pointer :: coupler_ptr !< Pointer of coupler component
  contains
    procedure, public :: setup => ocean_setup 
    procedure, public :: setup_vars => ocean_setup_vars   
    procedure, public :: set_coupler => ocean_set_coupler 
    procedure, public :: calc_tendency => ocean_calc_tendency
    procedure, public :: update => ocean_update
    procedure, public :: set_surface => ocean_set_surface
    procedure, public :: get_surface => ocean_get_surface
    procedure, public :: finalize => ocean_finalize
  end type oceancomponent
 
  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private procedure
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
  !-----------------------------------------------------------------------------
contains
 
  !> Setup an object to mange oceanic component
!OCL SERIAL
  subroutine ocean_setup( this )
    use scale_const, only: &
      undef8 => const_undef8
    
    use scale_file_history_meshfield, only: &
      file_history_meshfield_setup  
    use scale_time_manager, only: &
      time_manager_regist_component
 
    implicit none
    
    class(oceancomponent), intent(inout), target :: this
 
    logical :: ACTIVATE_FLAG = .false. !< Flag whether oceanic component is activated
 
    real(DP) :: TIME_DT                             = undef8  !< Timestep value of oceanic component
    character(len=H_SHORT) :: TIME_DT_UNIT          = 'SEC'   !< Timestep unit of oceanic component
    real(DP) :: TIME_DT_RESTART                     = undef8  !< Timestep value when outputting restart file for oceanic component
    character(len=H_SHORT) :: TIME_DT_RESTART_UNIT  = 'SEC'   !< Timestep unit when outputting restart file for oceanic component
 
    logical :: OCEAN_DYN_DO    = .true.  !< Flag whether dynamics process is considered
    character(len=H_SHORT) :: OCEAN_MESH_TYPE = 'REGIONAL'  !< Name of mesh type for oceanic component ('REGIONAL' or 'GLOBAL')
 
 
    namelist / param_ocean / &
      activate_flag,         &
      time_dt,               &
      time_dt_unit,          &
      time_dt_restart,       &
      time_dt_restart_unit,  &
      ocean_mesh_type,       &
      ocean_dyn_do
    
    integer :: ierr
    !--------------------------------------------------
 
    log_info('OceanComponent_setup',*) 'Oceanic model components '
    
    !--- read namelist
    rewind(io_fid_conf)
    read(io_fid_conf,nml=param_ocean,iostat=ierr)
    if( ierr < 0 ) then !--- missing
      log_info("Ocean_setup",*) 'Not found namelist. Default used.'
    elseif( ierr > 0 ) then !--- fatal error
      log_error("Ocean_setup",*) 'Not appropriate names in namelist PARAM_OCEAN. Check!'
      call prc_abort
    endif
    log_nml(param_ocean)
 
    !************************************************
    call this%ModelComponent_Init('OCEAN', activate_flag )
 
    if ( .not. activate_flag ) return
    
    call prof_rapstart( 'Ocean_setup', 1)
 
    !- Setup time manager
 
    call this%time_manager%Init( this%GetComponentName(), &
      time_dt, time_dt_unit,                              &
      time_dt_restart, time_dt_restart_unit               ) 
    
    call time_manager_regist_component( this%time_manager )
    
    !- Setup mesh & file I/O for oceanic component
    
    this%mesh_type = ocean_mesh_type
    select case( this%mesh_type )
    case('REGIONAL')
      call this%mesh_rm%Init()
      call file_history_meshfield_setup( mesh3d_=this%mesh_rm%mesh, & ! (in)
        dim_name_postfix_='_O',                                     & ! (in)
        registered_comp_id=this%vars%hist_comp_id                   ) ! (out)
      this%mesh => this%mesh_rm
    case('GLOBAL')
      call this%mesh_gm%Init()
      call file_history_meshfield_setup( meshcubedsphere3d_=this%mesh_gm%mesh, & ! (in)
        dim_name_postfix_='_O',                                                & ! (in)      
        registered_comp_id=this%vars%hist_comp_id )                              ! (out)
      this%mesh => this%mesh_gm
    case default
      log_error("Ocean_setup",*) 'Unsupported type of mesh is specified. Check!', this%mesh_type
      call prc_abort    
    end select
    
    !- setup common tools for oceanic model
 
    !- Setup each processes in ocean model ------------------------------------
 
    !- Setup the module for ocean / dynamics 
    call this%dyn_proc%ModelComponentProc_Init( 'OceanDyn', ocean_dyn_do )
    call this%dyn_proc%setup( this%mesh, this%time_manager )
    
    !- Setup
 
    log_newline
    log_info('OceanComponent_setup',*) 'Finish setup of each oceanic component.'
 
    call prof_rapend( 'Ocean_setup', 1)
 
    return
  end subroutine ocean_setup
 
  !> Setup variables with the oceanic component
!OCL SERIAL
  subroutine ocean_setup_vars( this )
    implicit none
    class(oceancomponent), intent(inout) :: this
    !----------------------------------------------------------
 
    call prof_rapstart( 'Ocean_setup_vars', 1)
    call this%vars%Init( this%mesh )
    call prof_rapend( 'Ocean_setup_vars', 1)
    return
  end subroutine ocean_setup_vars
 
!> Set coupler component to the oceanic component
!OCL SERIAL
  subroutine ocean_set_coupler( this, coupler )
    implicit none
    class(oceancomponent), intent(inout) :: this
    class(couplercomponent), target, intent(inout) :: coupler
    !----------------------------------------------------------
    this%coupler_ptr => coupler
    return
  end subroutine ocean_set_coupler  
 
!> Calculate tendencies with the oceanic component
!OCL SERIAL
  subroutine ocean_calc_tendency( this, force )
    use mod_ocean_vars, only: &
      alb_vis_dir_id => auxvar2d_sfc_alb_vis_dir_id, &
      sflx_mw_id => ocn_sflx_mw_id, &
      sflx_mu_id => ocn_sflx_mu_id, &
      sflx_mv_id => ocn_sflx_mv_id, &
      sflx_sh_id => ocn_sflx_sh_id, &
      sflx_lh_id => ocn_sflx_lh_id, &
      sflx_qv_id => ocn_sflx_qv_id, &
      sflx_rd_sw_dir_id => atmvar2d_sflx_rd_sw_dir_id, &
      sflx_rd_lw_dif_id => atmvar2d_sflx_rd_lw_dif_id, &
      rhoh_id => phytend_rhoh_id
    implicit none
    class(oceancomponent), intent(inout) :: this
    logical, intent(in) :: force
 
    class(meshbase3d), pointer :: mesh3D
    class(localmesh2d), pointer :: lmesh2D
    class(localmesh3d), pointer :: lmesh3D
    integer :: idom
 
    real(RP), allocatable :: SFLX_GH(:,:)
    !------------------------------------------------------------------
    
    !- Get surface data from coupler 
    call this%get_surface()
 
    call prof_rapstart( 'OCN_tendency', 1)
 
    mesh3d => this%mesh%ptr_mesh
 
    do idom=1, mesh3d%LOCAL_MESH_NUM
      lmesh3d => mesh3d%lcmesh_list(idom)
      lmesh2d => lmesh3d%lcmesh2D
 
      allocate( sflx_gh(lmesh2d%refElem2D%Np,lmesh2d%Ne) )
 
      !- Calculate surface fluxes
 
      call calculate_surface_flux( this, &
        this%vars%OCN_SFLX(sflx_mw_id)%local(idom)%val, &
        this%vars%OCN_SFLX(sflx_mu_id)%local(idom)%val, &
        this%vars%OCN_SFLX(sflx_mv_id)%local(idom)%val, &
        this%vars%OCN_SFLX(sflx_sh_id)%local(idom)%val, &
        this%vars%OCN_SFLX(sflx_lh_id)%local(idom)%val, &
        this%vars%OCN_SFLX(sflx_qv_id)%local(idom)%val, &
        !-
        this%vars%AUX_VARS2D(sfc_temp_id)%local(idom)%val, &
        this%vars%ATM_VARS2D(atm_sfc_dens_id)%local(idom)%val, &
        this%vars%ATM_VARS2D(atm_sfc_pres_id)%local(idom)%val, &
        this%vars%ATM_VARS2D(atm_temp_id)%local(idom)%val, &
        this%vars%ATM_VARS2D(atm_dens_id)%local(idom)%val, &
        this%vars%ATM_VARS2D(atm_pres_id)%local(idom)%val, &
        this%vars%ATM_VARS2D(atm_w_id)%local(idom)%val,  &
        this%vars%ATM_VARS2D(atm_u_id)%local(idom)%val,  &
        this%vars%ATM_VARS2D(atm_v_id)%local(idom)%val,  &
        this%vars%ATM_VARS2D(atm_qv_id)%local(idom)%val, &
        this%vars%ATM_VARS2D(dz_a_id)%local(idom)%val,   &
        lmesh2d, lmesh2d%refElem2D )
 
      !- Calculate ocean surface heat flux
      
      call calculate_sfc_heat_flux( sflx_gh,                     & ! (out)
        this%vars%AUX_VARS2D(sfc_temp_id)%local(idom)%val,       & ! (in)
        this%vars%ATM_VARS2D(sflx_rd_sw_dir_id)%local(idom)%val, & ! (in)
        this%vars%ATM_VARS2D(sflx_rd_lw_dif_id)%local(idom)%val, & ! (in)
        this%vars%AUX_VARS2D(alb_vis_dir_id)%local(idom)%val,    & ! (in)
        this%vars%OCN_SFLX(sflx_sh_id)%local(idom)%val,          & ! (in)
        this%vars%OCN_SFLX(sflx_lh_id)%local(idom)%val,          & ! (in)
        lmesh2d, lmesh2d%refElem2D )                               ! (in)
 
      !- Calculate tendencies with physics
 
      call calculate_phys_tendency( this%vars%PHY_TEND(rhoh_id)%local(idom)%val, & ! (out)
        sflx_gh, lmesh3d%zlev,                                                   & ! (in)
        lmesh3d, lmesh3d%refElem3D, lmesh2d, lmesh2d%refElem2D )                   ! (in)
      
      ! write(*,*) "------"
      ! write(*,*) "SFC_TEMP = ", this%vars%AUX_VARS2D(SFC_TEMP_ID)%local(idom)%val(1,1)
      ! write(*,*) "SFC_DENS = ", this%vars%ATM_VARS2D(ATM_SFC_DENS_ID)%local(idom)%val(1,1)
      ! write(*,*) "SFC_PRES = ", this%vars%ATM_VARS2D(ATM_SFC_PRES_ID)%local(idom)%val(1,1)
      ! write(*,*) "ATM_TEMP = ", this%vars%ATM_VARS2D(ATM_TEMP_ID)%local(idom)%val(1,1)
      ! write(*,*) "ATM_PRES = ", this%vars%ATM_VARS2D(ATM_PRES_ID)%local(idom)%val(1,1)
      ! write(*,*) "ATM_QV= ", this%vars%ATM_VARS2D(ATM_QV_ID)%local(idom)%val(1,1)
      ! write(*,*) "RD_SW_DIR = ", this%vars%ATM_VARS2D(SFLX_RD_SW_DIR_ID)%local(idom)%val(1,1)
      ! write(*,*) "RD_LW_DIF = ", this%vars%ATM_VARS2D(SFLX_RD_LW_DIF_ID)%local(idom)%val(1,1)
      ! write(*,*) "ZLEV_A = ", this%vars%ATM_VARS2D(ZLEV_A_ID)%local(idom)%val(1,1)
      ! write(*,*) "SLFX_MW=", this%vars%OCN_SFLX(SFLX_MW_ID)%local(idom)%val(1,1)
      ! write(*,*) "SLFX_MU=", this%vars%OCN_SFLX(SFLX_MU_ID)%local(idom)%val(1,1)
      ! write(*,*) "SLFX_MV=", this%vars%OCN_SFLX(SFLX_MV_ID)%local(idom)%val(1,1)
      ! write(*,*) "SH=", this%vars%OCN_SFLX(SFLX_SH_ID)%local(idom)%val(1,1)
      ! write(*,*) "LH=", this%vars%OCN_SFLX(SFLX_LH_ID)%local(idom)%val(1,1)
      ! write(*,*) 'SFLX_GH = ', SFLX_GH(1,1)
 
      deallocate( sflx_gh )
    end do
 
    call prof_rapend( 'OCN_tendency', 1)
 
    !-
    call this%set_surface( countup=.true. )
 
    return  
  end subroutine ocean_calc_tendency
 
!> Update variables with the oceanic component
!OCL SERIAL
  subroutine ocean_update( this )
    implicit none
    class(oceancomponent), intent(inout) :: this
    
    integer :: tm_process_id
    logical :: is_update
    integer :: inner_itr
 
    integer :: idom
    class(localmesh2d), pointer :: lmesh2D
    class(localmesh3d), pointer :: lmesh
    !--------------------------------------------------
    call prof_rapstart( 'OCN_update', 1)
 
    if ( this%dyn_proc%IsActivated() ) then
      call prof_rapstart('OCN_Dynamics', 1)
      tm_process_id = this%dyn_proc%tm_process_id
      is_update = this%time_manager%Do_process( tm_process_id )
 
      log_progress(*) 'ocean / dynamics'
      do inner_itr=1, this%time_manager%Get_process_inner_itr_num( tm_process_id )
        call this%dyn_proc%update( &
          this%mesh, this%vars%PROGVARS_manager, this%vars%QTRCVARS_manager, &
          this%vars%AUXVARS_manager, this%vars%PHYTENDS_manager, is_update   )
      end do
      call prof_rapend('OCN_Dynamics', 1)
    end if
 
    !- Set surface temperature
    do idom=1, this%mesh%ptr_mesh%LOCAL_MESH_NUM
      lmesh => this%mesh%ptr_mesh%lcmesh_list(idom)
      lmesh2d => lmesh%lcmesh2D
      call set_sfctemp_lc( this%vars%AUX_VARS2D(sfc_temp_id)%local(idom)%val, & ! (out)
        this%vars%PROG_VARS(prgvar_therm_id)%local(idom)%val,                 & ! (in)
        lmesh, lmesh%refElem3D, lmesh2d, lmesh2d%refElem2D )                    ! (in)
      ! write(*,*) "OCN_update: SFC_TEMP = ", this%vars%AUX_VARS2D(SFC_TEMP_ID)%local(idom)%val(1,1)
    end do
 
    call prof_rapend('OCN_update', 1)
    return  
  end subroutine ocean_update
!OCL SERIAL
  subroutine set_sfctemp_lc( SFC_TEMP, &
    THERM, lmesh, elem, lmesh2D, elem2D )
    implicit none
    class(localmesh3d), intent(in) :: lmesh
    class(elementbase3d), intent(in) :: elem
    class(localmesh2d), intent(in) :: lmesh2D
    class(elementbase2d), intent(in) :: elem2D
    real(RP), intent(out) :: SFC_TEMP(elem2D%Np,lmesh2D%NeA)
    real(RP), intent(in) :: THERM(elem%Np,lmesh%NeA)
 
    integer :: ke, ke2D
    integer :: hSlice(elem2D%Np)
    !------------------------------------------------------------------
 
    hslice(:) = elem%Hslice(:,elem%Nnode_v)
    !$omp parallel do private(ke)
    do ke2d=lmesh2d%NeS, lmesh2d%NeE
      ke = ke2d + (lmesh%NeZ-1)*lmesh%Ne2D
      sfc_temp(:,ke2d) = therm(hslice(:),ke)
    end do
    return
  end subroutine set_sfctemp_lc
 
  !> Set ocean quantities to coupler component
!OCL SERIAL
  subroutine ocean_set_surface( this, countup )
    implicit none
    class(oceancomponent), intent(inout) :: this
    logical, intent(in) :: countup
    !--------------------------------------------------
 
    call prof_rapstart( 'OCN_sfc_exch', 1)
 
    call this%coupler_ptr%vars%PutOCN( this%vars, countup )
 
    call prof_rapend( 'OCN_sfc_exch', 1)
    return
  end subroutine ocean_set_surface
 
  !> Get atmospheric quantities from coupler component
!OCL SERIAL
  subroutine ocean_get_surface( this )
    implicit none
    class(oceancomponent), intent(inout) :: this
    !--------------------------------------------------
 
    call prof_rapstart( 'OCN_sfc_exch', 1)
 
    call this%coupler_ptr%vars%Get_ATM_OCN( &
      this%vars%ATM_VARS2D(atm_sfc_dens_id),   &
      this%vars%ATM_VARS2D(atm_sfc_pres_id),   &
      this%vars%ATM_VARS2D(atm_temp_id),       &
      this%vars%ATM_VARS2D(atm_pres_id),       &
      this%vars%ATM_VARS2D(atm_w_id),          &
      this%vars%ATM_VARS2D(atm_u_id),          &
      this%vars%ATM_VARS2D(atm_v_id),          &
      this%vars%ATM_VARS2D(atm_qv_id),         &
      this%vars%ATM_VARS2D(sflx_rd_sw_dir_id), &
      this%vars%ATM_VARS2D(sflx_rd_lw_dif_id), &
      this%vars%ATM_VARS2D(dz_a_id) )
    
    call prof_rapend( 'OCN_sfc_exch', 1)
    return
  end subroutine ocean_get_surface
 
!> Finalize an object to manage the ocean component
!OCL SERIAL
  subroutine ocean_finalize( this )
    implicit none
    class(oceancomponent), intent(inout) :: this
    !--------------------------------------------------
 
    log_info('OceanComponent_finalize',*)
 
    if ( .not. this%IsActivated() ) return
 
    ! call this%dyn_proc%finalize()
    call this%vars%Final()
 
    select case( this%mesh_type )
    case('REGIONAL')
      call this%mesh_rm%Final()
    case('GLOBAL')
      call this%mesh_gm%Final()
    end select  
    this%mesh => null()
 
    call this%time_manager%Final()
 
    return  
  end subroutine ocean_finalize
 
!- Private subroutines -------------------------------------------------------------
 
  !> Calculate momentum and heat flux at the surface (tentative)
!OCL SERIAL
  subroutine calculate_surface_flux( this, &
    SFLX_MW, SFLX_MU, SFLX_MV, SFLX_SH, SFLX_LH, SFLX_QV, &
    SFC_TEMP, SFC_DENS, SFC_PRES,                                       &
    ATM_TEMP, ATM_DENS, ATM_PRES, ATM_W, ATM_U, ATM_V, ATM_QV, zlev_a, &
    lmesh, elem )
    use scale_const, only: &
      rplanet => const_radius    
    use scale_atm_phy_sf_bulk_simple, only: &
       atmos_phy_sf_simple_flux
    implicit none
    class(oceancomponent), intent(in) :: this
    class(localmesh2d), intent(in) :: lmesh
    class(elementbase2d), intent(in) :: elem
    real(RP), intent(out) :: SFLX_MW(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: SFLX_MU(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: SFLX_MV(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: SFLX_SH(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: SFLX_LH(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: SFLX_QV(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: SFC_TEMP(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: SFC_DENS(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: SFC_PRES(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: ATM_TEMP(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: ATM_DENS(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: ATM_PRES(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: ATM_W(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: ATM_U(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: ATM_V(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: ATM_QV(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: zlev_a(elem%Np,lmesh%NeA)
    
    ! dummy
    real(RP) :: U10(elem%Np,lmesh%NeA)
    real(RP) :: V10(elem%Np,lmesh%NeA)
 
    integer :: ke
    integer :: ke2D, ij
 
    real(RP) :: DZ1     (elem%Np,lmesh%NeA)
    real(RP) :: Z1      (elem%Np,lmesh%Ne)
    real(RP) :: ATM_W_lc(elem%Np,lmesh%NeA)
    real(RP) :: ATM_U_lc(elem%Np,lmesh%NeA)
    real(RP) :: ATM_V_lc(elem%Np,lmesh%NeA)
    !--------------------------------------------------
 
    !$omp parallel do collapse(2)
    do ke2d=lmesh%NeS, lmesh%NeE
    do ij=1, elem%Np
 
      atm_w_lc(ij,ke2d) = atm_w(ij,ke2d)
      atm_u_lc(ij,ke2d) = atm_u(ij,ke2d)
      atm_v_lc(ij,ke2d) = atm_v(ij,ke2d)
 
      z1(ij,ke2d) = zlev_a(ij,ke2d)
      dz1(ij,ke2d) = z1(ij,ke2d) - zlev_a(ij,ke2d)
      z1(ij,ke2d) = z1(ij,ke2d) + rplanet
    end do
    end do
 
    call convert_uv2localorthvec( &
      this%mesh%ptr_mesh, lmesh%pos_en(:,:,1), lmesh%pos_en(:,:,2), z1(:,:), elem%Np*lmesh%Ne, &
      atm_u_lc(:,lmesh%NeS:lmesh%NeE), atm_v_lc(:,lmesh%NeS:lmesh%NeE) ) ! (inout)
 
    call atmos_phy_sf_simple_flux( &
      elem%Np, 1, elem%Np, lmesh%NeA, 1, lmesh%Ne,                               & ! (in)
      atm_w_lc(:,:), atm_u_lc(:,:), atm_v_lc(:,:), atm_temp(:,:), atm_pres(:,:), & ! (in) Note: ATM_PRES is not used
      atm_qv(:,:),                                                               & ! (in)
      sfc_dens(:,:), sfc_temp(:,:), sfc_pres(:,:),                               & ! (in)
      dz1(:,:),                                                                  & ! (in)
      sflx_mw(:,:), sflx_mu(:,:), sflx_mv(:,:),                                  & ! (out)
      sflx_sh(:,:), sflx_lh(:,:), sflx_qv(:,:),                                  & ! (out)
      u10(:,:), v10(:,:)                                                         ) ! (out)
 
    call convert_localorth2uvvec( &
      this%mesh%ptr_mesh, lmesh%pos_en(:,:,1), lmesh%pos_en(:,:,2), z1(:,:), elem%Np*lmesh%Ne, &
      sflx_mu(:,lmesh%NeS:lmesh%NeE), sflx_mv(:,lmesh%NeS:lmesh%NeE) ) ! (inout)
    
    return
  end subroutine calculate_surface_flux
 
  !> Calculate ground heat flux (tentative)
  !!
!OCL SERIAL
  subroutine calculate_sfc_heat_flux( sflx_GH, &
    SFC_TEMP,                             & 
    SFLX_RD_SW_dn_dir, SFLX_RD_LW_dn_dif, &
    SFC_ALB_dir_vis,                      &
    SFLX_SH, SFLX_LH,                     &
    lmesh, elem )
    use scale_const, only: &
       stb   => const_stb
    implicit none
    class(localmesh2d), intent(in) :: lmesh
    class(elementbase2d), intent(in) :: elem
    real(RP), intent(out) :: sflx_GH(elem%Np,lmesh%Ne)
    real(RP), intent(in) :: SFC_TEMP(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: SFLX_RD_SW_dn_dir(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: SFLX_RD_LW_dn_dif(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: SFC_ALB_dir_vis(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: SFLX_SH(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: SFLX_LH(elem%Np,lmesh%NeA)
 
    integer :: ke
    real(RP) :: emis(elem%Np)
    real(RP) :: SWU(elem%Np), SWD(elem%Np), LWU(elem%Np), LWD(elem%Np)
    !--------------------------------------------------
 
    !$omp parallel do private(SWU,SWD,LWU,LWD,emis)
    do ke=lmesh%NeS, lmesh%NeE
      emis(:) = stb * sfc_temp(:,ke)**4
 
      lwd(:) = sflx_rd_lw_dn_dif(:,ke)
      lwu(:) = emis(:)
      swd(:) = sflx_rd_sw_dn_dir(:,ke)
      swu(:) = sfc_alb_dir_vis(:,ke) * swd(:)
 
      sflx_gh(:,ke) = swd(:) - swu(:) + lwd(:) - lwu(:) - sflx_sh(:,ke) - sflx_lh(:,ke)
    end do
    return
  end subroutine calculate_sfc_heat_flux
 
!OCL SERIAL
  subroutine calculate_phys_tendency( RHOH, &
    sflx_GH, zlev, lmesh, elem, lmesh2D, elem2D )
    implicit none
    class(localmesh3d), intent(in) :: lmesh
    class(elementbase3d), intent(in) :: elem
    class(localmesh2d), intent(in) :: lmesh2D
    class(elementbase2d), intent(in) :: elem2D
    real(RP), intent(out) :: RHOH(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: sflx_GH(elem2D%Np,lmesh2D%Ne)
    real(RP), intent(in) :: zlev(elem%Np,lmesh%Ne2D,lmesh%NeZ)
 
    integer :: ke, ke_z, ke2D
    integer :: ph, pz, p
    real(RP) :: r_ocn_depth(elem2D%Np,lmesh2D%Ne)
    integer :: hSlice_t(elem2D%Np), hSlice_b(elem2D%Np)
    !---------------------------------------
 
    hslice_t(:) = elem%Hslice(:,elem%Nnode_v)
    hslice_b(:) = elem%Hslice(:,1)
 
    !$omp parallel private(ke)
    !$omp do
    do ke2d=1, lmesh%Ne2D
      r_ocn_depth(:,ke2d) = 1.0_rp / ( zlev(hslice_t(:),ke2d,lmesh%NeZ) - zlev(hslice_b(:),ke2d,1) )
    end do
    !$omp do collapse(2)
    do ke_z=1, lmesh%NeZ
    do ke2d=1, lmesh%Ne2D
      ke = ke2d + (ke_z-1)*lmesh%Ne2D
      do pz=1, elem%Nnode_v
      do ph=1, elem%Nnode_h1D**2
        p = ph + (pz-1)*elem%Nnode_h1D**2
        rhoh(p,ke) = sflx_gh(ph,ke2d) * r_ocn_depth(ph,ke2d)
      end do
      end do
    end do
    end do
    !$omp end parallel
    return
  end subroutine calculate_phys_tendency
 
!OCL SERIAL
  subroutine convert_uv2localorthvec( mesh, x1, x2, r, N, & ! (in)
    u, v ) ! (inout)
    use scale_cubedsphere_coord_cnv, only: &
      cubedspherecoordcnv_cs2localorthvec_alpha
    use scale_mesh_cubedspheredom3d, only: meshcubedspheredom3d
 
    implicit none
    integer, intent(in) :: N
    class(meshbase3d), intent(in) :: mesh
    real(RP), intent(in) :: x1(N)
    real(RP), intent(in) :: x2(N)
    real(RP), intent(in) :: r(N)
    real(RP), intent(inout) :: U(N)
    real(RP), intent(inout) :: V(N)
 
    select type(mesh)
    type is (meshcubedspheredom3d)
      call cubedspherecoordcnv_cs2localorthvec_alpha( &
        x1, x2, r, n, & ! (in)
        u, v          ) ! (inout)
    end select
    return
  end subroutine convert_uv2localorthvec
 
!OCL SERIAL
  subroutine convert_localorth2uvvec( mesh, x1, x2, r, N, & ! (in)
    u, v ) ! (inout)
    use scale_cubedsphere_coord_cnv, only: &
      cubedspherecoordcnv_localorth2csvec_alpha
    use scale_mesh_cubedspheredom3d, only: meshcubedspheredom3d
 
    implicit none
    integer, intent(in) :: N
    class(meshbase3d), intent(in) :: mesh
    real(RP), intent(in) :: x1(N)
    real(RP), intent(in) :: x2(N)
    real(RP), intent(in) :: r(N)
    real(RP), intent(inout) :: U(N)
    real(RP), intent(inout) :: V(N)
 
    select type(mesh)
    type is (meshcubedspheredom3d)
      call cubedspherecoordcnv_localorth2csvec_alpha( &
        x1, x2, r, n, & ! (in)
        u, v          ) ! (inout)
    end select
    return
  end subroutine convert_localorth2uvvec 
 
end module mod_ocean_component