advance_local Subroutine

  subroutine advance_local(istep, remove_zonal_flows_switch)
    use run_parameters, only: has_phi, has_apar, has_bpar, reset
    use fields_implicit, only: remove_zonal_flows
    use fields_arrays, only: phi, apar, bpar, phinew, time_field
    use fields_arrays, only: aparnew, bparnew, apar_ext, time_field_mpi
    use dist_fn, only: timeadv, exb_shear, g_exb, g_exbfac
    use dist_fn, only: collisions_advance
    use dist_fn_arrays, only: g, gnew
    use antenna, only: antenna_amplitudes, no_driver
    use gs2_layouts, only: g_lo
    use unit_tests, only: debug_message
    use mp, only: get_mp_times
    use job_manage, only: time_message
    use array_utils, only: copy
    implicit none
    integer, intent(in) :: istep
    logical, intent(in) :: remove_zonal_flows_switch
    integer, parameter :: diagnostics = 1
    logical, parameter :: do_gather=.true.
    logical :: do_update
    integer, parameter :: verb = 4
    real :: mp_total, mp_total_after
    !do_gather=.true. => fields are collected from all procs and distributed to everyone (required)
    !do_update=.true. => "Smart" update routines are used which only update the local parts of the
    !field arrays. This could help when xy are strongly parallelised but is likely to be slower
    !when all local, hence we should turn it off if xy all local.
    do_update=do_smart_update

    if(g_lo%x_local.and.g_lo%y_local) do_update=.false.

    !GGH NOTE: apar_ext is initialized in this call
    if(.not.no_driver) call antenna_amplitudes (apar_ext)

    call debug_message(verb, &
      'fields_local::advance_local calling g_exb')
    !Apply flow shear if active
    if(abs(g_exb*g_exbfac)>epsilon(0.)) call exb_shear(gnew,phinew,aparnew,bparnew,istep,field_local_allreduce_sub)

    !Update g and fields
    call copy(gnew,g)
    if(do_update)then !Here we tie the "smart" update to do_update
       call fieldmat%update_fields_newstep
    else
       if (has_phi) call copy(phinew, phi)
       if (has_apar) call copy(aparnew, apar)
       if (has_bpar) call copy(bparnew, bpar)
    endif

    call debug_message(verb, &
      'fields_local::advance_local calling timeadv 1')

    !Find gnew given fields at time step midpoint
    call timeadv (phi, apar, bpar, phinew, &
         aparnew, bparnew, istep)
    if(reset) return !Return is resetting

    !Add in antenna driving if present
    !<DD>TAGGED: Should we only this is fapar>0 as well?
    if(.not.no_driver) aparnew=aparnew+apar_ext

    call debug_message(verb, &
         'fields_local::advance_local calling getfield_local')

    !For timing
    !if(debug)call barrier !!/These barriers influence the reported time
    call time_message(.false.,time_field,' Field Solver')
    call get_mp_times(total_time = mp_total)

    !Calculate the fields at the next time point
    call getfield_local(phinew,aparnew,bparnew,do_gather,do_update)

    !If we do the update in getfield_local don't do it here
    if(.not.do_update)then
       if (has_phi) phinew = phinew + phi
       if (has_apar) aparnew = aparnew + apar
       if (has_bpar) bparnew = bparnew + bpar
    endif

    !For timing
    !Barrier usage ensures that proc0 measures the longest time taken on
    !any proc.
    !if(debug)call barrier !!/These barriers influence the reported time
    call time_message(.false.,time_field,' Field Solver')
    call get_mp_times(total_time = mp_total_after)
    time_field_mpi = time_field_mpi + (mp_total_after - mp_total)

    !Remove zonal component if requested
    if(remove_zonal_flows_switch) call remove_zonal_flows

    call debug_message(verb, &
      'fields_local::advance_local calling timeadv 2')
    !Evolve gnew to next half time point
    call timeadv (phi, apar, bpar, phinew, &
         aparnew, bparnew, istep, diagnostics) 

    ! Advance collisions, if separate from timeadv
    call collisions_advance (phi, bpar, phinew, aparnew, bparnew, diagnostics)

  end subroutine advance_local