7 !! Definition of the components g11 and g22 of the metric tensor of the
8 !! applied coordinates.
12 !! Copyright 2009-2013 Ralf Greve
16 !! This file is part of SICOPOLIS.
18 !! SICOPOLIS is free software: you can redistribute it and/or modify
19 !! it under the terms of the GNU General Public License as published by
20 !! the Free Software Foundation, either version 3 of the License, or
21 !! (at your option) any later version.
23 !! SICOPOLIS is distributed in the hope that it will be useful,
24 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
25 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
26 !! GNU General Public License for more details.
28 !! You should have received a copy of the GNU General Public License
29 !! along with SICOPOLIS. If not, see <http://www.gnu.org/licenses/>.
31 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
34 !> Definition of the components g11 and g22 of the metric tensor of the
35 !! applied coordinates.
36 !<------------------------------------------------------------------------------
45 real(dp) :: g11_g(0:jmax,0:imax), g22_g(0:jmax,0:imax), &
46 g11_sgx(0:jmax,0:imax), g11_sgy(0:jmax,0:imax), &
47 g22_sgx(0:jmax,0:imax), g22_sgy(0:jmax,0:imax)
48 real(dp) :: xhelp, yhelp
71 k = (cos(0.25_dp*pi-0.5_dp*phi0))**2
72 else if (phi0.lt.(-eps))
then
73 k = (cos(0.25_dp*pi+0.5_dp*phi0))**2
75 stop
' metric: PHI0 must be different from zero!'
80 g11_g(j,i) = 1.0_dp/ &
81 ( k**2*(1.0_dp+(xi(i)**2+eta(j)**2)/(2.0_dp*r*k)**2)**2 )
82 g22_g(j,i) = g11_g(j,i)
90 xhelp = 0.5_dp*(xi(i)+xi(i+1))
91 g11_sgx(j,i) = 1.0_dp/ &
92 ( k**2*(1.0_dp+(xhelp**2+eta(j)**2)/(2.0_dp*r*k)**2)**2 )
93 g22_sgx(j,i) = g11_sgx(j,i)
99 yhelp = 0.5_dp*(eta(j)+eta(j+1))
100 g22_sgy(j,i) = 1.0_dp/ &
101 ( k**2*(1.0_dp+(xi(i)**2+yhelp**2)/(2.0_dp*r*k)**2)**2 )
102 g11_sgy(j,i) = g22_sgy(j,i)
112 g11_g(j,i) = r**2*(cos(eta(j)))**2
121 g11_sgx(j,i) = r**2*(cos(eta(j)))**2
128 yhelp = 0.5_dp*(eta(j)+eta(j+1))
130 g11_sgy(j,i) = r**2*(cos(yhelp))**2
141 sq_g11_g(j,i) = sqrt(g11_g(j,i))
142 sq_g22_g(j,i) = sqrt(g22_g(j,i))
143 insq_g11_g(j,i) = 1.0_dp/sq_g11_g(j,i)
144 insq_g22_g(j,i) = 1.0_dp/sq_g22_g(j,i)
150 sq_g11_sgx(j,i) = sqrt(g11_sgx(j,i))
151 sq_g22_sgx(j,i) = sqrt(g22_sgx(j,i))
152 insq_g11_sgx(j,i) = 1.0_dp/sq_g11_sgx(j,i)
158 sq_g22_sgy(j,i) = sqrt(g22_sgy(j,i))
159 sq_g11_sgy(j,i) = sqrt(g11_sgy(j,i))
160 insq_g22_sgy(j,i) = 1.0_dp/sq_g22_sgy(j,i)