submodule (parallel_map_m) parallel_map_build_s !! Rotines related with build of map implicit none contains module subroutine build_parallel_map(self, comm_handler, equi, mesh_cano, mesh_stag, dbgout) class(parallel_map_t), intent(inout), target :: self type(comm_handler_t), intent(in) :: comm_handler class(equilibrium_t) :: equi type(mesh_cart_t), intent(in) :: mesh_cano type(mesh_cart_t), intent(in) :: mesh_stag integer, intent(in), optional :: dbgout integer :: l, iz, n_points, n_points_stag, outi real(GP) :: phi_cano, phi_stag, dphi, dphi_half, x, y real(GP) :: fluxexp_fwd, fluxexp_bwd, vb type(mesh_cart_t), allocatable :: mesh_target self%nplanes = comm_handler%get_nplanes() self%iplane = comm_handler%get_plane() self%phi_cano = mesh_cano%get_phi() self%phi_stag = mesh_stag%get_phi() self%dphi = 2.0_GP * (self%phi_stag - self%phi_cano) dphi = self%dphi dphi_half = dphi/2.0_GP phi_cano = self%phi_cano phi_stag = self%phi_stag ! set output level outi = 0 if (present(dbgout)) then if (is_rank_info_writer) then outi = dbgout endif if (dbgout >= 3) then outi = dbgout ! every rank writes endif endif ! Build Map matrices to stag from cano call perf_start('../map_build_matrices') allocate(mesh_target) call mesh_cano%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=1, mesh_received=mesh_target) allocate(self%map_stag_cano_fwd) call create_map_matrix(map_matrix = self%map_stag_cano_fwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_stag, & mesh_target = mesh_target, & phi_direction = +1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) allocate(mesh_target) call mesh_cano%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=0, mesh_received=mesh_target) allocate(self%map_stag_cano_bwd) call create_map_matrix(map_matrix = self%map_stag_cano_bwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_stag, & mesh_target = mesh_target, & phi_direction = -1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) ! Build map matrices to cano from stag if (equi%is_axisymmetric()) then self%map_cano_stag_fwd => self%map_stag_cano_fwd self%map_cano_stag_bwd => self%map_stag_cano_bwd else allocate(mesh_target) call mesh_stag%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=0, mesh_received=mesh_target) allocate(self%map_cano_stag_fwd_mem) self%map_cano_stag_fwd => self%map_cano_stag_fwd_mem call create_map_matrix(map_matrix = self%map_cano_stag_fwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_cano, & mesh_target = mesh_target, & phi_direction = +1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) allocate(mesh_target) call mesh_stag%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=-1, mesh_received=mesh_target) allocate(self%map_cano_stag_bwd_mem) self%map_cano_stag_bwd => self%map_cano_stag_bwd_mem call create_map_matrix(map_matrix = self%map_cano_stag_bwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_cano, & mesh_target = mesh_target, & phi_direction = -1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) endif ! Build parallel gradient matrices to stag from cano allocate(mesh_target) call mesh_cano%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=+1, mesh_received=mesh_target) allocate(self%grad_stag_cano_fwd) call create_map_matrix(map_matrix = self%grad_stag_cano_fwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_stag, & mesh_target = mesh_target, & phi_direction = +1, & intorder = self%intorder, & xorder = self%xorder, & flux_box_surface_integral = .true., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) allocate(mesh_target) call mesh_cano%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=0, mesh_received=mesh_target) allocate(self%grad_stag_cano_bwd) call create_map_matrix(map_matrix = self%grad_stag_cano_bwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_stag, & mesh_target = mesh_target, & phi_direction = -1, & intorder = self%intorder, & xorder = self%xorder, & flux_box_surface_integral = .true., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) ! Build map matrices to cano from cano allocate(mesh_target) call mesh_cano%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=+1, mesh_received=mesh_target) allocate(self%map_cano_cano_fwd) call create_map_matrix(map_matrix = self%map_cano_cano_fwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_cano, & mesh_target = mesh_target, & phi_direction = +1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) allocate(mesh_target) call mesh_cano%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=-1, mesh_received=mesh_target) allocate(self%map_cano_cano_bwd) call create_map_matrix(map_matrix = self%map_cano_cano_bwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_cano, & mesh_target = mesh_target, & phi_direction = -1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) ! Build map matrices to cano from next but one cano allocate(mesh_target) call mesh_cano%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=+2, mesh_received=mesh_target) allocate(self%map_cano_cano_fwd2) call create_map_matrix(map_matrix = self%map_cano_cano_fwd2, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_cano, & mesh_target = mesh_target, & phi_direction = +1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) allocate(mesh_target) call mesh_cano%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=-2, mesh_received=mesh_target) allocate(self%map_cano_cano_bwd2) call create_map_matrix(map_matrix = self%map_cano_cano_bwd2, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_cano, & mesh_target = mesh_target, & phi_direction = -1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) ! Build map matrices to stag from stag if (equi%is_axisymmetric()) then self%map_stag_stag_fwd => self%map_cano_cano_fwd self%map_stag_stag_bwd => self%map_cano_cano_bwd else allocate(mesh_target) call mesh_stag%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=+1, mesh_received=mesh_target) allocate(self%map_stag_stag_fwd_mem) self%map_stag_stag_fwd => self%map_stag_stag_fwd_mem call create_map_matrix(map_matrix = self%map_stag_stag_fwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_stag, & mesh_target = mesh_target, & phi_direction = +1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) allocate(mesh_target) call mesh_stag%exchange_mesh_mpi(comm_handler%get_comm_planes(), & step=-1, mesh_received=mesh_target) allocate(self%map_stag_stag_bwd_mem) self%map_stag_stag_bwd => self%map_stag_stag_bwd_mem call create_map_matrix(map_matrix = self%map_stag_stag_bwd, & comm = comm_handler%get_comm_planes(), & equi = equi, & mesh_base = mesh_stag, & mesh_target = mesh_target, & phi_direction = -1, & intorder = self%intorder, & xorder = 0, & flux_box_surface_integral = .false., & use_fixed_stencil = self%use_fixed_stencil, & use_gauss_quadrature = self%use_gauss_quadrature, & dbgout = outi ) deallocate(mesh_target) endif call perf_stop('../map_build_matrices') call perf_start('../map_build_meta') ! Compute parallel distances and flux box volumes for canonical mesh n_points = mesh_cano%get_n_points() allocate(self%dpar_cano_stag_fwd(n_points)) allocate(self%dpar_cano_stag_bwd(n_points)) allocate(self%dpar_cano_cano_fwd(n_points)) allocate(self%dpar_cano_cano_bwd(n_points)) allocate(self%fluxbox_vol_cano(n_points)) !$omp parallel & !$omp default(none) & !$omp shared(n_points, self, phi_cano, dphi_half, dphi, equi, mesh_cano) & !$omp private(l, x, y, fluxexp_fwd, fluxexp_bwd) !$omp do do l = 1, n_points x = mesh_cano%get_x(l) y = mesh_cano%get_y(l) call trace(x_init = x, & y_init = y, & phi_init = phi_cano, & dphi = dphi_half, & equi = equi, & arclen = self%dpar_cano_stag_fwd(l), & fluxexp = fluxexp_fwd) call trace(x_init = x, & y_init = y, & phi_init = phi_cano, & dphi = -dphi_half, & equi = equi, & arclen = self%dpar_cano_stag_bwd(l), & fluxexp = fluxexp_bwd) self%fluxbox_vol_cano(l) = abs( equi%btor(x, y, phi_cano) * & mesh_cano%get_spacing_f()**2 * (fluxexp_fwd + fluxexp_bwd) ) call trace(x_init = x, & y_init = y, & phi_init = phi_cano, & dphi = dphi, & equi = equi, & arclen = self%dpar_cano_cano_fwd(l)) call trace(x_init = x, & y_init = y, & phi_init = phi_cano, & dphi = -dphi, & equi = equi, & arclen = self%dpar_cano_cano_bwd(l)) enddo !$omp end do !$omp end parallel ! Parallel distances and flux box volumes for staggerd mesh n_points_stag = mesh_stag%get_n_points() if (equi%is_axisymmetric()) then ! For axisymmetric case quantities for staggered mesh ! are the same asfor canonical mesh self%dpar_stag_cano_fwd => self%dpar_cano_stag_fwd self%dpar_stag_cano_bwd => self%dpar_cano_stag_bwd self%dpar_stag_stag_fwd => self%dpar_cano_cano_fwd self%dpar_stag_stag_bwd => self%dpar_cano_cano_bwd self%fluxbox_vol_stag => self%fluxbox_vol_cano else ! For non axisymetric case compute quantities for staggered mesh allocate(self%dpar_stag_cano_fwd_mem(n_points_stag)) self%dpar_stag_cano_fwd => self%dpar_stag_cano_fwd_mem allocate(self%dpar_stag_cano_bwd_mem(n_points_stag)) self%dpar_stag_cano_bwd => self%dpar_stag_cano_bwd_mem allocate(self%dpar_stag_stag_fwd_mem(n_points_stag)) self%dpar_stag_stag_fwd => self%dpar_stag_stag_fwd_mem allocate(self%dpar_stag_stag_bwd_mem(n_points_stag)) self%dpar_stag_stag_bwd => self%dpar_stag_stag_bwd_mem allocate(self%fluxbox_vol_stag_mem(n_points_stag)) self%fluxbox_vol_stag => self%fluxbox_vol_stag_mem !$omp parallel & !$omp default(none) & !$omp shared(n_points_stag, self, phi_stag, dphi_half, dphi, equi, mesh_stag) & !$omp private(l, x, y, fluxexp_fwd, fluxexp_bwd) !$omp do do l = 1, n_points_stag x = mesh_stag%get_x(l) y = mesh_stag%get_y(l) call trace(x_init = x, & y_init = y, & phi_init = phi_stag, & dphi = dphi_half, & equi = equi, & arclen = self%dpar_stag_cano_fwd(l), & fluxexp = fluxexp_fwd) call trace(x_init = x, & y_init = y, & phi_init = phi_stag, & dphi = -dphi_half, & equi = equi, & arclen = self%dpar_stag_cano_bwd(l), & fluxexp = fluxexp_bwd) self%fluxbox_vol_stag(l) = abs( equi%btor(x, y, phi_stag) * & mesh_stag%get_spacing_f()**2 * (fluxexp_fwd + fluxexp_bwd) ) call trace(x_init = x, & y_init = y, & phi_init = phi_stag, & dphi = dphi, & equi = equi, & arclen = self%dpar_stag_stag_fwd(l)) call trace(x_init = x, & y_init = y, & phi_init = phi_stag, & dphi = -dphi, & equi = equi, & arclen = self%dpar_stag_stag_bwd(l)) enddo !$omp end do !$omp end parallel endif ! Divide matrices Q by vb !$omp parallel & !$omp default(none) & !$omp shared(n_points_stag, self, phi_stag, equi, mesh_stag) & !$omp private(l, x, y, vb, iz) !$omp do do l = 1, n_points_stag x = mesh_stag%get_x(l) y = mesh_stag%get_y(l) vb = (self%dpar_stag_cano_fwd(l) + self%dpar_stag_cano_bwd(l)) & * equi%btor(x, y, phi_stag) ! Length along field line do iz = self%grad_stag_cano_fwd%i(l), self%grad_stag_cano_fwd%i(l+1)-1 self%grad_stag_cano_fwd%val(iz) = self%grad_stag_cano_fwd%val(iz) / vb enddo do iz = self%grad_stag_cano_bwd%i(l), self%grad_stag_cano_bwd%i(l+1)-1 self%grad_stag_cano_bwd%val(iz) = self%grad_stag_cano_bwd%val(iz) / vb enddo enddo !$omp end do !$omp end parallel ! Build derived operators, i.e. parallel divergence matrices pdiv_cano_stag_fwd and pdiv_cano_stag_bwd call self%build_derived_operators(comm_handler) #if ENABLE_CUDA call self%build_cuda_struct() #endif call perf_stop('../map_build_meta') if (outi >0) then call self%display() endif end subroutine module subroutine build_derived_operators(self, comm_handler) class(parallel_map_t), intent(inout) :: self type(comm_handler_t), intent(in) :: comm_handler integer :: n_points, n_points_stag, n_points_stag_bwd, l, iz, j real(GP), allocatable, dimension(:) :: fluxbox_vol_stag_bwd type(csrmat_t), allocatable :: wrk_csr ! Compute parallel divergence matrices ! pdiv_cano_stag_fwd = -V^(-1)*grad_stag_cano_bwd^T*V_stag ! pdiv_cano_stag_bwd = V^(-1)*grad_stag_cano_fwd^T*V_stag_bwd ! Transposition call getdata_fwdbwdplane(comm_handler%get_comm_planes(), step=-1, & acsr_send=self%grad_stag_cano_fwd, & acsr_recv=wrk_csr) call csr_transpose(wrk_csr, self%pdiv_cano_stag_bwd) deallocate(wrk_csr) call csr_transpose(self%grad_stag_cano_bwd, self%pdiv_cano_stag_fwd) n_points = self%pdiv_cano_stag_fwd%ndim n_points_stag = self%pdiv_cano_stag_fwd%ncol ! Get flux box volumes from backward staggered mesh call getdata_fwdbwdplane(comm_handler%get_comm_planes(), step=-1, & nsend=n_points_stag, & usend=self%fluxbox_vol_stag, & nrecv=n_points_stag_bwd, & urecv=fluxbox_vol_stag_bwd) ! Correct with flux box volume factors !$omp parallel & !$omp default(none) & !$omp shared(n_points, self, fluxbox_vol_stag_bwd) & !$omp private(l, iz, j) !$omp do do l = 1, n_points do iz = self%pdiv_cano_stag_fwd%i(l), self%pdiv_cano_stag_fwd%i(l+1)-1 j = self%pdiv_cano_stag_fwd%j(iz) self%pdiv_cano_stag_fwd%val(iz) = -self%pdiv_cano_stag_fwd%val(iz) * self%fluxbox_vol_stag(j) / self%fluxbox_vol_cano(l) enddo do iz = self%pdiv_cano_stag_bwd%i(l), self%pdiv_cano_stag_bwd%i(l+1)-1 j = self%pdiv_cano_stag_bwd%j(iz) self%pdiv_cano_stag_bwd%val(iz) = -self%pdiv_cano_stag_bwd%val(iz) * fluxbox_vol_stag_bwd(j) / self%fluxbox_vol_cano(l) enddo enddo !$omp end do !$omp end parallel deallocate(fluxbox_vol_stag_bwd) end subroutine module subroutine destructor_parallel_map_t(self) type(parallel_map_t), intent(inout) :: self if (associated(self%map_stag_stag_bwd)) nullify(self%map_stag_stag_bwd) if (allocated(self%map_stag_stag_bwd_mem)) deallocate(self%map_stag_stag_bwd_mem) if (associated(self%map_stag_stag_fwd)) nullify(self%map_stag_stag_fwd) if (allocated(self%map_stag_stag_fwd_mem)) deallocate(self%map_stag_stag_fwd_mem) if (allocated(self%map_cano_cano_bwd)) deallocate(self%map_cano_cano_bwd) if (allocated(self%map_cano_cano_fwd)) deallocate(self%map_cano_cano_fwd) if (allocated(self%pdiv_cano_stag_bwd)) deallocate(self%pdiv_cano_stag_bwd) if (allocated(self%pdiv_cano_stag_fwd)) deallocate(self%pdiv_cano_stag_fwd) if (allocated(self%grad_stag_cano_bwd)) deallocate(self%grad_stag_cano_bwd) if (allocated(self%grad_stag_cano_fwd)) deallocate(self%grad_stag_cano_fwd) if (associated(self%map_cano_stag_bwd)) nullify(self%map_cano_stag_bwd) if (allocated(self%map_cano_stag_bwd_mem)) deallocate(self%map_cano_stag_bwd_mem) if (associated(self%map_cano_stag_fwd)) nullify(self%map_cano_stag_fwd) if (allocated(self%map_cano_stag_fwd_mem)) deallocate(self%map_cano_stag_fwd_mem) if (allocated(self%map_stag_cano_bwd)) deallocate(self%map_stag_cano_bwd) if (allocated(self%map_stag_cano_fwd)) deallocate(self%map_stag_cano_fwd) if (associated(self%fluxbox_vol_stag)) nullify(self%fluxbox_vol_stag) if (allocated(self%fluxbox_vol_stag_mem)) deallocate(self%fluxbox_vol_stag_mem) if (associated(self%dpar_stag_stag_bwd)) nullify(self%dpar_stag_stag_bwd) if (allocated(self%dpar_stag_stag_bwd_mem)) deallocate(self%dpar_stag_stag_bwd_mem) if (associated(self%dpar_stag_stag_fwd)) nullify(self%dpar_stag_stag_fwd) if (allocated(self%dpar_stag_stag_fwd_mem)) deallocate(self%dpar_stag_stag_fwd_mem) if (associated(self%dpar_stag_cano_bwd)) nullify(self%dpar_stag_cano_bwd) if (allocated(self%dpar_stag_cano_bwd_mem)) deallocate(self%dpar_stag_cano_bwd_mem) if (associated(self%dpar_stag_cano_fwd)) nullify(self%dpar_stag_cano_fwd) if (allocated(self%dpar_stag_cano_fwd_mem)) deallocate(self%dpar_stag_cano_fwd_mem) if (allocated(self%fluxbox_vol_cano)) deallocate(self%fluxbox_vol_cano) if (allocated(self%dpar_cano_cano_bwd)) deallocate(self%dpar_cano_cano_bwd) if (allocated(self%dpar_cano_cano_fwd)) deallocate(self%dpar_cano_cano_fwd) if (allocated(self%dpar_cano_stag_bwd)) deallocate(self%dpar_cano_stag_bwd) if (allocated(self%dpar_cano_stag_fwd)) deallocate(self%dpar_cano_stag_fwd) end subroutine destructor_parallel_map_t #if ENABLE_CUDA module subroutine build_cuda_struct(self) class(parallel_map_t), intent(inout), target :: self call allocate_cuda_struct(self%grad_stag_cano_fwd_cuda, self%grad_stag_cano_fwd%ndim, & self%grad_stag_cano_fwd%nnz, self%grad_stag_cano_fwd%i, self%grad_stag_cano_fwd%j, & self%grad_stag_cano_fwd%val) call allocate_cuda_struct(self%grad_stag_cano_bwd_cuda, self%grad_stag_cano_bwd%ndim, & self%grad_stag_cano_bwd%nnz, self%grad_stag_cano_bwd%i, self%grad_stag_cano_bwd%j, & self%grad_stag_cano_bwd%val) call allocate_cuda_struct(self%pdiv_cano_stag_bwd_cuda, self%pdiv_cano_stag_bwd%ndim, & self%pdiv_cano_stag_bwd%nnz, self%pdiv_cano_stag_bwd%i, self%pdiv_cano_stag_bwd%j, & self%pdiv_cano_stag_bwd%val) call allocate_cuda_struct(self%pdiv_cano_stag_fwd_cuda, self%pdiv_cano_stag_fwd%ndim, & self%pdiv_cano_stag_fwd%nnz, self%pdiv_cano_stag_fwd%i, self%pdiv_cano_stag_fwd%j, & self%pdiv_cano_stag_fwd%val) call allocate_cuda_struct(self%map_cano_stag_bwd_cuda, self%map_cano_stag_bwd%ndim, & self%map_cano_stag_bwd%nnz, self%map_cano_stag_bwd%i, self%map_cano_stag_bwd%j, & self%map_cano_stag_bwd%val) call allocate_cuda_struct(self%map_cano_stag_fwd_cuda, self%map_cano_stag_fwd%ndim, & self%map_cano_stag_fwd%nnz, self%map_cano_stag_fwd%i, self%map_cano_stag_fwd%j, & self%map_cano_stag_fwd%val) end subroutine subroutine allocate_cuda_struct(cuda_map, nrows, nnz, rows, cols, vals) !! Exposes CSR matrix to CUDA type(cuda_map_t), intent(inout) :: cuda_map !! Contains exposure information of CSR matrix integer, intent(in) :: nrows !! Number of rows of matrix integer, intent(in) :: nnz !! Number of non-zeros of matrix integer, dimension(nrows+1), intent(in) :: rows !! Row array of CSR matrix integer, dimension(nnz), intent(in) :: cols !! Column indices of CSR matrix real(GP), dimension(nnz), intent(in) :: vals !! Values of CSR matrix cuda_map%nrows = nrows cuda_map%nnz = nnz call allocate_memory_int(cuda_map%cols, 1, cuda_map%nnz) call allocate_memory_int(cuda_map%rows, 1, cuda_map%nrows+1) call allocate_memory_double(cuda_map%vals, 1, cuda_map%nnz) call allocate_memory_double(cuda_map%x, 1, cuda_map%nrows) call allocate_memory_double(cuda_map%y, 1, cuda_map%nrows) cuda_map%rows = rows - 1 cuda_map%cols = cols - 1 cuda_map%vals = vals call allocate_cusparse_struct(cuda_map%cuda_struct, cuda_map%nrows, & cuda_map%nrows, cuda_map%nnz, cuda_map%rows, & cuda_map%cols, cuda_map%vals, cuda_map%x, cuda_map%y) end subroutine module subroutine destructor_cuda_map_t(self) type(cuda_map_t), intent(inout) :: self call freeint(self%rows) call freeint(self%cols) call freedouble(self%vals) call freedouble(self%x) call freedouble(self%y) self%cuda_struct = c_null_ptr end subroutine #endif end submodule