fdurastante/psblas3
1
0
Fork 0
mirror of https://github.com/sfilippone/psblas3.git synced 2026-10-20 21:02:26 +00:00
Code Issues Packages Projects Releases Wiki Activity

[ADD] Binary matrix format with a parallel reader

Browse Source 
Matrix Market is ASCII with variable-length lines, so the byte offset of
nonzero k cannot be computed without scanning the file: reading it in parallel
is impossible without a redistribution pass. psb_mm2bin converts it once into a
row-sorted binary layout, and psb_binmat_mod reads that layout in parallel.

Layout, 8-byte native-endian throughout: a header carrying a magic number, the
index and real kinds the file was written with, and the dimensions; then
row_ptr(nrows+1), col_idx(nnz) and values(nnz). The kinds are recorded so a
mismatched reader fails loudly instead of reinterpreting the bytes.

Each rank reads the header, then only its own slice of row_ptr, and from that
seeks straight to the nonzeros of the rows it owns. Because the block-row
partition is implicit in the file order, psb_matdist is not needed at all --
that scatter alone accounted for 28 s of MPI time at 448 ranks.

Geo_1438, 1.44M rows, 63.2M nonzeros:
  Matrix Market read, serial        27.0 s
  conversion, once                  27.9 s
  binary read at 8 ranks, including
  assembly and one SpMV              0.98 s

Verified: byte-for-byte layout check on a 4x4 symmetric matrix with the lower
triangle expanded; A*1 identical at 1, 2, 3 and 4 ranks there, the 3-rank case
covering the uneven split; sum(A*1) identical to sixteen digits at 1, 4 and 8
ranks on Geo_1438.

