coord_transform/m_coord_transform_toroidal.f90
| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | !> @brief Module for the toroidal coordinate transformation | ||
| 2 | module m_coord_transform_toroidal | ||
| 3 | use m_common, only: wp, pi | ||
| 4 | use m_coord_transform_abstract, only: CoordTransformAbstract | ||
| 5 | |||
| 6 | implicit none | ||
| 7 | |||
| 8 | private | ||
| 9 | public :: ToroidalTransform | ||
| 10 | public :: toroidal_transform, toroidal_inverse_transform, toroidal_jacobian | ||
| 11 | public :: toroidal_jacobian_matrix, toroidal_jacobian_matrix_inv, toroidal_G_matrix, toroidal_G_matrix_inv | ||
| 12 | |||
| 13 | !> @brief Toroidal coordinates given by | ||
| 14 | !> \f{eqnarray*}{ | ||
| 15 | !> x &= (R0 + R1 xp \sin(2 \pi yp)) \cos(2 \pi zp)\\ | ||
| 16 | !> y &= (R0 + R1 xp \sin(2 \pi yp)) \sin(2 \pi zp)\\ | ||
| 17 | !> z &= R1 xp \cos(2 \pi yp) | ||
| 18 | !> \f} | ||
| 19 | type, extends(CoordTransformAbstract) :: ToroidalTransform | ||
| 20 | !> @brief The major radius of the toroidal coordinate system | ||
| 21 | real(wp) :: major_radius | ||
| 22 | |||
| 23 | !> @brief The minor radius of the toroidal coordinate system | ||
| 24 | real(wp) :: minor_radius | ||
| 25 | end type ToroidalTransform | ||
| 26 | |||
| 27 | interface ToroidalTransform | ||
| 28 | module procedure init_ToroidalTransform | ||
| 29 | end interface ToroidalTransform | ||
| 30 | |||
| 31 | interface | ||
| 32 | |||
| 33 | !> @brief Create a new ToroidalTransform object | ||
| 34 | !> @param[in] major_radius The major radius of the toroidal coordinate system | ||
| 35 | !> @param[in] minor_radius The minor radius of the toroidal coordinate system | ||
| 36 | !> @return A new ToroidalTransform object | ||
| 37 | module function init_ToroidalTransform(major_radius, minor_radius) result(this) | ||
| 38 | real(wp), intent(in) :: major_radius, minor_radius | ||
| 39 | type(ToroidalTransform) :: this | ||
| 40 | end function init_ToroidalTransform | ||
| 41 | |||
| 42 | !> @brief The i-th component of the transformation | ||
| 43 | !> @param[in] this The ToroidalTransform object | ||
| 44 | !> @param[in] i The index of the component | ||
| 45 | !> @param[in] xp The radial coordinate | ||
| 46 | !> @param[in] yp The poloidal angle | ||
| 47 | !> @param[in] zp The toroidal angle | ||
| 48 | !> @return The i-th component of the transformation | ||
| 49 | module pure real(wp) function toroidal_transform(this, i, xp, yp, zp) result(ans) | ||
| 50 | !$acc routine seq | ||
| 51 | type(ToroidalTransform), intent(in) :: this | ||
| 52 | integer, intent(in) :: i | ||
| 53 | real(wp), intent(in) :: xp, yp, zp | ||
| 54 | end function toroidal_transform | ||
| 55 | |||
| 56 | !> @brief The i-th component of the inverse transformation | ||
| 57 | !> @param[in] this The ToroidalTransform object | ||
| 58 | !> @param[in] i The index of the component | ||
| 59 | !> @param[in] x The x coordinate | ||
| 60 | !> @param[in] y The y coordinate | ||
| 61 | !> @param[in] z The z coordinate | ||
| 62 | !> @return The i-th coordinate | ||
| 63 | module pure real(wp) function toroidal_inverse_transform(this, i, x, y, z) result(ans) | ||
| 64 | !$acc routine seq | ||
| 65 | type(ToroidalTransform), intent(in) :: this | ||
| 66 | integer, intent(in) :: i | ||
| 67 | real(wp), intent(in) :: x, y, z | ||
| 68 | end function toroidal_inverse_transform | ||
| 69 | |||
| 70 | !> @brief The Jacobian determinant of the transformation | ||
| 71 | !> @param[in] this The ToroidalTransform object | ||
| 72 | !> @param[in] xp The radial coordinate | ||
| 73 | !> @param[in] yp The poloidal angle | ||
| 74 | !> @param[in] zp The toroidal angle | ||
| 75 | !> @return The Jacobian determinant | ||
| 76 | module pure real(wp) function toroidal_jacobian(this, xp, yp, zp) result(ans) | ||
| 77 | !$acc routine seq | ||
| 78 | type(ToroidalTransform), intent(in) :: this | ||
| 79 | real(wp), intent(in) :: xp, yp, zp | ||
| 80 | end function toroidal_jacobian | ||
| 81 | |||
| 82 | !> @brief The (i, j)-th component of the Jacobian matrix | ||
| 83 | !> @param[in] this The ToroidalTransform object | ||
| 84 | !> @param[in] i The row index | ||
| 85 | !> @param[in] j The column index | ||
| 86 | !> @param[in] xp The radial coordinate | ||
| 87 | !> @param[in] yp The poloidal angle | ||
| 88 | !> @param[in] zp The toroidal angle | ||
| 89 | !> @return The (i, j)-th component of the Jacobian matrix | ||
| 90 | module pure real(wp) function toroidal_jacobian_matrix(this, i, j, xp, yp, zp) result(ans) | ||
| 91 | !$acc routine seq | ||
| 92 | type(ToroidalTransform), intent(in) :: this | ||
| 93 | integer, intent(in) :: i, j | ||
| 94 | real(wp), intent(in) :: xp, yp, zp | ||
| 95 | end function toroidal_jacobian_matrix | ||
| 96 | |||
| 97 | !> @brief The (i, j)-th component of the inverse of the Jacobian matrix | ||
| 98 | !> @param[in] this The ToroidalTransform object | ||
| 99 | !> @param[in] i The row index | ||
| 100 | !> @param[in] j The column index | ||
| 101 | !> @param[in] xp The radial coordinate | ||
| 102 | !> @param[in] yp The poloidal angle | ||
| 103 | !> @param[in] zp The toroidal angle | ||
| 104 | !> @return The (i, j)-th component of the inverse of the Jacobian matrix | ||
| 105 | module pure real(wp) function toroidal_jacobian_matrix_inv(this, i, j, xp, yp, zp) result(ans) | ||
| 106 | !$acc routine seq | ||
| 107 | type(ToroidalTransform), intent(in) :: this | ||
| 108 | integer, intent(in) :: i, j | ||
| 109 | real(wp), intent(in) :: xp, yp, zp | ||
| 110 | end function toroidal_jacobian_matrix_inv | ||
| 111 | |||
| 112 | !> @brief The (i, j)-th component of the metric tensor | ||
| 113 | !> @param[in] this The ToroidalTransform object | ||
| 114 | !> @param[in] i The row index | ||
| 115 | !> @param[in] j The column index | ||
| 116 | !> @param[in] xp The radial coordinate | ||
| 117 | !> @param[in] yp The poloidal angle | ||
| 118 | !> @param[in] zp The toroidal angle | ||
| 119 | !> @return The (i, j)-th component of the metric tensor | ||
| 120 | module pure real(wp) function toroidal_G_matrix(this, i, j, xp, yp, zp) result(ans) | ||
| 121 | !$acc routine seq | ||
| 122 | type(ToroidalTransform), intent(in) :: this | ||
| 123 | integer, intent(in) :: i, j | ||
| 124 | real(wp), intent(in) :: xp, yp, zp | ||
| 125 | end function toroidal_G_matrix | ||
| 126 | |||
| 127 | !> @brief The (i, j)-th component of the inverse metric tensor | ||
| 128 | !> @param[in] this The ToroidalTransform object | ||
| 129 | !> @param[in] i The row index | ||
| 130 | !> @param[in] j The column index | ||
| 131 | !> @param[in] xp The radial coordinate | ||
| 132 | !> @param[in] yp The poloidal angle | ||
| 133 | !> @param[in] zp The toroidal angle | ||
| 134 | !> @return The (i, j)-th component of the inverse metric tensor | ||
| 135 | module pure real(wp) function toroidal_G_matrix_inv(this, i, j, xp, yp, zp) result(ans) | ||
| 136 | !$acc routine seq | ||
| 137 | type(ToroidalTransform), intent(in) :: this | ||
| 138 | integer, intent(in) :: i, j | ||
| 139 | real(wp), intent(in) :: xp, yp, zp | ||
| 140 | end function toroidal_G_matrix_inv | ||
| 141 | |||
| 142 | end interface | ||
| 143 | |||
| 144 | ✗ | end module m_coord_transform_toroidal | |
| 145 |