scalar_ic.f90

Go to the documentation of this file.

! Copyright (c) 2021-2026, The Neko Authors
! All rights reserved.
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
!
!   * Redistributions of source code must retain the above copyright
!     notice, this list of conditions and the following disclaimer.
!
!   * Redistributions in binary form must reproduce the above
!     copyright notice, this list of conditions and the following
!     disclaimer in the documentation and/or other materials provided
!     with the distribution.
!
!   * Neither the name of the authors nor the names of its
!     contributors may be used to endorse or promote products derived
!     from this software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
module scalar_ic
  use gather_scatter, only : gs_t, gs_op_add
  use neko_config, only : neko_bcknd_device
  use num_types, only : rp
  use device_math, only : device_col2
  use device, only : device_memcpy, host_to_device, device_to_host
  use field, only : field_t
  use utils, only : neko_error, filename_chsuffix, filename_suffix, &
       neko_warning, neko_fname_len, extract_fld_file_index
  use coefs, only : coef_t
  use math, only : col2, cfill, cfill_mask
  use user_intf, only : user_initial_conditions_intf
  use json_module, only : json_file
  use json_utils, only : json_get, json_get_or_default, &
       json_get_or_lookup, json_get_or_lookup_or_default, &
       json_get_subdict_or_empty
  use point_zone, only : point_zone_t
  use point_zone_registry, only : neko_point_zone_registry
  use logger, only : neko_log, log_size
  use fld_file_data, only : fld_file_data_t
  use fld_file, only : fld_file_t
  use checkpoint, only : chkp_t
  use file, only : file_t
  use global_interpolation, only : global_interpolation_t
  use interpolation, only : interpolator_t
  use space, only : space_t, gll
  use field_list, only : field_list_t
  use import_field_utils, only : import_fields
  use expression, only : expression_eval_static
  implicit none
  private
 
  interface set_scalar_ic
     module procedure set_scalar_ic_int, set_scalar_ic_usr
  end interface set_scalar_ic
 
  public :: set_scalar_ic
 
