Skip to content
Snippets Groups Projects
Commit 8dcc64de authored by Antoine Hoffmann's avatar Antoine Hoffmann
Browse files

comment s-alpha geometry (not ready)

parent f1f70613
No related branches found
No related tags found
No related merge requests found
...@@ -77,67 +77,67 @@ contains ...@@ -77,67 +77,67 @@ contains
! !
!-------------------------------------------------------------------------------- !--------------------------------------------------------------------------------
! UNUSED ! UNUSED
subroutine eval_salpha_geometry ! subroutine eval_salpha_geometry
! evaluate s-alpha geometry model ! ! evaluate s-alpha geometry model
implicit none ! implicit none
REAL(dp) :: z, kx, ky ! REAL(dp) :: z, kx, ky
!
zloop: DO iz = izs,ize ! zloop: DO iz = izs,ize
z = zarray(iz) ! z = zarray(iz)
gxx(iz) = 1._dp ! gxx(iz) = 1._dp
gxy(iz) = shear*z ! gxy(iz) = shear*z
gyy(iz) = 1._dp + (shear*z)**2 ! gyy(iz) = 1._dp + (shear*z)**2
gyz(iz) = 1._dp/eps ! gyz(iz) = 1._dp/eps
gxz(iz) = 0._dp ! gxz(iz) = 0._dp
!
! Relative strengh of radius ! ! Relative strengh of radius
hatR(iz) = 1._dp + eps*cos(z) ! hatR(iz) = 1._dp + eps*cos(z)
!
! Jacobian ! ! Jacobian
Jacobian(iz) = q0*hatR(iz) ! Jacobian(iz) = q0*hatR(iz)
!
! Relative strengh of modulus of B ! ! Relative strengh of modulus of B
hatB(iz) = 1._dp / hatR( iz) ! hatB(iz) = 1._dp / hatR( iz)
!
! Derivative of the magnetic field strenght ! ! Derivative of the magnetic field strenght
gradxB(iz) = - cos( z) / hatR(iz) ! gradxB(iz) = - cos( z) / hatR(iz)
gradzB(iz) = eps * sin(z) ! gradzB(iz) = eps * sin(z)
!
! Gemoetrical coefficients for the curvature operator ! ! Gemoetrical coefficients for the curvature operator
! Note: Gamma2 and Gamma3 are obtained directly form Gamma1 in the expression of the curvature operator implemented here ! ! Note: Gamma2 and Gamma3 are obtained directly form Gamma1 in the expression of the curvature operator implemented here
! ! !
Gamma1(iz) = gxy(iz) * gxy(iz) - gxx(iz) * gyy(iz) ! Gamma1(iz) = gxy(iz) * gxy(iz) - gxx(iz) * gyy(iz)
Gamma2(iz) = gxz(iz) * gxy(iz) - gxx(iz) * gyz(iz) ! Gamma2(iz) = gxz(iz) * gxy(iz) - gxx(iz) * gyz(iz)
Gamma3(iz) = gxz(iz) * gyy(iz) - gxy(iz) * gyz(iz) ! Gamma3(iz) = gxz(iz) * gyy(iz) - gxy(iz) * gyz(iz)
!
! Curvature operator ! ! Curvature operator
DO iky = ikys, ikye ! DO iky = ikys, ikye
ky = kyarray(iky) ! ky = kyarray(iky)
DO ikx= ikxl, ikxr ! DO ikx= ikxl, ikxr
kx = kxarray(ikx,iky) ! kx = kxarray(ikx,iky)
Cxy(ikx, iky, iz) = (-sin(z)*kx - (cos(z) + shear* z* sin(z))*ky) * hatB(iz) ! .. multiply by hatB to cancel the 1/ hatB factor in moments_eqs_rhs.f90 routine ! Cxy(ikx, iky, iz) = (-sin(z)*kx - (cos(z) + shear* z* sin(z))*ky) * hatB(iz) ! .. multiply by hatB to cancel the 1/ hatB factor in moments_eqs_rhs.f90 routine
ENDDO ! ENDDO
ENDDO ! ENDDO
!
! coefficient in the front of parallel derivative ! ! coefficient in the front of parallel derivative
gradz_coeff(iz) = 1._dp / Jacobian(iz) / hatB(iz) ! gradz_coeff(iz) = 1._dp / Jacobian(iz) / hatB(iz)
ENDDO ! ENDDO
!
! Evaluate perpendicular wavenumber ! ! Evaluate perpendicular wavenumber
! k_\perp^2 = g^{xx} k_x^2 + 2 g^{xy}k_x k_y + k_y^2 g^{yy} ! ! k_\perp^2 = g^{xx} k_x^2 + 2 g^{xy}k_x k_y + k_y^2 g^{yy}
! normalized to rhos_ ! ! normalized to rhos_
DO iky = ikys, ikye ! DO iky = ikys, ikye
ky = kyarray(iky) ! ky = kyarray(iky)
DO ikx = ikxl, ikxr ! DO ikx = ikxl, ikxr
kx = kxarray(ikx,iky) ! kx = kxarray(ikx,iky)
kperp_array(ikx, iky, iz) = sqrt( gxx(iz)*kx**2 + 2._dp* gxy(iz) * kx*ky + gyy(iz)* ky**2) !! / hatB( iz) ! there is a factor 1/B from the normalization; important to match GENE ! kperp_array(ikx, iky, iz) = sqrt( gxx(iz)*kx**2 + 2._dp* gxy(iz) * kx*ky + gyy(iz)* ky**2) !! / hatB( iz) ! there is a factor 1/B from the normalization; important to match GENE
ENDDO ! ENDDO
ENDDO ! ENDDO
ENDDO zloop ! ENDDO zloop
! ! !
IF( me .eq. 0 ) PRINT*, 'max kperp = ', maxval( kperp_array) ! IF( me .eq. 0 ) PRINT*, 'max kperp = ', maxval( kperp_array)
!
END SUBROUTINE eval_salpha_geometry ! END SUBROUTINE eval_salpha_geometry
! !
!-------------------------------------------------------------------------------- !--------------------------------------------------------------------------------
! !
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment