solfp1_standard_layout Subroutine

private subroutine solfp1_standard_layout(g, g1, gc1, gc2, diagnostics, gtoc, ctog)

FIXME : Add documentation

Arguments

Type IntentOptional Attributes Name
complex, intent(inout), dimension (-ntgrid:,:,g_lo%llim_proc:) :: g
complex, intent(inout), dimension (-ntgrid:,:,g_lo%llim_proc:) :: g1
complex, intent(inout), dimension (-ntgrid:,:,g_lo%llim_proc:) :: gc1
complex, intent(inout), dimension (-ntgrid:,:,g_lo%llim_proc:) :: gc2
integer, intent(in), optional :: diagnostics
logical, intent(in), optional :: gtoc
logical, intent(in), optional :: ctog

Contents


Source Code

  subroutine solfp1_standard_layout (g, g1, gc1, gc2, diagnostics, gtoc, ctog)
    use gs2_layouts, only: g_lo, it_idx, ik_idx, ie_idx, is_idx
    use theta_grid, only: ntgrid
    use run_parameters, only: beta
    use gs2_time, only: code_dt
    use le_grids, only: energy => energy_maxw
    use species, only: spec, is_electron_species
    use dist_fn_arrays, only: vpa, aj0
    use fields_arrays, only: aparnew
    use run_parameters, only: ieqzip
    use kt_grids, only: kwork_filter, kperp2
    use optionals, only: get_option_with_default
    implicit none

    complex, dimension (-ntgrid:,:,g_lo%llim_proc:), intent (in out) :: g, g1, gc1, gc2
    integer, optional, intent (in) :: diagnostics
    logical, optional, intent (in) :: gtoc, ctog
    integer :: isgn, it, ik, ie, is, iglo
!CMR, 12/9/2013: 
!CMR   New logical optional input parameters gtoc, ctog used to set
!CMR   flags (g_to_c and c_to_g) to control whether redistributes required
!CMR   to map g_lo to collision_lo, and collision_lo to g_lo.
!CMR   All redistributes are performed by default.
!CMR  
    logical :: g_to_c, c_to_g

    g_to_c = get_option_with_default(gtoc, .true.)
    c_to_g = get_option_with_default(ctog, .true.)

    if (has_diffuse) then
       call solfp_ediffuse_standard_layout (g)
       if (conserve_moments) call conserve_diffuse (g, g1)
    end if

    if (has_lorentz) then
       if (drag) then
          !$OMP PARALLEL DO DEFAULT(none) &
          !$OMP PRIVATE(iglo, is, it, ik, ie, isgn) &
          !$OMP SHARED(g_lo, spec, kwork_filter, g, ieqzip, vnmult, code_dt, vpa, &
          !$OMP kperp2, aparnew, aj0, beta, energy) &
          !$OMP SCHEDULE(static)
          do iglo = g_lo%llim_proc, g_lo%ulim_proc
             is = is_idx(g_lo,iglo)
             if (.not. is_electron_species(spec(is))) cycle
             it = it_idx(g_lo,iglo)
             ik = ik_idx(g_lo,iglo)
             if(kwork_filter(it,ik)) cycle
             if(ieqzip(it,ik)) cycle
             ie = ie_idx(g_lo,iglo)
             do isgn = 1, 2
                g(:, isgn, iglo) = g(:, isgn, iglo) + vnmult(1)*spec(is)%vnewk*code_dt &
                     * vpa(:,isgn,iglo)*kperp2(:,it,ik)*aparnew(:,it,ik)*aj0(:,iglo) &
                     / ((-spec(is)%z*spec(is)%dens)*beta*spec(is)%stm*energy(ie)**1.5)
                ! probably need 1/(spec(is_ion)%z*spec(is_ion)%dens) above
                ! This has been implemented as 1/-spec(electron)%z*spec(electron)%dens
                ! in an attempt handle the multi-ion species case.
             end do
          end do
          !$OMP END PARALLEL DO
       end if

       call solfp_lorentz (g, gc1, gc2, diagnostics)
       if (conserve_moments) call conserve_lorentz (g, g1)
    end if
  end subroutine solfp1_standard_layout