s_spline_interpolator_1d__init Subroutine

private subroutine s_spline_interpolator_1d__init(self, bspl, bc_xmin, bc_xmax)

Type Bound

sll_t_spline_interpolator_1d

Arguments

Type IntentOptional Attributes Name
class(sll_t_spline_interpolator_1d), intent(out) :: self

1D spline interpolator
class(sll_c_bsplines), intent(in), target :: bspl

B-splines (basis)
integer, intent(in) :: bc_xmin

boundary condition at xmin
integer, intent(in) :: bc_xmax

boundary condition at xmax

Calls

proc~~s_spline_interpolator_1d__init~~CallsGraph proc~s_spline_interpolator_1d__init sll_t_spline_interpolator_1d%s_spline_interpolator_1d__init factorize factorize proc~s_spline_interpolator_1d__init->factorize proc~s_build_system s_build_system proc~s_spline_interpolator_1d__init->proc~s_build_system proc~s_compute_interpolation_points_non_uniform s_compute_interpolation_points_non_uniform proc~s_spline_interpolator_1d__init->proc~s_compute_interpolation_points_non_uniform proc~s_compute_interpolation_points_uniform s_compute_interpolation_points_uniform proc~s_spline_interpolator_1d__init->proc~s_compute_interpolation_points_uniform proc~s_compute_num_diags_non_uniform s_compute_num_diags_non_uniform proc~s_spline_interpolator_1d__init->proc~s_compute_num_diags_non_uniform proc~s_compute_num_diags_uniform s_compute_num_diags_uniform proc~s_spline_interpolator_1d__init->proc~s_compute_num_diags_uniform proc~sll_s_spline_matrix_new sll_s_spline_matrix_new proc~s_spline_interpolator_1d__init->proc~sll_s_spline_matrix_new [] [] proc~s_build_system->[] eval_basis eval_basis proc~s_build_system->eval_basis eval_basis_and_n_derivs eval_basis_and_n_derivs proc~s_build_system->eval_basis_and_n_derivs set_element set_element proc~s_build_system->set_element find_cell find_cell proc~s_compute_num_diags_non_uniform->find_cell init init proc~sll_s_spline_matrix_new->init proc~sll_s_error_handler sll_s_error_handler proc~sll_s_spline_matrix_new->proc~sll_s_error_handler interface~c_abort~2 c_abort proc~sll_s_error_handler->interface~c_abort~2 proc~errout errout proc~sll_s_error_handler->proc~errout

Source Code

  subroutine s_spline_interpolator_1d__init( self, bspl, bc_xmin, bc_xmax )

    class(sll_t_spline_interpolator_1d), intent(  out) :: self
    class(sll_c_bsplines),       target, intent(in   ) :: bspl
    integer                            , intent(in   ) :: bc_xmin
    integer                            , intent(in   ) :: bc_xmax

    integer :: kl, ku
    character(len=32) :: matrix_type

    ! Sanity checks
    SLL_ASSERT( any( bc_xmin == allowed_bcs ) )
    SLL_ASSERT( any( bc_xmax == allowed_bcs ) )
    if (bspl % periodic) then
      SLL_ASSERT( bc_xmin == sll_p_periodic )
      SLL_ASSERT( bc_xmax == sll_p_periodic )
    end if

    ! Save pointer to B-splines
    ! (later needed to verify 1D spline input to 'compute_interpolant')
    self % bspl => bspl

    ! Save other data
    self % bc_xmin  = bc_xmin
    self % bc_xmax  = bc_xmax
    self % nbc_xmin = merge ( bspl%degree/2, 0, bc_xmin == sll_p_hermite )
    self % nbc_xmax = merge ( bspl%degree/2, 0, bc_xmax == sll_p_hermite )
    self % odd      = modulo( bspl%degree  , 2 )

    ! Save average cell size for normalization of derivatives
    self%dx = (bspl%xmax - bspl%xmin) / bspl%ncells

    ! Compute interpolation points and number of diagonals in linear system
    if (bspl % uniform) then
      call s_compute_interpolation_points_uniform( self, self%tau )
      call s_compute_num_diags_uniform( self, kl, ku )
    else
      call s_compute_interpolation_points_non_uniform( self, self%tau )
      call s_compute_num_diags_non_uniform( self, kl, ku )
    end if

    ! Special case: linear spline
    ! No need for matrix assembly
    if (self % bspl % degree == 1) return

    ! Determine matrix storage type and allocate matrix
    associate( nbasis => self % bspl % nbasis )

      if (bc_xmin == sll_p_periodic) then
        if (kl+1+ku >= nbasis) then
          matrix_type = "dense"
        else
          matrix_type = "periodic_banded"
        end if
      else
        matrix_type = "banded"
      end if

      call sll_s_spline_matrix_new( self%matrix, matrix_type, nbasis, kl, ku )

    end associate ! nbasis

    ! Fill in entries A_ij of linear system A*x = b for interpolation
    call s_build_system( self, self%matrix )

    ! Factorize matrix A to speed-up solution for any given b
    call self % matrix % factorize()

  end subroutine s_spline_interpolator_1d__init