contains
 
  subroutine set_scalar_ic_int(s, coef, gs, type, params, i)
    type(field_t), intent(inout) :: s
    type(coef_t), intent(in) :: coef
    type(gs_t), intent(inout) :: gs
    character(len=*) :: type
    type(json_file), intent(inout) :: params
    integer, intent(in) :: i
 
    ! Variables for retrieving JSON parameters
    real(kind=rp) :: ic_value
    character(len=:), allocatable :: read_str
    real(kind=rp) :: zone_value
 
    if (trim(type) .eq. 'uniform') then
 
       call json_get_or_lookup(params, 'value', ic_value)
       call set_scalar_ic_uniform(s, ic_value)
 
    else if (trim(type) .eq. 'expression') then
 
       call json_get(params, 'value', read_str)
       call set_scalar_ic_expression(s, read_str)
       if (allocated(read_str)) deallocate(read_str)
 
    else if (trim(type) .eq. 'point_zone') then
 
       call json_get_or_lookup(params, 'base_value', ic_value)
       call json_get(params, 'zone_name', read_str)
       call json_get_or_lookup(params, 'zone_value', zone_value)
 
       call set_scalar_ic_point_zone(s, ic_value, read_str, zone_value)
 
    else if (trim(type) .eq. 'field') then
 
       block
         character(len=NEKO_FNAME_LEN) :: fname, mesh_fname
         logical :: interpolate
         type(json_file) :: interp_subdict
         integer :: tgt_scal_idx
 
         call json_get(params, 'file_name', read_str)
         fname = trim(read_str)
 
         call json_get_or_default(params, 'interpolate', interpolate, &
              .false.)
 
         call json_get_or_default(params, 'mesh_file_name', read_str, &
              "none")
         mesh_fname = trim(read_str)
 
         ! Give the user the option to select which scalar they want to import
         ! the values from, in the fld file. 0 corresponds to temperature.
         call json_get_or_default(params, 'target_index', tgt_scal_idx, i)
 
         call json_get_subdict_or_empty(params, "interpolation", &
              interp_subdict)
 
         call set_scalar_ic_fld(s, fname, interpolate, mesh_fname, i, &
              tgt_scal_idx, interp_subdict)
 
       end block
 
    else
       call neko_error('Invalid initial condition')
    end if
 
    call set_scalar_ic_common(s, coef, gs)
 
  end subroutine set_scalar_ic_int
 
  subroutine set_scalar_ic_usr(scheme_name, s, coef, gs, user_proc)
    character(len=*), intent(in) :: scheme_name
    type(field_t), target, intent(inout) :: s
    type(coef_t), intent(in) :: coef
    type(gs_t), intent(inout) :: gs
    procedure(user_initial_conditions_intf) :: user_proc
    type(field_list_t) :: fields
 
    call neko_log%message("Type: user")
 
    call fields%init(1)
    call fields%assign_to_field(1, s)
 
    call user_proc(scheme_name, fields)
    call set_scalar_ic_common(s, coef, gs)
 
  end subroutine set_scalar_ic_usr
 
  subroutine set_scalar_ic_common(s, coef, gs)
    type(field_t), intent(inout) :: s
    type(coef_t), intent(in) :: coef
    type(gs_t), intent(inout) :: gs
    integer :: n
 
    n = s%dof%size()
    if (neko_bcknd_device .eq. 1) then
       call device_memcpy(s%x, s%x_d, n, host_to_device, sync = .false.)
    end if
 
    ! Ensure continuity across elements for initial conditions
    call gs%op(s%x, n, gs_op_add)
 
    if (neko_bcknd_device .eq. 1) then
       call device_col2(s%x_d, coef%mult_d, n)
    else
       call col2(s%x, coef%mult, n)
    end if
 
  end subroutine set_scalar_ic_common
 
  subroutine set_scalar_ic_uniform(s, ic_value)
    type(field_t), intent(inout) :: s
    real(kind=rp), intent(in) :: ic_value
    integer :: n
    character(len=LOG_SIZE) :: log_buf
 
    call neko_log%message("Type : uniform")
    write (log_buf, '(A,ES12.6)') "Value: ", ic_value
    call neko_log%message(log_buf)
 
    s = ic_value
    n = s%dof%size()
    if (neko_bcknd_device .eq. 1) then
       call cfill(s%x, ic_value, n)
    end if
 
  end subroutine set_scalar_ic_uniform
 
  subroutine set_scalar_ic_expression(s, expr)
    type(field_t), intent(inout) :: s
    character(len=*), intent(in) :: expr
 
    call neko_log%message("Type : expression")
    call neko_log%message("Value: " // trim(expr))
 
    call expression_eval_static(expr, s%x, s%dof%size(), &
         s%dof%x, s%dof%y, s%dof%z, 'scalar initial condition')
 
  end subroutine set_scalar_ic_expression
 
  subroutine set_scalar_ic_point_zone(s, base_value, zone_name, zone_value)
    type(field_t), intent(inout) :: s
    real(kind=rp), intent(in) :: base_value
    character(len=*), intent(in) :: zone_name
    real(kind=rp), intent(in) :: zone_value
 
    ! Internal variables
    character(len=LOG_SIZE) :: log_buf
    class(point_zone_t), pointer :: zone
    integer :: size
 
    call neko_log%message("Type       : point_zone")
    write (log_buf, '(A,ES12.6)') "Base value: ", base_value
    call neko_log%message(log_buf)
    call neko_log%message("Zone name : " // trim(zone_name))
    write (log_buf, '(A,ES12.6)') "Zone value: ", zone_value
    call neko_log%message(log_buf)
 
    size = s%dof%size()
    zone => neko_point_zone_registry%get_point_zone(trim(zone_name))
 
    call set_scalar_ic_uniform(s, base_value)
    call cfill_mask(s%x, zone_value, size, zone%mask%get(), zone%size)
 
  end subroutine set_scalar_ic_point_zone
 
  subroutine set_scalar_ic_fld(s, file_name, &
       interpolate, mesh_file_name, i, target_idx, global_interp_subdict)
    type(field_t), target, intent(inout) :: s
    character(len=*), intent(in) :: file_name
    logical, intent(in) :: interpolate
    character(len=*), intent(inout) :: mesh_file_name
    integer, intent(in) :: i
    integer, intent(in) :: target_idx
    type(json_file), intent(inout) :: global_interp_subdict
 
    character(len=LOG_SIZE) :: log_buf
    type(field_t), pointer :: ss
    type(field_list_t) :: s_tgt_list
 
    if (i .ne. target_idx) then
       write (log_buf, '(A,I0,A,I0)') "Loading scalar #", target_idx, &
            " into scalar #", i
       call neko_log%message(log_buf)
    end if
 
    ! use a pointer since import_fields needs a pointer as input
    ss => s
 
    ! Put ss in a field list of 1 element
    call s_tgt_list%init(1)
    call s_tgt_list%assign(1, ss)
 
    call import_fields(file_name, global_interp_subdict, mesh_file_name, &
         s_target_list = s_tgt_list, & ! The target field
         s_index_list = [target_idx], & ! Take values from target scalar
         interpolate = interpolate)
 
    call s_tgt_list%free()
 
    nullify(ss)
 
    ! If we are on GPU we need to move s back to the host
    ! since set_scalar_ic_common copies it again to the device.
    call s%copy_from(device_to_host, .true.)
 
  end subroutine set_scalar_ic_fld
 
end module scalar_ic