Fix a bug that manifests while running the test suite of 'python-igraph': https://github.com/opencollab/arpack-ng/issues/401 https://github.com/opencollab/arpack-ng/pull/414 From d885b7be4ecdc9c1496f2d6f256f6c0d34962459 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Szabolcs=20Horva=CC=81t?= Date: Sun, 9 Apr 2023 16:36:35 +0200 Subject: [PATCH] fix: ensure that LAPACK RNG state is propagated - fixes #401, #410, #411 - restores 'inits' variable removed in ce2e69a849da1d10dad5d6d3ec4db6120b3ecf50, ensuring that the RNG state is propagated - reverts e0d67054f573da351f12a226f7c7cc65a690ef3d to ensure that seed is different on each parallel thread - updates seed initialization of parallel pdgetv0/psgetv0 so that they match that of pzgetv0/pcgetv0 --- PARPACK/SRC/MPI/pcgetv0.f | 48 +++++++++++++++++++++++---------------- PARPACK/SRC/MPI/pdgetv0.f | 40 ++++++++++++++++++++++++++------ PARPACK/SRC/MPI/psgetv0.f | 43 ++++++++++++++++++++++++++--------- PARPACK/SRC/MPI/pzgetv0.f | 48 +++++++++++++++++++++++---------------- SRC/cgetv0.f | 21 ++++++++++++----- SRC/dgetv0.f | 21 ++++++++++++----- SRC/sgetv0.f | 21 ++++++++++++----- SRC/zgetv0.f | 21 ++++++++++++----- 8 files changed, 183 insertions(+), 80 deletions(-) diff --git a/PARPACK/SRC/MPI/pcgetv0.f b/PARPACK/SRC/MPI/pcgetv0.f index 59e3d1658..24fe8a0f1 100644 --- a/PARPACK/SRC/MPI/pcgetv0.f +++ b/PARPACK/SRC/MPI/pcgetv0.f @@ -176,13 +176,13 @@ subroutine pcgetv0 c | Local Scalars & Arrays | c %------------------------% c - logical first, orth + logical first, inits, orth integer idist, iseed(4), iter, msglvl, jj, myid, igen Real & rnorm0 Complex & cnorm, cnorm2 - save first, iseed, iter, msglvl, orth, rnorm0 + save first, iseed, inits, iter, msglvl, orth, rnorm0 c Complex & cnorm_buf, buf2(1) @@ -203,6 +203,12 @@ subroutine pcgetv0 & ccdotc external ccdotc, pscnorm2, slapy2 c +c %-----------------% +c | Data Statements | +c %-----------------% +c + data inits /.true./ +c c %-----------------------% c | Executable Statements | c %-----------------------% @@ -213,26 +219,30 @@ subroutine pcgetv0 c | random number generator | c %-----------------------------------% c + if (inits) then c -c %-----------------------------------% -c | Generate a seed on each processor | -c | using process id (myid). | -c | Note: the seed must be between 1 | -c | and 4095. iseed(4) must be odd. | -c %-----------------------------------% +c %-----------------------------------% +c | Generate a seed on each processor | +c | using process id (myid). | +c | Note: the seed must be between 1 | +c | and 4095. iseed(4) must be odd. | +c %-----------------------------------% c - call MPI_COMM_RANK(comm, myid, ierr) - igen = 1000 + 2*myid + 1 - if (igen .gt. 4095) then - write(0,*) 'Error in p_getv0: seed exceeds 4095!' - end if + call MPI_COMM_RANK(comm, myid, ierr) + igen = 1000 + 2*myid + 1 + if (igen .gt. 4095) then + write(0,*) 'Error in p_getv0: seed exceeds 4095!' + end if +c + iseed(1) = igen/1000 + igen = mod(igen,1000) + iseed(2) = igen/100 + igen = mod(igen,100) + iseed(3) = igen/10 + iseed(4) = mod(igen,10) c - iseed(1) = igen/1000 - igen = mod(igen,1000) - iseed(2) = igen/100 - igen = mod(igen,100) - iseed(3) = igen/10 - iseed(4) = 7 + inits = .false. + end if c if (ido .eq. 0) then c diff --git a/PARPACK/SRC/MPI/pdgetv0.f b/PARPACK/SRC/MPI/pdgetv0.f index 0f348b820..5a1956997 100644 --- a/PARPACK/SRC/MPI/pdgetv0.f +++ b/PARPACK/SRC/MPI/pdgetv0.f @@ -177,11 +177,11 @@ subroutine pdgetv0 c | Local Scalars & Arrays | c %------------------------% c - logical first, orth - integer idist, iseed(4), iter, msglvl, jj + logical first, inits, orth + integer idist, iseed(4), iter, msglvl, jj, myid, igen Double precision & rnorm0, buf2(1) - save first, iseed, iter, msglvl, orth, rnorm0 + save first, iseed, inits, iter, msglvl, orth, rnorm0 c Double precision & rnorm_buf @@ -206,6 +206,12 @@ subroutine pdgetv0 c intrinsic abs, sqrt c +c %-----------------% +c | Data Statements | +c %-----------------% +c + data inits /.true./ +c c %-----------------------% c | Executable Statements | c %-----------------------% @@ -216,10 +222,30 @@ subroutine pdgetv0 c | random number generator | c %-----------------------------------% c - iseed(1) = 1 - iseed(2) = 3 - iseed(3) = 5 - iseed(4) = 7 + if (inits) then +c +c %-----------------------------------% +c | Generate a seed on each processor | +c | using process id (myid). | +c | Note: the seed must be between 1 | +c | and 4095. iseed(4) must be odd. | +c %-----------------------------------% +c + call MPI_COMM_RANK(comm, myid, ierr) + igen = 1000 + 2*myid + 1 + if (igen .gt. 4095) then + write(0,*) 'Error in p_getv0: seed exceeds 4095!' + end if +c + iseed(1) = igen/1000 + igen = mod(igen,1000) + iseed(2) = igen/100 + igen = mod(igen,100) + iseed(3) = igen/10 + iseed(4) = mod(igen,10) +c + inits = .false. + end if c if (ido .eq. 0) then c diff --git a/PARPACK/SRC/MPI/psgetv0.f b/PARPACK/SRC/MPI/psgetv0.f index d79a513b2..078e4fa8c 100644 --- a/PARPACK/SRC/MPI/psgetv0.f +++ b/PARPACK/SRC/MPI/psgetv0.f @@ -177,11 +177,11 @@ subroutine psgetv0 c | Local Scalars & Arrays | c %------------------------% c - logical first, orth - integer idist, iseed(4), iter, msglvl, jj + logical first, inits, orth + integer idist, iseed(4), iter, msglvl, jj, myid, igen Real & rnorm0 - save first, iseed, iter, msglvl, orth, rnorm0 + save first, iseed, inits, iter, msglvl, orth, rnorm0 c Real & rnorm_buf @@ -206,20 +206,41 @@ subroutine psgetv0 c intrinsic abs, sqrt c +c %-----------------% +c | Data Statements | +c %-----------------% +c + data inits /.true./ +c c %-----------------------% c | Executable Statements | c %-----------------------% c c -c %-----------------------------------% -c | Initialize the seed of the LAPACK | -c | random number generator | -c %-----------------------------------% + if (inits) then c - iseed(1) = 1 - iseed(2) = 3 - iseed(3) = 5 - iseed(4) = 7 +c %-----------------------------------% +c | Generate a seed on each processor | +c | using process id (myid). | +c | Note: the seed must be between 1 | +c | and 4095. iseed(4) must be odd. | +c %-----------------------------------% +c + call MPI_COMM_RANK(comm, myid, ierr) + igen = 1000 + 2*myid + 1 + if (igen .gt. 4095) then + write(0,*) 'Error in p_getv0: seed exceeds 4095!' + end if +c + iseed(1) = igen/1000 + igen = mod(igen,1000) + iseed(2) = igen/100 + igen = mod(igen,100) + iseed(3) = igen/10 + iseed(4) = mod(igen,10) +c + inits = .false. + end if c if (ido .eq. 0) then c diff --git a/PARPACK/SRC/MPI/pzgetv0.f b/PARPACK/SRC/MPI/pzgetv0.f index 731fb319f..94fb705f3 100644 --- a/PARPACK/SRC/MPI/pzgetv0.f +++ b/PARPACK/SRC/MPI/pzgetv0.f @@ -176,13 +176,13 @@ subroutine pzgetv0 c | Local Scalars & Arrays | c %------------------------% c - logical first, orth + logical first, inits, orth integer idist, iseed(4), iter, msglvl, jj, myid, igen Double precision & rnorm0 Complex*16 & cnorm, cnorm2 - save first, iseed, iter, msglvl, orth, rnorm0 + save first, iseed, inits, iter, msglvl, orth, rnorm0 c Complex*16 & cnorm_buf, buf2(1) @@ -203,6 +203,12 @@ subroutine pzgetv0 & zzdotc external zzdotc , pdznorm2 , dlapy2 c +c %-----------------% +c | Data Statements | +c %-----------------% +c + data inits /.true./ +c c %-----------------------% c | Executable Statements | c %-----------------------% @@ -213,26 +219,30 @@ subroutine pzgetv0 c | random number generator | c %-----------------------------------% c + if (inits) then c -c %-----------------------------------% -c | Generate a seed on each processor | -c | using process id (myid). | -c | Note: the seed must be between 1 | -c | and 4095. iseed(4) must be odd. | -c %-----------------------------------% +c %-----------------------------------% +c | Generate a seed on each processor | +c | using process id (myid). | +c | Note: the seed must be between 1 | +c | and 4095. iseed(4) must be odd. | +c %-----------------------------------% c - call MPI_COMM_RANK(comm, myid, ierr) - igen = 1000 + 2*myid + 1 - if (igen .gt. 4095) then - write(0,*) 'Error in p_getv0: seed exceeds 4095!' - end if + call MPI_COMM_RANK(comm, myid, ierr) + igen = 1000 + 2*myid + 1 + if (igen .gt. 4095) then + write(0,*) 'Error in p_getv0: seed exceeds 4095!' + end if +c + iseed(1) = igen/1000 + igen = mod(igen,1000) + iseed(2) = igen/100 + igen = mod(igen,100) + iseed(3) = igen/10 + iseed(4) = mod(igen,10) c - iseed(1) = igen/1000 - igen = mod(igen,1000) - iseed(2) = igen/100 - igen = mod(igen,100) - iseed(3) = igen/10 - iseed(4) = 7 + inits = .false. + end if c if (ido .eq. 0) then c diff --git a/SRC/cgetv0.f b/SRC/cgetv0.f index b49e66708..c231eadcb 100644 --- a/SRC/cgetv0.f +++ b/SRC/cgetv0.f @@ -156,13 +156,13 @@ subroutine cgetv0 c | Local Scalars & Arrays | c %------------------------% c - logical first, orth + logical first, inits, orth integer idist, iseed(4), iter, msglvl, jj Real & rnorm0 Complex & cnorm - save first, iseed, iter, msglvl, orth, rnorm0 + save first, iseed, inits, iter, msglvl, orth, rnorm0 c c %----------------------% c | External Subroutines | @@ -180,6 +180,12 @@ subroutine cgetv0 & ccdotc external ccdotc, scnrm2, slapy2 c +c %-----------------% +c | Data Statements | +c %-----------------% +c + data inits /.true./ +c c %-----------------------% c | Executable Statements | c %-----------------------% @@ -190,10 +196,13 @@ subroutine cgetv0 c | random number generator | c %-----------------------------------% c - iseed(1) = 1 - iseed(2) = 3 - iseed(3) = 5 - iseed(4) = 7 + if (inits) then + iseed(1) = 1 + iseed(2) = 3 + iseed(3) = 5 + iseed(4) = 7 + inits = .false. + end if c if (ido .eq. 0) then c diff --git a/SRC/dgetv0.f b/SRC/dgetv0.f index 8be4fa26d..1d6dc01bd 100644 --- a/SRC/dgetv0.f +++ b/SRC/dgetv0.f @@ -157,11 +157,11 @@ subroutine dgetv0 c | Local Scalars & Arrays | c %------------------------% c - logical first, orth + logical first, inits, orth integer idist, iseed(4), iter, msglvl, jj Double precision & rnorm0 - save first, iseed, iter, msglvl, orth, rnorm0 + save first, iseed, inits, iter, msglvl, orth, rnorm0 c c %----------------------% c | External Subroutines | @@ -183,6 +183,12 @@ subroutine dgetv0 c intrinsic abs, sqrt c +c %-----------------% +c | Data Statements | +c %-----------------% +c + data inits /.true./ +c c %-----------------------% c | Executable Statements | c %-----------------------% @@ -193,10 +199,13 @@ subroutine dgetv0 c | random number generator | c %-----------------------------------% c - iseed(1) = 1 - iseed(2) = 3 - iseed(3) = 5 - iseed(4) = 7 + if (inits) then + iseed(1) = 1 + iseed(2) = 3 + iseed(3) = 5 + iseed(4) = 7 + inits = .false. + end if c if (ido .eq. 0) then c diff --git a/SRC/sgetv0.f b/SRC/sgetv0.f index 26130a014..d861b2d6d 100644 --- a/SRC/sgetv0.f +++ b/SRC/sgetv0.f @@ -157,11 +157,11 @@ subroutine sgetv0 c | Local Scalars & Arrays | c %------------------------% c - logical first, orth + logical first, inits, orth integer idist, iseed(4), iter, msglvl, jj Real & rnorm0 - save first, iseed, iter, msglvl, orth, rnorm0 + save first, iseed, inits, iter, msglvl, orth, rnorm0 c c %----------------------% c | External Subroutines | @@ -183,6 +183,12 @@ subroutine sgetv0 c intrinsic abs, sqrt c +c %-----------------% +c | Data Statements | +c %-----------------% +c + data inits /.true./ +c c %-----------------------% c | Executable Statements | c %-----------------------% @@ -193,10 +199,13 @@ subroutine sgetv0 c | random number generator | c %-----------------------------------% c - iseed(1) = 1 - iseed(2) = 3 - iseed(3) = 5 - iseed(4) = 7 + if (inits) then + iseed(1) = 1 + iseed(2) = 3 + iseed(3) = 5 + iseed(4) = 7 + inits = .false. + end if c if (ido .eq. 0) then c diff --git a/SRC/zgetv0.f b/SRC/zgetv0.f index cc13c3cfb..1fbd50851 100644 --- a/SRC/zgetv0.f +++ b/SRC/zgetv0.f @@ -156,13 +156,13 @@ subroutine zgetv0 c | Local Scalars & Arrays | c %------------------------% c - logical first, orth + logical first, inits, orth integer idist, iseed(4), iter, msglvl, jj Double precision & rnorm0 Complex*16 & cnorm - save first, iseed, iter, msglvl, orth, rnorm0 + save first, iseed, inits, iter, msglvl, orth, rnorm0 c c %----------------------% c | External Subroutines | @@ -180,6 +180,12 @@ subroutine zgetv0 & zzdotc external zzdotc, dznrm2, dlapy2 c +c %-----------------% +c | Data Statements | +c %-----------------% +c + data inits /.true./ +c c %-----------------------% c | Executable Statements | c %-----------------------% @@ -190,10 +196,13 @@ subroutine zgetv0 c | random number generator | c %-----------------------------------% c - iseed(1) = 1 - iseed(2) = 3 - iseed(3) = 5 - iseed(4) = 7 + if (inits) then + iseed(1) = 1 + iseed(2) = 3 + iseed(3) = 5 + iseed(4) = 7 + inits = .false. + end if c if (ido .eq. 0) then c