submodule(splitting_m) splitting_gauss_seidel_redblack_cpu_s !! Splitting solver with Red-Black Gauss-Seidel splitting method on CPU use screen_io_m, only : get_stderr implicit none contains module subroutine create_gauss_seidel_redblack_cpu(self, mesh) class(splitting_gauss_seidel_redblack_cpu_t), intent(inout) :: self type(mesh_cart_t), intent(in), target :: mesh self%ndim = mesh%get_n_points() self%n_points_red = mesh%get_n_points_red() self%n_points_black = mesh%get_n_points_black() self%redblack_indices => mesh%redblack_indices end subroutine module subroutine apply_gauss_seidel_redblack_cpu(self, a, dinv, x, b) class(splitting_gauss_seidel_redblack_cpu_t), intent(inout) :: self type(csrmat_t), intent(in) :: a real(FP), intent(in), dimension(self%ndim) :: dinv real(FP), intent(inout), dimension(self%ndim) :: x real(FP), intent(in), dimension(self%ndim) :: b if (.not.((a%ndim) == (self%ndim))) then write(get_stderr(), *) & 'error splitting_gauss_seidel_redblack_cpu%apply): & &called for wrong-size matrix' write(get_stderr(), *) & 'matrix ndim is', a%ndim, ', gsrb%ndim is', self%ndim error stop endif call compute_subloop(1, self%n_points_red) call compute_subloop(self%n_points_red + 1, & self%n_points_red + self%n_points_black) call compute_subloop(self%n_points_red + self%n_points_black + 1, & self%ndim) contains subroutine compute_subloop(kinit, kfinal) !! Contained subroutine for the copy/pasted inner loop lines integer, intent(in) :: kinit integer, intent(in) :: kfinal integer :: kk, ll ! Sum is supposed to avoid the diagonal term (which would be ! x(l) / d_inverse(l)). Adding x(l) to the entire right-hand ! side effectively removes the diagonal term. Also this avoids ! the divisions/multiplications required if the diagonal term ! would be dealt with directly by the dot_product. !$omp do private(kk, ll) do kk = kinit, kfinal ll = self%redblack_indices(kk) x(ll) = x(ll) + (b(ll) - & dot_product(a%val(a%i(ll):a%i(ll+1)-1), & x(a%j(a%i(ll):a%i(ll+1)-1)))) * dinv(ll) enddo !$omp end do end subroutine compute_subloop end subroutine module subroutine destructor_gauss_seidel_redblack_cpu(self) type(splitting_gauss_seidel_redblack_cpu_t), intent(inout) :: self self%redblack_indices => null() self%ndim = 0 self%n_points_red = 0 self%n_points_black = 0 end subroutine end submodule