coriolis.F90

Go to the documentation of this file.

module coriolis_mod
  use datadefn_mod, only : default_precision
  use monc_component_mod, only : component_descriptor_type
  use optionsdatabase_mod, only : options_get_logical, options_get_real
  use state_mod, only : model_state_type
  use grids_mod, only : z_index, y_index, x_index
  implicit none

#ifndef TEST_MODE
  private
#endif

  logical :: baroclinicity_use_geostrophic_shear
  real(kind=DEFAULT_PRECISION), dimension(:), allocatable :: geostrophic_wind_x, geostrophic_wind_y
  real(kind=DEFAULT_PRECISION) :: fcoriol

  public coriolis_get_descriptor

  contains

  type(component_descriptor_type) function coriolis_get_descriptor()
    coriolis_get_descriptor%name="coriolis"
    coriolis_get_descriptor%version=0.1
    coriolis_get_descriptor%initialisation=>initialisation_callback
    coriolis_get_descriptor%timestep=>timestep_callback
  end function coriolis_get_descriptor

  subroutine initialisation_callback(current_state)
    type(model_state_type), target, intent(inout) :: current_state

    integer :: k
    
    baroclinicity_use_geostrophic_shear=options_get_logical(current_state%options_database, "baroclinicity_use_geostrophic_shear")
    fcoriol=options_get_real(current_state%options_database, "fcoriol")
    current_state%geostrophic_wind_rate_of_change_in_x=options_get_real(current_state%options_database, &
         "geostrophic_wind_rate_of_change_in_x")
    current_state%geostrophic_wind_rate_of_change_in_y=options_get_real(current_state%options_database, &
         "geostrophic_wind_rate_of_change_in_y")
    current_state%surface_geostrophic_wind_x=options_get_real(current_state%options_database, "surface_geostrophic_wind_x")
    current_state%surface_geostrophic_wind_y=options_get_real(current_state%options_database, "surface_geostrophic_wind_y")

    allocate(geostrophic_wind_x(current_state%local_grid%size(z_index)), &
         geostrophic_wind_y(current_state%local_grid%size(z_index)))

    do k=1,current_state%local_grid%size(z_index)
      geostrophic_wind_x(k)=current_state%surface_geostrophic_wind_x
      geostrophic_wind_y(k)=current_state%surface_geostrophic_wind_y
      if (baroclinicity_use_geostrophic_shear) then
        geostrophic_wind_x(k)=geostrophic_wind_x(k)+current_state%geostrophic_wind_rate_of_change_in_x*&
             current_state%global_grid%configuration%vertical%zn(k)
        geostrophic_wind_y(k)=geostrophic_wind_y(k)+current_state%geostrophic_wind_rate_of_change_in_y*&
             current_state%global_grid%configuration%vertical%zn(k)
      end if
    end do
  end subroutine initialisation_callback

  subroutine timestep_callback(current_state)
    type(model_state_type), target, intent(inout) :: current_state

    integer :: k

    if (current_state%halo_column) then
      if (.not. ((current_state%column_local_y == current_state%local_grid%halo_size(y_index) .and. &
           current_state%column_local_x .le. current_state%local_grid%local_domain_end_index(x_index) .and. &
           current_state%column_local_x .ge. current_state%local_grid%local_domain_start_index(x_index)-1) .or. &
           (current_state%column_local_x == current_state%local_grid%halo_size(x_index) .and. &
           current_state%column_local_y .ge. current_state%local_grid%local_domain_start_index(y_index) &
           .and. current_state%column_local_y .le. current_state%local_grid%local_domain_end_index(y_index)) )) return
    end if
    
    do k=2,current_state%local_grid%size(z_index)
#if defined(U_ACTIVE) && defined(V_ACTIVE)
      current_state%su%data(k, current_state%column_local_y, current_state%column_local_x)=&
           current_state%su%data(k, current_state%column_local_y, current_state%column_local_x)+fcoriol*&
           (0.25_default_precision*(current_state%v%data(k, current_state%column_local_y, current_state%column_local_x)+&
           current_state%v%data(k, current_state%column_local_y, current_state%column_local_x+1)+&
           current_state%v%data(k, current_state%column_local_y-1, current_state%column_local_x)+&
           current_state%v%data(k, current_state%column_local_y-1, current_state%column_local_x+1))+current_state%vgal-&
           geostrophic_wind_y(k))

      current_state%sv%data(k, current_state%column_local_y, current_state%column_local_x)=&
           current_state%sv%data(k, current_state%column_local_y, current_state%column_local_x)-fcoriol*&
           (0.25_default_precision*(current_state%u%data(k, current_state%column_local_y, current_state%column_local_x)+&
           current_state%u%data(k, current_state%column_local_y, current_state%column_local_x-1)+&
           current_state%u%data(k, current_state%column_local_y+1, current_state%column_local_x)+&
           current_state%u%data(k, current_state%column_local_y+1, current_state%column_local_x-1))+current_state%ugal-&
           geostrophic_wind_x(k))
#endif
    end do
  end subroutine timestep_callback
end module coriolis_mod