phynexis-deps/armadillo-code
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

initialise arrays to zeros

Browse Source
This commit is contained in:
conrad
2022-03-14 16:10:00 +10:00
parent 392fd04ad4
commit 44bc0686e8
1 changed files with 22 additions and 13 deletions
Download Patch File Download Diff File Expand all files Collapse all files
include/armadillo_bits/sp_auxlib_meat.hpp
+22 -13
View File
@@ -482,6 +482,8 @@ sp_auxlib::eigs_sym_arpack(Col<eT>& eigval, Mat<eT>& eigvec, const SpMat<eT>& X,
podarray<blas_int> select(ncv); // Logical array of dimension NCV.
blas_int ldz = n;
select.zeros();
// seupd() will output directly into the eigval and eigvec objects.
eigval.zeros(n_eigvals);
eigvec.zeros(n, n_eigvals);
@@ -843,7 +845,7 @@ sp_auxlib::eigs_gen_arpack(Col< std::complex<T> >& eigval, Mat< std::complex<T>
}
if(info != 0) { return false; }
// The process has converged, and now we need to recover the actual eigenvectors using neupd().
blas_int rvec = 1; // .TRUE
blas_int nev = blas_int(n_eigvals);
@@ -858,9 +860,11 @@ sp_auxlib::eigs_gen_arpack(Col< std::complex<T> >& eigval, Mat< std::complex<T>
blas_int ldz = n;
podarray<T> workev(3 * ncv);
select.zeros();
dr.zeros();
di.zeros();
z.zeros();
workev.zeros();
arpack::neupd(&rvec, &howmny, select.memptr(), dr.memptr(), di.memptr(), z.memptr(), &ldz, (T*) &sigmar, (T*) &sigmai, workev.memptr(), &bmat, &n, which, &nev, &tol, resid.memptr(), &ncv, v.memptr(), &ldv, iparam.memptr(), ipntr.memptr(), workd.memptr(), workl.memptr(), &lworkl, rwork.memptr(), &info);
@@ -1090,6 +1094,11 @@ sp_auxlib::eigs_gen(Col< std::complex<T> >& eigval, Mat< std::complex<T> >& eigv
blas_int ldz = n;
podarray<std::complex<T>> workev(2 * ncv);
select.zeros();
d.zeros();
z.zeros();
workev.zeros();
// Prepare the outputs; neupd() will write directly to them.
eigval.zeros(n_eigvals);
eigvec.zeros(n, n_eigvals);
@@ -1952,7 +1961,7 @@ sp_auxlib::run_aupd_plain
n = X.n_rows; // The size of the matrix (should already be set outside).
blas_int nev = n_eigvals;
resid.set_size(n);
resid.zeros(n);
// Two contraints on NCV: (NCV > NEV) for sym problems or
// (NCV > NEV + 2) for gen problems and (NCV <= N)
@@ -1965,8 +1974,8 @@ sp_auxlib::run_aupd_plain
if(ncv < (nev + (sym ? 1 : 3))) { ncv = (nev + (sym ? 1 : 3)); }
if(ncv > n ) { ncv = n; }
v.set_size(n * ncv); // Array N by NCV (output).
rwork.set_size(ncv); // Work array of size NCV for complex calls.
v.zeros(n * ncv); // Array N by NCV (output).
rwork.zeros(ncv); // Work array of size NCV for complex calls.
ldv = n; // "Leading dimension of V exactly as declared in the calling program."
// IPARAM: integer array of length 11.
@@ -1976,16 +1985,16 @@ sp_auxlib::run_aupd_plain
iparam(6) = 1; // Mode 1: A * x = lambda * x.
// IPNTR: integer array of length 14 (output).
ipntr.set_size(14);
ipntr.zeros(14);
// Real work array used in the basic Arnoldi iteration for reverse communication.
workd.set_size(3 * n);
workd.zeros(3 * n);
// lworkl must be at least 3 * NCV^2 + 6 * NCV.
lworkl = 3 * (ncv * ncv) + 6 * ncv;
// Real work array of length lworkl.
workl.set_size(lworkl);
workl.zeros(lworkl);
info = 0; // Set to 0 initially to use random initial vector.
@@ -2124,7 +2133,7 @@ sp_auxlib::run_aupd_shiftinvert
n = X.n_rows; // The size of the matrix (should already be set outside).
blas_int nev = n_eigvals;
resid.set_size(n);
resid.zeros(n);
// Two contraints on NCV: (NCV > NEV) for sym problems or
// (NCV > NEV + 2) for gen problems and (NCV <= N)
@@ -2137,8 +2146,8 @@ sp_auxlib::run_aupd_shiftinvert
if(ncv < (nev + (sym ? 1 : 3))) { ncv = (nev + (sym ? 1 : 3)); }
if(ncv > n ) { ncv = n; }
v.set_size(n * ncv); // Array N by NCV (output).
rwork.set_size(ncv); // Work array of size NCV for complex calls.
v.zeros(n * ncv); // Array N by NCV (output).
rwork.zeros(ncv); // Work array of size NCV for complex calls.
ldv = n; // "Leading dimension of V exactly as declared in the calling program."
// IPARAM: integer array of length 11.
@@ -2151,16 +2160,16 @@ sp_auxlib::run_aupd_shiftinvert
iparam(6) = 3; // Mode 3: A * x = lambda * M * x, M symmetric semi-definite. OP = inv[A - sigma*M]*M (A complex) or Real_Part{ inv[A - sigma*M]*M } (A real) and B = M.
// IPNTR: integer array of length 14 (output).
ipntr.set_size(14);
ipntr.zeros(14);
// Real work array used in the basic Arnoldi iteration for reverse communication.
workd.set_size(3 * n);
workd.zeros(3 * n);
// lworkl must be at least 3 * NCV^2 + 6 * NCV.
lworkl = 3 * (ncv * ncv) + 6 * ncv;
// Real work array of length lworkl.
workl.set_size(lworkl);
workl.zeros(lworkl);
info = 0; // Set to 0 initially to use random initial vector.