csr_build_from_coo Subroutine

   subroutine csr_build_from_coo(self, n_rows, n_cols, rows, cols, vals, err)
      !! Build the matrix from coordinate (COO) triplets.
      !!
      !! `rows`, `cols` and `vals` must have the same length and may be given
      !! in any order; entries are bucketed by row and then sorted by column.
      !! **Repeated `(i, j)` coordinates are summed.** Explicitly stored zeros
      !! are kept, so a triplet with value zero still occupies a slot.
      !!
      !! Reports `ERROR_VALIDATION` when the three arrays differ in length,
      !! when `n_rows` or `n_cols` is negative, or when any coordinate falls
      !! outside `1:n_rows` / `1:n_cols`. On any error the matrix is left
      !! empty rather than half built.
      class(csr_t), intent(inout) :: self
      integer(default_int), intent(in) :: n_rows
         !! Number of rows of the assembled matrix
      integer(default_int), intent(in) :: n_cols
         !! Number of columns of the assembled matrix
      integer(default_int), intent(in) :: rows(:)
         !! Row coordinate of each triplet
      integer(default_int), intent(in) :: cols(:)
         !! Column coordinate of each triplet
      real(dp), intent(in) :: vals(:)
         !! Value of each triplet
      type(error_t), intent(out), optional :: err

      integer(default_int) :: n_in, k, i, j, lo, hi, write_pos
      integer(default_int), allocatable :: start(:), next_free(:), tmp_col(:)
      real(dp), allocatable :: tmp_val(:)

      call self%destroy()

      n_in = int(size(rows), default_int)
      if (int(size(cols), default_int) /= n_in .or. int(size(vals), default_int) /= n_in) then
         call fail(err, "pic_csr: build_from_coo requires rows, cols and vals of equal length")
         return
      end if

      if (n_rows < 0 .or. n_cols < 0) then
         call fail(err, "pic_csr: build_from_coo requires non-negative matrix dimensions")
         return
      end if

      do k = 1, n_in
         if (rows(k) < 1 .or. rows(k) > n_rows) then
            call fail(err, "pic_csr: build_from_coo got a row coordinate outside 1:n_rows")
            return
         end if
         if (cols(k) < 1 .or. cols(k) > n_cols) then
            call fail(err, "pic_csr: build_from_coo got a column coordinate outside 1:n_cols")
            return
         end if
      end do

      self%nrow = n_rows
      self%ncol = n_cols

      allocate (self%row_ptr(n_rows + 1))
      allocate (start(n_rows + 1))
      allocate (next_free(n_rows))
      allocate (tmp_col(n_in))
      allocate (tmp_val(n_in))

      ! bucket by row: count, prefix sum, then scatter
      start = 0
      do k = 1, n_in
         start(rows(k) + 1) = start(rows(k) + 1) + 1
      end do
      start(1) = 1
      do i = 1, n_rows
         start(i + 1) = start(i + 1) + start(i)
         next_free(i) = start(i)
      end do
      do k = 1, n_in
         i = rows(k)
         tmp_col(next_free(i)) = cols(k)
         tmp_val(next_free(i)) = vals(k)
         next_free(i) = next_free(i) + 1
      end do

      ! sort each row by column index, then merge duplicates by summing
      write_pos = 0
      self%row_ptr(1) = 1
      do i = 1, n_rows
         lo = start(i)
         hi = start(i + 1) - 1
         call sort_row(tmp_col, tmp_val, lo, hi)
         k = lo
         do while (k <= hi)
            write_pos = write_pos + 1
            tmp_col(write_pos) = tmp_col(k)
            tmp_val(write_pos) = tmp_val(k)
            j = k + 1
            do while (j <= hi)
               if (tmp_col(j) /= tmp_col(write_pos)) exit
               tmp_val(write_pos) = tmp_val(write_pos) + tmp_val(j)
               j = j + 1
            end do
            k = j
         end do
         self%row_ptr(i + 1) = write_pos + 1
      end do

      allocate (self%col_idx(write_pos))
      allocate (self%values(write_pos))
      do k = 1, write_pos
         self%col_idx(k) = tmp_col(k)
         self%values(k) = tmp_val(k)
      end do
   end subroutine csr_build_from_coo