coordinate_system.F90

Go to the documentation of this file.

!! Copyright (C) 2021 M. Oliveira
!! Copyright (C) 2025 S. Ohlmann
!!
!! This program is free software; you can redistribute it and/or modify
!! it under the terms of the GNU General Public License as published by
!! the Free Software Foundation; either version 2, or (at your option)
!! any later version.
!!
!! This program is distributed in the hope that it will be useful,
!! but WITHOUT ANY WARRANTY; without even the implied warranty of
!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!! GNU General Public License for more details.
!!
!! You should have received a copy of the GNU General Public License
!! along with this program; if not, write to the Free Software
!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
!! 02110-1301, USA.
!!
 
#include "global.h"
 
module coordinate_system_oct_m
  use, intrinsic :: iso_fortran_env
  use global_oct_m
  use lalg_adv_oct_m
  use messages_oct_m
  use namespace_oct_m
  implicit none
 
  private
  public :: coordinate_system_t
 
  type, abstract :: coordinate_system_t
    logical :: local_basis
    logical :: orthogonal
    integer :: dim
    real(real64) :: min_mesh_scaling_product
  contains
    generic   :: vector_from_cartesian => dvector_from_cartesian, zvector_from_cartesian
    
    procedure :: dvector_from_cartesian => dcoordinate_system_vector_from_cartesian
    
    procedure :: zvector_from_cartesian => zcoordinate_system_vector_from_cartesian
    
    generic   :: covector_to_cartesian => dcovector_to_cartesian, zcovector_to_cartesian
    
    procedure :: dcovector_to_cartesian => dcoordinate_system_covector_to_cartesian
    
    procedure :: zcovector_to_cartesian => zcoordinate_system_covector_to_cartesian
    
    procedure :: metric => coordinate_system_metric
    
    procedure :: metric_inverse => coordinate_system_metric_inverse
    
    procedure :: jacobian_determinant => coordinate_system_jacobian_determinant
    
    procedure :: trace_hessian => coordinate_system_trace_hessian
    
    procedure(coordinate_system_to_cartesian),           deferred :: to_cartesian
    
    procedure(coordinate_system_from_cartesian),         deferred :: from_cartesian
    
    procedure(coordinate_system_write_info),             deferred :: write_info
    
    procedure(coordinates_surface_element),              deferred :: surface_element
    
    procedure(coordinate_system_jacobian),               deferred :: jacobian
    
    procedure(coordinate_system_jacobian_inverse),       deferred :: jacobian_inverse
    
  end type coordinate_system_t
 
  abstract interface
    ! ---------------------------------------------------------
    function coordinate_system_to_cartesian(this, chi) result(xx)
      import coordinate_system_t
      import real64
      class(coordinate_system_t), target, intent(in)  :: this
      real(real64),                       intent(in)  :: chi(:)
      real(real64) :: xx(1:this%dim)
    end function coordinate_system_to_cartesian
 
    ! ---------------------------------------------------------
    function coordinate_system_from_cartesian(this, xx) result(chi)
      import coordinate_system_t
      import real64
      class(coordinate_system_t), target, intent(in)  :: this
      real(real64),                       intent(in)  :: xx(:)
      real(real64) :: chi(1:this%dim)
    end function coordinate_system_from_cartesian
 
    ! ---------------------------------------------------------
    subroutine coordinate_system_write_info(this, iunit, namespace)
      import coordinate_system_t
      import namespace_t
      class(coordinate_system_t),           intent(in) :: this
      integer,                    optional, intent(in) :: iunit
      type(namespace_t),          optional, intent(in) :: namespace
    end subroutine coordinate_system_write_info
 
    ! ---------------------------------------------------------
    real(real64) function coordinates_surface_element(this, idir) result(ds)
      import coordinate_system_t
      import real64
      class(coordinate_system_t), intent(in) :: this
      integer,                    intent(in) :: idir
    end function coordinates_surface_element
 
    ! ---------------------------------------------------------
    function coordinate_system_jacobian(this, chi) result(jacobian)
      import coordinate_system_t
      import real64
      class(coordinate_system_t), intent(in)  :: this
      real(real64),               intent(in)  :: chi(:)
      real(real64) :: jacobian(1:this%dim, 1:this%dim)
    end function coordinate_system_jacobian
 
    ! ---------------------------------------------------------
    function coordinate_system_jacobian_inverse(this, chi) result(jacobian_inverse)
      import coordinate_system_t
      import real64
      class(coordinate_system_t), intent(in)  :: this
      real(real64),               intent(in)  :: chi(:)
      real(real64) :: jacobian_inverse(1:this%dim, 1:this%dim)
    end function coordinate_system_jacobian_inverse
  end interface
 
contains
  ! ---------------------------------------------------------
  function coordinate_system_metric(this, chi) result(metric)
    class(coordinate_system_t), intent(in)  :: this
    real(real64),               intent(in)  :: chi(:)
    real(real64) :: metric(1:this%dim, 1:this%dim)
 
    real(real64) :: jacobian(1:this%dim, 1:this%dim)
 
    jacobian = this%jacobian(chi)
    metric = matmul(transpose(jacobian), jacobian)
  end function coordinate_system_metric
 
  ! ---------------------------------------------------------
  function coordinate_system_metric_inverse(this, chi) result(metric_inverse)
    class(coordinate_system_t), intent(in)  :: this
    real(real64),               intent(in)  :: chi(:)
    real(real64) :: metric_inverse(1:this%dim, 1:this%dim)
 
    real(real64) :: jacobian_inverse(1:this%dim, 1:this%dim)
 
    jacobian_inverse = this%jacobian_inverse(chi)
    metric_inverse = matmul(jacobian_inverse, transpose(jacobian_inverse))
  end function coordinate_system_metric_inverse
 
  ! ---------------------------------------------------------
  function coordinate_system_jacobian_determinant(this, chi) result(jacobian_determinant)
    class(coordinate_system_t), intent(in)  :: this
    real(real64),               intent(in)  :: chi(:)
    real(real64) :: jacobian_determinant
 
    real(real64) :: jacobian(1:this%dim, 1:this%dim)
 
    jacobian = this%jacobian(chi)
    jacobian_determinant = lalg_determinant(this%dim, jacobian, preserve_mat=.false.)
  end function coordinate_system_jacobian_determinant
 
  ! ---------------------------------------------------------
  function coordinate_system_trace_hessian(this, chi) result(trace_hessian)
    class(coordinate_system_t), intent(in)  :: this
    real(real64),               intent(in)  :: chi(:)
    real(real64) :: trace_hessian(1:this%dim)
 
    ! the Hessian is non-zero only if we have a local basis, i.e. for general curvilinear coordinates
    if (this%local_basis) then
      ! needs to be implemented in derived classes
      call messages_not_implemented("Trace of Hessian for the coordinate system in use")
    else
      trace_hessian(:) = m_zero
    end if
  end function coordinate_system_trace_hessian
 
#include "undef.F90"
#include "real.F90"
#include "coordinate_system_inc.F90"
 
#include "undef.F90"
#include "complex.F90"
#include "coordinate_system_inc.F90"
 
end module coordinate_system_oct_m
 
!! Local Variables:
!! mode: f90
!! coding: utf-8
!! End: