metric.F90

Go to the documentation of this file.

!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!
!  Subroutine :  m e t r i c
!
!> @file
!!
!! Definition of the components g11 and g22 of the metric tensor of the
!! applied coordinates.
!!
!! @section Copyright
!!
!! Copyright 2009-2016 Ralf Greve
!!
!! @section License
!!
!! This file is part of SICOPOLIS.
!!
!! SICOPOLIS 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 3 of the License, or
!! (at your option) any later version.
!!
!! SICOPOLIS 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 SICOPOLIS.  If not, see <http://www.gnu.org/licenses/>.
!<
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

!-------------------------------------------------------------------------------
!> Definition of the components g11 and g22 of the metric tensor of the
!! applied coordinates.
!<------------------------------------------------------------------------------
subroutine metric()

use sico_types
use sico_variables
use sico_vars

implicit none
integer(i4b) :: i, j
real(dp) :: k
real(dp) :: g11_g(0:jmax,0:imax), g22_g(0:jmax,0:imax), &
            g11_sgx(0:jmax,0:imax), g11_sgy(0:jmax,0:imax), &
            g22_sgx(0:jmax,0:imax), g22_sgy(0:jmax,0:imax)
real(dp) :: xhelp, yhelp

#if GRID==0

!-------- Components g11, g22 on the grid points (_g) and between
!         the grid points (_sg) --------

do i=0, imax
do j=0, jmax
   g11_g(j,i) = 1
   g22_g(j,i) = 1
   g11_sgx(j,i) = 1
   g11_sgy(j,i) = 1
   g22_sgx(j,i) = 1
   g22_sgy(j,i) = 1
end do
end do

#elif GRID==1

!-------- Components g11, g22 on the grid points (_g) --------

if (phi0 > eps) then   ! for northern hemisphere
   k = (cos(0.25_dp*pi-0.5_dp*phi0))**2
else if (phi0 < (-eps)) then   ! for southern hemisphere
   k = (cos(0.25_dp*pi+0.5_dp*phi0))**2
else
   stop ' metric: PHI0 must be different from zero!'
end if

do i=0, imax
do j=0, jmax
   g11_g(j,i) = 1.0_dp/ &
      ( k**2*(1.0_dp+(xi(i)**2+eta(j)**2)/(2.0_dp*r*k)**2)**2 )
   g22_g(j,i) = g11_g(j,i)
end do
end do

!-------- Components g11, g22 between the grid points (_sg) --------

do i=0, imax-1
do j=0, jmax
   xhelp = 0.5_dp*(xi(i)+xi(i+1))
   g11_sgx(j,i) = 1.0_dp/ &
      ( k**2*(1.0_dp+(xhelp**2+eta(j)**2)/(2.0_dp*r*k)**2)**2 )
   g22_sgx(j,i) = g11_sgx(j,i)
end do
end do

do i=0, imax
do j=0, jmax-1
   yhelp = 0.5_dp*(eta(j)+eta(j+1))
   g22_sgy(j,i) = 1.0_dp/ &
      ( k**2*(1.0_dp+(xi(i)**2+yhelp**2)/(2.0_dp*r*k)**2)**2 )
   g11_sgy(j,i) = g22_sgy(j,i)
end do
end do

#elif GRID==2

!-------- Components g11, g22 on the grid points (_g) --------

do i=0, imax
do j=0, jmax
   g11_g(j,i) = r**2*(cos(eta(j)))**2
   g22_g(j,i) = r**2
end do
end do

!-------- Components g11, g22 between the grid points (_sg) --------

do i=0, imax-1
do j=0, jmax
   g11_sgx(j,i) = r**2*(cos(eta(j)))**2
   g22_sgx(j,i) = r**2
end do
end do

do i=0, imax
do j=0, jmax-1
   yhelp = 0.5_dp*(eta(j)+eta(j+1))
   g22_sgy(j,i) = r**2
   g11_sgy(j,i) = r**2*(cos(yhelp))**2
end do
end do

#endif

!-------- Square roots (sq_) and inverse square roots (insq_) of
!         g11 and g22 --------

do i=0, imax
do j=0, jmax
   sq_g11_g(j,i)   = sqrt(g11_g(j,i))
   sq_g22_g(j,i)   = sqrt(g22_g(j,i))
   insq_g11_g(j,i) = 1.0_dp/sq_g11_g(j,i)
   insq_g22_g(j,i) = 1.0_dp/sq_g22_g(j,i)
end do
end do

do i=0, imax-1
do j=0, jmax
   sq_g11_sgx(j,i)   = sqrt(g11_sgx(j,i))
   sq_g22_sgx(j,i)   = sqrt(g22_sgx(j,i))
   insq_g11_sgx(j,i) = 1.0_dp/sq_g11_sgx(j,i)
end do
end do

do i=0, imax
do j=0, jmax-1
   sq_g22_sgy(j,i)   = sqrt(g22_sgy(j,i))
   sq_g11_sgy(j,i)   = sqrt(g11_sgy(j,i))
   insq_g22_sgy(j,i) = 1.0_dp/sq_g22_sgy(j,i)
end do
end do

end subroutine metric
!