enforce_parity Subroutine

private subroutine enforce_parity(redist_obj, ik_ind)

Enforce requested parity

Arguments

Type IntentOptional Attributes Name
type(redist_type), intent(in), optional :: redist_obj
integer, intent(in), optional :: ik_ind

Contents

Source Code


Source Code

  subroutine enforce_parity(redist_obj,ik_ind)
    use theta_grid, only:ntgrid
    use dist_fn_arrays, only: gnew
    use redistribute, only: scatter,redist_type
    use gs2_layouts, only: g_lo,ik_idx
    implicit none
    type(redist_type),intent(in),optional :: redist_obj
    integer, intent(in),optional :: ik_ind
    type(redist_type) :: redist_local
    integer :: ik_local
    integer :: iglo,mult

    !If enforced parity not requested then exit
    if(.not.(def_parity))return
    
    !Set multiplier
    if(even) then
       mult=1
    else
       mult=-1
    endif
    
    !Behaviour depends upon if we're in flux tube or ballooning space
    if(boundary_option_switch.eq.boundary_option_linked) then !Flux-tube
       !Ensure a redist object is present, if not default to parity_redist
       if(present(redist_obj))then
          redist_local=redist_obj
       else
          redist_local=parity_redist
       endif

       !Redistribute the data
       call scatter(redist_local,gnew,gnew)

       !Multiply by factor if required
       if(mult.ne.1) gnew(:,1,:)=mult*gnew(:,1,:)
    else !Ballooning/extended space
       !Ensure ik_local is specified
       if(present(ik_ind))then
          ik_local=ik_ind
       else
          ik_local=-1
       endif

       !Loop over all local iglo
       do iglo=g_lo%llim_proc,g_lo%ulim_alloc
          !Skip if needed
          if(ik_local.gt.0) then
             if(ik_idx(g_lo,iglo).ne.ik_local) cycle
          endif

          !Apply parity filter
          gnew(-ntgrid:-1,1,iglo)=mult*gnew(ntgrid:1:-1,2,iglo)
          gnew(1:ntgrid,1,iglo)=mult*gnew(-1:-ntgrid:-1,2,iglo)
       enddo
    endif
    
  end subroutine enforce_parity