From 20abc9d147fbd0e7bc7d1d6fb69f557dd6740965 Mon Sep 17 00:00:00 2001 From: Jeff Candy Date: Thu, 24 Oct 2024 15:24:08 -0700 Subject: [PATCH] Streamlining global code (no SHEAR_METHOD=3) --- cgyro/src/cgyro_globalshear.F90 | 15 +++------------ cgyro/src/cgyro_make_profiles.F90 | 7 ++----- 2 files changed, 5 insertions(+), 17 deletions(-) diff --git a/cgyro/src/cgyro_globalshear.F90 b/cgyro/src/cgyro_globalshear.F90 index 632d0317e..19f376dce 100644 --- a/cgyro/src/cgyro_globalshear.F90 +++ b/cgyro/src/cgyro_globalshear.F90 @@ -28,8 +28,7 @@ subroutine cgyro_globalshear(ij) !$acc private(ivc,ir,l,iccj,j,ll,rl,llnt,h1,h2) & !$acc present(rhs(:,:,:,ij),omega_ss,omega_sbeta,field,cap_h_c,h_x,c_wave) #else -!$omp parallel do collapse(3) private(ivc,ir,l,iccj,j,ll,rl,llnt,h1,h2) & -!$omp firstprivate(shear_method) +!$omp parallel do collapse(3) private(ivc,ir,l,iccj,j,ll,rl,llnt,h1,h2) #endif do itor=nt1,nt2 do ivc=1,nv_loc @@ -49,11 +48,7 @@ subroutine cgyro_globalshear(ij) if ( (ir+ll) <= n_radial ) then ! ExB shear - if (shear_method == 3) then - h1 = omega_eb_base*itor*cap_h_c(iccj+llnt,ivc,itor) - else - h1 = omega_eb_base*itor*h_x(iccj+llnt,ivc,itor) - endif + 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) ! omega_star shear @@ -64,11 +59,7 @@ subroutine cgyro_globalshear(ij) if ( (ir-ll) >= 1 ) then ! ExB shear - if (shear_method == 3) then - h2 = omega_eb_base*itor*cap_h_c(iccj-llnt,ivc,itor) - else - h2 = omega_eb_base*itor*h_x(iccj-llnt,ivc,itor) - endif + 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) ! omega_star shear diff --git a/cgyro/src/cgyro_make_profiles.F90 b/cgyro/src/cgyro_make_profiles.F90 index 011bc2bfe..3b2a806f0 100644 --- a/cgyro/src/cgyro_make_profiles.F90 +++ b/cgyro/src/cgyro_make_profiles.F90 @@ -416,12 +416,9 @@ subroutine cgyro_make_profiles omega_eb_base = k_theta_base*length*gamma_e/(2*pi) select case (shear_method) case (1) - call cgyro_info('ExB shear: Hammett discrete shift (do not use for production simulations)') + call cgyro_info('ExB shear: Hammett discrete shift (WARNING: TESTING PURPOSES ONLY)') case (2) - call cgyro_info('ExB shear: Wavenumber advection (h)') - source_flag = 1 - case (3) - call cgyro_info('ExB shear: Wavenumber advection (H)') + call cgyro_info('ExB shear: Wavenumber advection') source_flag = 1 case default call cgyro_error('Unknown ExB shear method')