advance_gs2_gryfx Subroutine

public subroutine advance_gs2_gryfx(dens_ky0, upar_ky0, tpar_ky0, tprp_ky0, qpar_ky0, qprp_ky0, phi_ky0, first_half_step)

FIXME : Add documentation

Arguments

Type IntentOptional Attributes Name
complex(kind=real32), intent(inout), dimension (:) :: dens_ky0
complex(kind=real32), intent(inout), dimension (:) :: upar_ky0
complex(kind=real32), intent(inout), dimension (:) :: tpar_ky0
complex(kind=real32), intent(inout), dimension (:) :: tprp_ky0
complex(kind=real32), intent(inout), dimension (:) :: qpar_ky0
complex(kind=real32), intent(inout), dimension (:) :: qprp_ky0
complex(kind=real32), intent(inout), dimension (:) :: phi_ky0
logical, intent(in) :: first_half_step

Contents

Source Code


Source Code

  subroutine advance_gs2_gryfx(dens_ky0, upar_ky0, tpar_ky0, tprp_ky0, &
       qpar_ky0, qprp_ky0, phi_ky0, first_half_step)

    use nonlinear_terms, only: gryfx_zonal
    use gs2_main, only: evolve_equations
    use dist_fn, only: getmoms_gryfx_dist 
    use mp, only: broadcast
    use unit_tests, only: debug_message

    implicit none
    complex(real32), dimension (:), intent (inout) :: &
         dens_ky0, upar_ky0, tpar_ky0, tprp_ky0, qpar_ky0, qprp_ky0
    complex(real32), dimension (:), intent (inout) :: phi_ky0
    logical, intent (in) :: first_half_step

    !state%istep_end = istep
    state%print_times = .false.

    call debug_message(verb+1, 'gs2_gryfx_zonal::advance_gs2_gryfx starting')

    ! only proc0 knows values of arrays, so broadcast
    call broadcast(dens_ky0)
    call broadcast(upar_ky0)
    call broadcast(tpar_ky0)
    call broadcast(tprp_ky0)
    call broadcast(qpar_ky0)
    call broadcast(qprp_ky0)
    ! all procs know value of first_half_step already

    call debug_message(verb+1, &
         'gs2_gryfx_zonal::advance_gs2_gryfx broadcasted hybrid arrays')

    call interpolate_theta(gryfx_2_gs2_grid, dens_ky0, .false.)
    call interpolate_theta(gryfx_2_gs2_grid, upar_ky0, .false.)
    call interpolate_theta(gryfx_2_gs2_grid, tpar_ky0, .false.)
    call interpolate_theta(gryfx_2_gs2_grid, tprp_ky0, .false.)
    call interpolate_theta(gryfx_2_gs2_grid, qpar_ky0, .false.)
    call interpolate_theta(gryfx_2_gs2_grid, qprp_ky0, .false.)

    call debug_message(verb+1, &
         'gs2_gryfx_zonal::advance_gs2_gryfx interpolated hybrid arrays')

    ! now that all procs know values, copy into gryfx_zonal structure
    gryfx_zonal%first_half_step = first_half_step 
    gryfx_zonal%NLdens_ky0 = dens_ky0
    gryfx_zonal%NLupar_ky0 = upar_ky0
    gryfx_zonal%NLtpar_ky0 = tpar_ky0
    gryfx_zonal%NLtprp_ky0 = tprp_ky0
    gryfx_zonal%NLqpar_ky0 = qpar_ky0
    gryfx_zonal%NLqprp_ky0 = qprp_ky0

    call debug_message(verb+1, &
         'gs2_gryfx_zonal::advance_gs2_gryfx set hybrid arrays')

    if(first_half_step) then
       state%dont_change_timestep = .true.
       state%skip_diagnostics = .true.
    else 
       state%dont_change_timestep = .false.
       state%skip_diagnostics = .false.
    endif

    call debug_message(verb+1, &
         'gs2_gryfx_zonal::advance_gs2_gryfx evolving equations')
    ! only advance one timestep
    call evolve_equations(state, 1)
    ! getmoms_gryfx takes moments of g and puts results into arguments.
    ! since these are the same arguments that are passed into 
    ! advance_gs2_gryfx, we don't need to copy them like on the way in

    !write (*,*) 'Calling getmoms_gryfx'

    call debug_message(verb+1, &
         'gs2_gryfx_zonal::advance_gs2_gryfx getting gryfx moments')

    call getmoms_gryfx_dist(dens_ky0, upar_ky0, tpar_ky0, &
         tprp_ky0, qpar_ky0, qprp_ky0, phi_ky0)

    call debug_message(verb+1, &
         'gs2_gryfx_zonal::advance_gs2_gryfx interpolating back')

    call interpolate_theta(gs2_2_gryfx_grid, dens_ky0, .false.)
    call interpolate_theta(gs2_2_gryfx_grid, upar_ky0, .false.)
    call interpolate_theta(gs2_2_gryfx_grid, tpar_ky0, .false.)
    call interpolate_theta(gs2_2_gryfx_grid, tprp_ky0, .false.)
    call interpolate_theta(gs2_2_gryfx_grid, qpar_ky0, .false.)
    call interpolate_theta(gs2_2_gryfx_grid, qprp_ky0, .false.)
    call interpolate_theta(gs2_2_gryfx_grid, phi_ky0, .true.)

    call debug_message(verb+1, &
         'gs2_gryfx_zonal::advance_gs2_gryfx finished')

  end subroutine advance_gs2_gryfx