Positioned stream reads rather than MPI-IO: the ranks touch disjoint byte
ranges, so there is nothing to coordinate. If collective buffering turns out to
matter at high rank counts, MPI_File_read_at_all is a change local to
psb_binmat_mod.
This commit is contained in:
Stack-1
2026-08-17 16:12:19 +02:00
parent b85532c0ce
commit ea273e8e93
2 changed files with 231 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files
test/comm/psb_binmat_mod.f90
+125
View File
@@ -0,0 +1,125 @@
!
! Parallel read of the binary layout written by psb_mm2bin.
!
! Every rank opens the file, reads the header and only its own slice of row_ptr,
! then seeks straight to the nonzeros of the rows it owns. Nothing is read twice
! and nothing is redistributed: the block-row partition is implicit in the file
! order, so psb_matdist -- 28 s of MPI time at 448 ranks -- is not needed at all.
!
! Positioned stream reads are used rather than MPI-IO. The ranks touch disjoint
! byte ranges, so there is nothing to coordinate; if collective buffering ever
! turns out to matter at high rank counts, swapping these for
! MPI_File_read_at_all is a change local to this file.
!
module psb_binmat_mod
use psb_base_mod
implicit none
integer(psb_lpk_), parameter :: psb_bin_magic = 6510534141251461121_psb_lpk_
integer(psb_lpk_), parameter :: psb_bin_version = 1_psb_lpk_
contains
subroutine psb_bin_mat_read(ctxt, fname, afmt, a, desc_a, info)
type(psb_ctxt_type), intent(in) :: ctxt
character(len=*), intent(in) :: fname, afmt
type(psb_dspmat_type), intent(out) :: a
type(psb_desc_type), intent(out) :: desc_a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_) :: iam, np, iu, nloc, ist
integer(psb_lpk_) :: magic, vers, ipkb, rpkb, nr, nc, nz, res
integer(psb_lpk_) :: r1, r2, k1, k2, lnz, i, k, base, cbase, vbase, nb, rem
integer(psb_lpk_), allocatable :: rp(:), ia(:), ja(:)
real(psb_dpk_), allocatable :: val(:)
info = psb_success_
call psb_info(ctxt, iam, np)
iu = 71 + iam
open(unit=iu, file=trim(fname), access='stream', form='unformatted', &
& status='old', action='read', iostat=ist)
if (ist /= 0) then
if (iam == psb_root_) write(psb_err_unit,*) 'cannot open ', trim(fname)
info = psb_err_internal_error_
return
end if
read(iu, pos=1, iostat=ist) magic, vers, ipkb, rpkb, nr, nc, nz, res
if ((ist /= 0).or.(magic /= psb_bin_magic)) then
if (iam == psb_root_) write(psb_err_unit,*) 'not a psb binary matrix: ', trim(fname)
info = psb_err_internal_error_ ; close(iu) ; return
end if
! The header records the kinds the file was written with: a mismatch has to
! fail here rather than silently reinterpret the bytes.
if ((ipkb /= 8).or.(rpkb /= 8)) then
if (iam == psb_root_) write(psb_err_unit,*) 'kind mismatch in ', trim(fname), ipkb, rpkb
info = psb_err_internal_error_ ; close(iu) ; return
end if
if (nr /= nc) then
if (iam == psb_root_) write(psb_err_unit,*) 'matrix is not square'
info = psb_err_internal_error_ ; close(iu) ; return
end if
! Block-row split, in rank order: the same assignment psb_cdall(nl=) makes,
! which is what lets the file order stand in for the partition.
nb = nr/np
rem = mod(nr,np)
if (iam < rem) then
nloc = int(nb+1, psb_ipk_)
r1 = iam*(nb+1) + 1
else
nloc = int(nb, psb_ipk_)
r1 = rem*(nb+1) + (iam-rem)*nb + 1
end if
r2 = r1 + nloc - 1
base = 65 ! row_ptr starts here (1-based stream pos)
cbase = base + (nr+1)*8 ! col_idx
vbase = cbase + nz*8 ! values
allocate(rp(nloc+1), stat=ist)
if (ist /= 0) then
info = psb_err_alloc_dealloc_ ; close(iu) ; return
end if
read(iu, pos=base + (r1-1)*8, iostat=ist) rp
if (ist /= 0) then
info = psb_err_internal_error_ ; close(iu) ; return
end if
k1 = rp(1)
k2 = rp(nloc+1) - 1
lnz = k2 - k1 + 1
if (lnz < 0) lnz = 0
allocate(ia(max(lnz,1_psb_lpk_)), ja(max(lnz,1_psb_lpk_)), &
& val(max(lnz,1_psb_lpk_)), stat=ist)
if (ist /= 0) then
info = psb_err_alloc_dealloc_ ; close(iu) ; return
end if
if (lnz > 0) then
read(iu, pos=cbase + (k1-1)*8, iostat=ist) ja(1:lnz)
if (ist == 0) read(iu, pos=vbase + (k1-1)*8, iostat=ist) val(1:lnz)
if (ist /= 0) then
info = psb_err_internal_error_ ; close(iu) ; return
end if
! Global row index for every nonzero, from the row_ptr slice.
do i = 1, nloc
do k = rp(i), rp(i+1)-1
ia(k-k1+1) = r1 + i - 1
end do
end do
end if
close(iu)
call psb_cdall(ctxt, desc_a, info, nl=nloc)
if (info == psb_success_) call psb_spall(a, desc_a, info, nnz=int(lnz,psb_ipk_))
if ((info == psb_success_).and.(lnz > 0)) &
& call psb_spins(int(lnz,psb_ipk_), ia, ja, val, a, desc_a, info)
if (info == psb_success_) call psb_cdasb(desc_a, info)
if (info == psb_success_) call psb_spasb(a, desc_a, info, afmt=afmt)
end subroutine psb_bin_mat_read
end module psb_binmat_mod
test/comm/psb_mm2bin.f90
+106
View File
@@ -0,0 +1,106 @@
!
! Convert a Matrix Market file into a row-sorted binary CSR-like layout that can
! be read in parallel.
!
! Why this exists: Matrix Market is ASCII with variable-length lines, so the byte
! offset of nonzero k cannot be computed without scanning the file. Reading it in
! parallel is therefore impossible without a redistribution step. Converting once
! to a fixed-record binary layout lets every rank seek straight to its own rows,
! which also removes the need for psb_matdist -- at 448 ranks that scatter cost
! 28 s of MPI time on its own.
!
! Layout, all 8-byte native-endian:
!
! header magic, version, ipk_bytes, rpk_bytes, nrows, ncols, nnz, reserved
! row_ptr nrows+1 entries, 1-based, row i occupies row_ptr(i)..row_ptr(i+1)-1
! col_idx nnz entries, 1-based global column indices
! values nnz doubles
!
! Native endianness on purpose: this is an intermediate for local runs, and the
! header records the kinds so a mismatched reader fails loudly instead of
! producing garbage.
!
program psb_mm2bin
use psb_base_mod
use psb_util_mod
implicit none
integer(psb_lpk_), parameter :: bin_magic = 6510534141251461121_psb_lpk_
integer(psb_lpk_), parameter :: bin_version = 1_psb_lpk_
type(psb_ldspmat_type) :: aux
type(psb_ld_coo_sparse_mat) :: coo
character(len=1024) :: infile, outfile
integer(psb_ipk_) :: info, iout
integer(psb_lpk_) :: nr, nc, nz, k, irow
integer(psb_lpk_), allocatable :: row_ptr(:)
real(psb_dpk_) :: t0, t1
if (command_argument_count() < 2) then
write(psb_err_unit,*) 'usage: psb_mm2bin <in.mtx> <out.psb>'
stop 1
end if
call get_command_argument(1, infile)
call get_command_argument(2, outfile)
call cpu_time(t0)
call mm_mat_read(aux, info, filename=trim(infile))
if (info /= psb_success_) then
write(psb_err_unit,*) 'cannot read ', trim(infile), ' info=', info
stop 2
end if
call cpu_time(t1)
write(psb_out_unit,'("read : ",f8.2," s")') t1-t0
nr = aux%get_nrows()
nc = aux%get_ncols()
! Move the storage into a COO object and sort it row-major: fix() is what
! guarantees the entries of a row are contiguous, which is the whole premise
! of the layout below.
call cpu_time(t0)
call aux%mv_to(coo)
call coo%fix(info)
if (info /= psb_success_) then
write(psb_err_unit,*) 'fix failed, info=', info
stop 3
end if
nz = coo%get_nzeros()
call cpu_time(t1)
write(psb_out_unit,'("sort : ",f8.2," s")') t1-t0
write(psb_out_unit,'("rows ",i0," cols ",i0," nonzeros ",i0)') nr, nc, nz
! row_ptr from the sorted row indices
allocate(row_ptr(nr+1), stat=info)
if (info /= 0) stop 4
row_ptr = 0
do k = 1, nz
irow = coo%ia(k)
row_ptr(irow+1) = row_ptr(irow+1) + 1
end do
row_ptr(1) = 1
do k = 1, nr
row_ptr(k+1) = row_ptr(k+1) + row_ptr(k)
end do
if (row_ptr(nr+1) /= nz+1) then
write(psb_err_unit,*) 'row_ptr inconsistent: ', row_ptr(nr+1), nz+1
stop 5
end if
call cpu_time(t0)
iout = 91
open(unit=iout, file=trim(outfile), access='stream', form='unformatted', &
& status='replace', action='write', iostat=info)
if (info /= 0) stop 6
write(iout) bin_magic, bin_version, &
& int(storage_size(nz)/8, psb_lpk_), int(storage_size(coo%val(1))/8, psb_lpk_), &
& nr, nc, nz, 0_psb_lpk_
write(iout) row_ptr
write(iout) coo%ja(1:nz)
write(iout) coo%val(1:nz)
close(iout)
call cpu_time(t1)
write(psb_out_unit,'("write : ",f8.2," s")') t1-t0
write(psb_out_unit,'("output : ",a)') trim(outfile)
end program psb_mm2bin