Skip to content

Commit

Permalink
SBETA testing parameters
Browse files Browse the repository at this point in the history
jcandy committed Jan 2, 2025

Verified

This commit was signed with the committer’s verified signature. The key has expired.
1 parent 582ffe3 commit da4dd74
Showing 5 changed files with 17 additions and 6 deletions.
3 changes: 3 additions & 0 deletions cgyro/bin/cgyro_parse.py
Original file line number Diff line number Diff line change
@@ -153,6 +153,9 @@
x.add('NU_EE_SCALE','1.0')
x.add('ZF_SCALE','1.0')

x.add('SBETA_CONST_FLAG','0')
x.add('SBETA_H','0.0')

# Perform the parsing
x.read_input('input.cgyro')

4 changes: 3 additions & 1 deletion cgyro/src/cgyro_globals.F90
Original file line number Diff line number Diff line change
@@ -14,7 +14,9 @@ module cgyro_globals
#endif

use, intrinsic :: iso_fortran_env


integer :: sbeta_const_flag
real :: sbeta_h

!---------------------------------------------------------------
! Input parameters:
6 changes: 4 additions & 2 deletions cgyro/src/cgyro_globalshear.F90
Original file line number Diff line number Diff line change
@@ -50,7 +50,8 @@ subroutine cgyro_globalshear(ij)
! ExB shear
h1 = omega_eb_base*itor*h_x(iccj+llnt,ivc,itor)
! beta_star shear
h1 = h1+omega_sbeta(iccj+llnt,ivc,itor)*cap_h_c(iccj+llnt,ivc,itor)
h1 = h1+omega_sbeta(iccj+llnt,ivc,itor)*( &
(1-sbeta_h)*cap_h_c(iccj+llnt,ivc,itor)+sbeta_h*h_x(iccj+llnt,ivc,itor))
! omega_star shear
h1 = h1+sum(omega_ss(:,iccj+llnt,ivc,itor)*field(:,iccj+llnt,itor))
else
@@ -61,7 +62,8 @@ subroutine cgyro_globalshear(ij)
! ExB shear
h2 = omega_eb_base*itor*h_x(iccj-llnt,ivc,itor)
! beta_star shear
h2 = h2+omega_sbeta(iccj-llnt,ivc,itor)*cap_h_c(iccj-llnt,ivc,itor)
h2 = h2+omega_sbeta(iccj-llnt,ivc,itor)*( &
(1-sbeta_h)*cap_h_c(iccj-llnt,ivc,itor)+sbeta_h*h_x(iccj-llnt,ivc,itor))
! omega_star shear
h2 = h2+sum(omega_ss(:,iccj-llnt,ivc,itor)*field(:,iccj-llnt,itor))
else
6 changes: 3 additions & 3 deletions cgyro/src/cgyro_init_arrays.F90
Original file line number Diff line number Diff line change
@@ -465,10 +465,10 @@ subroutine cgyro_init_arrays
sm = sdlnndr(is)+sdlntdr(is)*(energy(ie)-1.5)

! generalized beta/drift shear (acts on H)
if (nup_alpha == 3) then
sb = -sbeta(is)*energy(ie)*xi(ix)**2/bmag(it)**3
else
if (sbeta_const_flag == 1) then
sb = -sbeta(is)
else
sb = -sbeta(is)*energy(ie)*xi(ix)**2/bmag(it)**3
endif

arg = k_theta_base*itor*length/(2*pi)
4 changes: 4 additions & 0 deletions cgyro/src/cgyro_read_input.f90
Original file line number Diff line number Diff line change
@@ -143,6 +143,10 @@ subroutine cgyro_read_input
call cgyro_readbc_real(nu_ee_scale)
call cgyro_readbc_real(zf_scale)

call cgyro_readbc_int(sbeta_const_flag,'SBETA_CONST_FLAG')
call cgyro_readbc_real(sbeta_h)


if (i_proc == 0) close(1)

end subroutine cgyro_read_input

0 comments on commit da4dd74

Please sign in to comment.