!> FIXME : Add documentation module le_derivatives implicit none private public :: vspace_derivatives, time_vspace_derivatives, time_vspace_derivatives_mpi real :: time_vspace_derivatives(2) = 0., time_vspace_derivatives_mpi = 0. contains !> FIXME : Add documentation subroutine vspace_derivatives (g, gold, g1, phi, bpar, phinew, bparnew, diagnostics, gtoc, ctog) use array_utils, only: copy, zero_array use redistribute, only: gather, scatter use dist_fn_arrays, only: g_adjust, to_g_gs2, from_g_gs2 use gs2_layouts, only: g_lo, le_lo use gs2_time, only: code_dt use theta_grid, only: ntgrid use le_grids, only: integrate_moment, negrid, g2le use collisions, only: solfp1, colls, adjust, heating, hyper_colls, use_le_layout, nxi_lim, c_rate use job_manage, only: time_message use mp, only: get_mp_times use optionals, only: get_option_with_default implicit none complex, dimension (-ntgrid:,:,g_lo%llim_proc:), intent (in out) :: g, gold, g1 complex, dimension (-ntgrid:,:,:), intent (in) :: phi, bpar, phinew, bparnew integer, optional, intent (in) :: diagnostics !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 !CMR, 3/10/2013: !CMR Just realised we need to be careful with g_adjust if we avoid !CMR mapping g_lo <=> le_lo Mmmm, a little more to think about here ;-( !CMR logical, intent(in), optional :: gtoc, ctog logical :: g_to_c, c_to_g complex, dimension (:,:,:), allocatable :: gle complex, dimension (:,:,:), allocatable :: gc1, gc2, gc3 logical :: heating_flag real :: mp_total, mp_total_after call time_message(.false., time_vspace_derivatives, "Vspace derivative") call get_mp_times(total_time = mp_total) g_to_c = get_option_with_default(gtoc, .true.) c_to_g = get_option_with_default(ctog, .true.) heating_flag = heating .and. present(diagnostics) if(use_le_layout) then if (colls) then if (adjust) then call g_adjust (g, phinew, bparnew, direction = from_g_gs2) if (heating_flag) call g_adjust (gold, phi, bpar, direction = from_g_gs2) end if if (heating_flag) then allocate (gc3(-ntgrid:ntgrid,2,g_lo%llim_proc:g_lo%ulim_alloc)) call copy(g, gc3) end if allocate (gle(nxi_lim,negrid+1,le_lo%llim_proc:le_lo%ulim_alloc)) call zero_array(gle) ! map data from g_layout to le_layout if (g_to_c) call gather (g2le, g, gle) if (colls) then ! update distribution function to take into account collisions call solfp1 (gle, diagnostics) end if ! remap from le_layout to g_layout if (c_to_g) call scatter (g2le, gle, g) deallocate (gle) if (heating_flag) then ! form (h_i+1 + h_i)/2 * C(h_i+1) and integrate. gc3 = 0.5*conjg(g+gold)*(g-gc3)/code_dt ! We only want the real part of this result so we could ! either provide a real argument overload or reduce to ! a complex temporary before storing just the real part. ! This would allow c_rate to be declared real, halving the ! memory requirements. call integrate_moment (gc3, c_rate(:,:,:,:,3)) deallocate (gc3) endif if (adjust) then call g_adjust (g, phinew, bparnew, direction = to_g_gs2) if (heating_flag) call g_adjust (gold, phi, bpar, direction = to_g_gs2) end if end if else if (colls) then if (heating_flag) then allocate (gc1(-ntgrid:ntgrid,2,g_lo%llim_proc:g_lo%ulim_alloc)) if (hyper_colls) then allocate (gc2(-ntgrid:ntgrid,2,g_lo%llim_proc:g_lo%ulim_alloc)) else allocate (gc2(1,1,1)) end if allocate (gc3(-ntgrid:ntgrid,2,g_lo%llim_proc:g_lo%ulim_alloc)) else allocate (gc1(1,1,1)) allocate (gc2(1,1,1)) allocate (gc3(1,1,1)) end if call zero_array(gc1) call zero_array(gc2) call zero_array(gc3) if (adjust) then call g_adjust (g, phinew, bparnew, direction = from_g_gs2) if (heating_flag) call g_adjust (gold, phi, bpar, direction = from_g_gs2) end if if (heating_flag) call copy(g, gc3) ! update distribution function to take into account collisions call solfp1 (g, g1, gc1, gc2, diagnostics) if (heating_flag) then call integrate_moment (gc1, c_rate(:,:,:,:,1)) deallocate (gc1) ! We only want the real part of this result so we could ! either provide a real argument overload or reduce to ! a complex temporary before storing just the real part. ! This would allow c_rate to be declared real, halving the ! memory requirements. if (hyper_colls) call integrate_moment (gc2, c_rate(:,:,:,:,2)) deallocate (gc2) ! form (h_i+1 + h_i)/2 * C(h_i+1) and integrate. gc3 = 0.5*conjg(g+gold)*(g-gc3)/code_dt call integrate_moment (gc3, c_rate(:,:,:,:,3)) deallocate (gc3) end if if (adjust) then call g_adjust (g, phinew, bparnew, direction = to_g_gs2) if (heating_flag) call g_adjust (gold, phi, bpar, direction = to_g_gs2) end if end if end if call time_message(.false., time_vspace_derivatives, "Vspace derivative") call get_mp_times(total_time = mp_total_after) time_vspace_derivatives_mpi = time_vspace_derivatives_mpi + (mp_total_after - mp_total) end subroutine vspace_derivatives end module le_derivatives