suggestion for matrix algebra library?

I'm currently looking for a free software matrix algebra library written in C (or with C interfaces),
preferably if it supports indefinitely large matrices. Currently I'm considering the GNU scientific library
but as it was the only library I found and it appears that it's matrix data type is only able to provide
matrices with at most (size_t)-1 elements then I believe I may be missing something.

So, what else is there? What libraries do you recommend?


Thanks in advance,
Rui Maciel

Rui Maciel

Rui Maciel wrote:
> I'm currently looking for a free software matrix algebra library written in C (or with C interfaces),
> preferably if it supports indefinitely large matrices. Currently I'm considering the GNU scientific library
> but as it was the only library I found and it appears that it's matrix data type is only able to provide
> matrices with at most (size_t)-1 elements then I believe I may be missing something.

Do you know of a 64 bit machine with more than (size_t)-1 bytes of memory?

--
Ian Collins

Ian Collins

Ian Collins wrote:

> Do you know of a 64 bit machine with more than (size_t)-1 bytes of memory?

You have a point. I was trying to avoid any arbitrary limits on the matrices' size and at the time I didn't
knew exactly how big a size_t could get and the 2^32 value seemed a bit low for some uses. But it turns out
x86-64 saved the day.


Thanks for the help,
Rui Maciel

Rui Maciel

Malcolm McLean wrote:

> If your matrices are so big they won't fit into memory you need to code
> routines to access from disk specially. I would look at obtaining a
> Beowulf cluster and doing the calculations in parallel.

Thanks for the tip, Malcolm. You have a point. As soon as it is necessary to partition the matrix the maximum
matrix size ceases to be an issue.


> However the matrices really do have to be very big for that to happen.

Unfortunately they have to be. I'm currently designing a small finite element method program to solve certain
3D structural analysis problems and it would be nice if it was possible to avoid having to set a limit on how
many elements the models may have.


Rui Maciel

Rui Maciel

On 2009-11-01, Rui Maciel <> wrote:
> Malcolm McLean wrote:
>
>> If your matrices are so big they won't fit into memory you need to code
>> routines to access from disk specially.
>
> Unfortunately they have to be. I'm currently designing a small finite
> element method program to solve certain 3D structural analysis
> problems and it would be nice if it was possible to avoid having to
> set a limit on how many elements the models may have.

A [ (size_t) times (size_t) ] matrix is probably much larger
than any conceivable need in finite elements. If your
matrix reaches that limit, then you will have so many other
problems that the matrix storage itself will be the least of
your worries.

The main reason for the popularity of the finite element
methods is that the coefficient matrices obtained from it
are _sparse_. This is because each elements interacts only
with its neighboring elements and none other. As a results,
the coefficient matrix consists mostly of zeros. Storing all
those zeros is a waste --- all you need are the (relatively few)
nonzero entries.

Look up the UMFPACK package (it's written in standard C) for
storing and manipulating sparse matrices.

--
Rouben Rostamian

rouben

Rui Maciel a écrit :
> I'm currently looking for a free software matrix algebra library written in C (or with C interfaces),
> preferably if it supports indefinitely large matrices. Currently I'm considering the GNU scientific library
> but as it was the only library I found and it appears that it's matrix data type is only able to provide
> matrices with at most (size_t)-1 elements then I believe I may be missing something.
>
> So, what else is there? What libraries do you recommend?
>
>
> Thanks in advance,
> Rui Maciel

lcc-win ships with the Meschach matrix library. It is included in the standard distribution.
If you want the source code please drop me an email.

Here is the index. I think it is fairly complete

Linear algebra functions

In the descriptions below, matrices are represented by capital letters, vectors by lower case
letters and scalars by alpha.

Function Description

band2mat Convert band matrix to dense matrix
bd_free Deallocate (destroy) band matrix
bd_get Allocate and initialize band matrix
bd_transp Transpose band matrix
bd_resize Resize band matrix
bdLDLfactor Band LDL^T factorization
bdLDLsolve Solve Ax=b using band LDL^T factors
bdLUfactor Band LU factorization
bdLUsolve Solve Ax=b using band LU factors
bisvd SVD of bi-diagonal matrix
BKPfactor Bunch-Kaufman-Parlett factorization
BKPsolve Bunch-Kaufman-Parlett solver
catch() Catch a raised error (macro)
catchall() Catch any raised error (macro)
catch_FPE() Catch floating point error (sets flag)
CHfactor Dense Cholesky factorization
CHsolve Cholesky solver
d_save Save real in MATLAB format
Dsolve Solve Dx=y , D diagonal
ERRABORT() Abort on error (sets flag, macro)
ERREXIT() Exit on error (sets flag, macro)
error() Raise an error (macro, see ev_err())
err_list_attach Attach new list of errors
err_list_free Discard list of errors
err_is_list_attached Checks for an error list
ev_err Raise an error (function)
fft Computes Fast Fourier Transform
finput Input a simple data item from a stream
fprompter Print prompt to stderr
get_col Extract a column from a matrix
get_row Extract a row from a matrix
givens Compute Givens parameters
hhtrcols Compute AP^T where P is a Householder matrix
hhtrrows Compute PA where P is a Householder matrix
hhtrvec Compute Px where P is a Householder matrix
hhvec Compute parameters for a Householder matrix
ifft Computes inverse FFT
in_prod Inner product of vectors
input Input a simple data item from stdin (macro)
iter_arnoldi Arnoldi iterative method
iter_arnoldi_iref Arnoldi iterative method with refinement
iter_ATx Set A^T in ITER structure
iter_Ax Set A in ITER structure
iter_Bx Set preconditioner in ITER structure
iter_cg Conjugate gradients iterative method
iter_cgne Conjugate gradients for normal equations
iter_cgs CGS iterative method
iter_copy Copy ITER data structures
iter_copy2 Shallow copy of ITER data structures
iter_dump Dump ITER data structure to a stream
iter_free Free (deallocate) ITER structure
iter_get Allocate ITER structure
iter_gmres GMRES iterative method
iter_lsqr LSQR iterative method
iter_mgcr MGCR iterative method
iter_resize Resize vectors in an ITER data structure
iter_spcg Sparse matrix CG method
iter_spcgne Sparse matrix CG method for normal equations
iter_spcgs Sparse matrix CGS method
iter_spgmres Sparse matrix GMRES method
iter_splsqr Sparse matrix LSQR method
iter_spmgcr Sparse matrix MGCR method
iv_add Add integer vectors
iv_copy Copy integer vector
iv_dump Dump integer vector to a stream
iv_finput Input integer vector from a stream
iv_foutput Output integer vector to a stream
IV_FREE Free (deallocate) an integer vector (macro)
iv_free Free (deallocate) integer vector (function)
iv_free_vars Free a list of integer vectors
iv_get Allocate and initialise an integer vector
iv_get_vars Allocate list of integer vectors
iv_input Input integer vector from stdin (macro)
iv_output Output integer vector to stdout (macro)
iv_resize Resize an integer vector
iv_resize_vars Resize a list of integer vectors
iv_sub Subtract integer vectors
LDLfactor LDL^T factorization
LDLsolve LDL^T solver
LDLupdate Update LDL^T factorization
Lsolve Solve Lx=y , L lower triangular
LTsolve Solve L^Tx=y , L lower triangular
LUcondest Estimate a condition number using LU factors
LUfactor Compute LU factors with implicit scaled partial pivoting
LUsolve Solve Ax=b using LU factors
LUTsolve Solve A^Tx=b usng LU factors
m_add Add matrices
makeQ Form Q matrix for QR factorization
makeR Form R matrix for QR factorization
mat2band Extract band matrix from dense matrix
MCHfactor Modified Cholesky factorization
m_copy Copy dense matrix
m_dump Dump matrix data structure to a stream
m_exp Computes matrix exponential
_m_exp Matrix exponential
m_finput Input matrix from a stream
m_foutput Output matrix to a stream
M_FREE Free (deallocate) a matrix (macro)
m_free Free (deallocate) matrix (function)
m_free_vars Free a list of matrices
m_get Allocate and initialize a matrix
m_get_vars Allocate list of matrices
m_ident Sets matrix to identity matrix
m_input Input matrix from stdin (macro)
m_inverse Invert matrix
m_load Load matrix in MATLAB format
m_mlt Multiplies matrices
mmtr_mlt Computes AB^T
m_norm1 Computes ||A||_1 of a matrix
m_norm_frob Computes the Frobenius norm of a matrix
m_norm_inf Computes ||A||_inf of a matrix
m_ones Set matrix to all 1's
m_output Output matrix to stdout (macro)
m_poly Computes a matrix polynomial
m_pow Computes integer power of a matrix
mrand Generates pseudo-random real number
m_rand Randomise entries of a matrix
mrandlist Generates array of pseudo-random numbers
m_resize Resize matrix
m_resize_vars Resize a list of matrices
m_save Save matrix in MATLAB format
m_sub Subtract matrices
m_transp Transpose matrix
mtrm_mlt Computes A^TB
mv_mlt Computes Ax
mv_mltadd Computes y <- Ax+y
m_zero Zero a matrix
ON_ERROR() Error handler (macro)
prompter Print prompt message to stdout
px_cols Permute the columns of a matrix
px_copy Copy permutation
px_dump Dump permutation data structure to a stream
px_finput Input permutation from a stream
px_foutput Output permutation to a stream
PX_FREE Free (deallocate) a permutation (macro)
px_free Free (deallocate) permutation (function)
px_free_vars Free a list of permutations
px_get Allocate and initialize a permutation
px_get_vars Allocate a list of permutations
px_ident Sets permutation to identity
px_input Input permutation from stdin (macro)
px_inv Invert permutation
pxinv_vec Computes P^Tx where P is a permutation matrix
pxinv_zvec Computes P^Tx where P is a permutation matrix (complex)
px_mlt Multiply permutations
px_output Output permutation to stdout (macro)
px_resize Resize a permutation
px_resize_vars Resize a list of permutations
px_rows Permute the rows of a matrix
px_sign Returns the sign of the permutation
px_transp Transpose a pair of entries
px_vec Computes Px where P is a permutation matrix
px_zvec Computes Px where P is a permutation matrix (complex)
QRCPfactor QR factorization with column pivoting
QRfactor QR factorization
QRCPsolve Solve Ax=b
QRsolve Solve Ax=b using QR factorization
QRTsolve Solve A^Tx=b using QR factorization
QRupdate Update explicit QR factors
rot_cols Apply Givens rotation to the columns of a matrix
rot_rows Apply Givens rotation to the rows of a matrix
rot_vec Apply Givens rotation to a vector
rot_zvec Apply complex Givens rotation to a vector
schur Compute real Schur form
schur_evals Compute eigenvalues from the real Schur form
schur_vecs Compute eigenvectors from the real Schur form
set_col Set the column of a matrix to a given vector
set_err_flag Control behavior of ev_err()
set_row Set the row of a matrix to a given vector
sm_mlt Scalar-matrix multiplication
smrand Set seed for mrand()
spBKPfactor Sparse symmetric indefinite factorization
spBKPsolve Sparse symmetric indefinite solver
spCHfactor Sparse Cholesky factorization
spCHsolve Sparse Cholesky solver
spCHsymb Symbolic sparse Cholesky factorization (no floating point operations)
sp_col_access Sets up column access paths for a sparse matrix
sp_compact Eliminates zero entries in a sparse matrix
sp_copy Copies a sparse matrix
sp_copy2 Copies a sparse matrix into another
sp_diag_access Sets up diagonal access paths for a sparse matrix
sp_dump Dump sparse matrix data structure to a stream
sp_finput Input sparse matrix from a stream
sp_foutput Output a sparse matrix to a stream
sp_free Free (deallocate) a sparse matrix
sp_get Allocate and initialize a sparse matrix
sp_get_val Get the (i,j) entry of a sparse matrix
spICHfactor Sparse incomplete Cholesky factorization
sp_input Input a sparse matrix form stdin
spLUfactor Sparse LU factorization using partial pivoting
spLUsolve Solves Ax=b using sparse LU factors
spLUTsolve Solves A^Tx=b using sparse LU factors
sp_mv_mlt Computes Ax for sparse A
sp_output Outputs a sparse matrix to a stream (macro)
sp_resize Resize a sparse matrix
sprow_add Adds a pair of sparse rows
sprow_foutput Output sparse row to a stream
sprow_get Allocate and initialize a sparse row
sprow_idx Get location of an entry in a sparse row
sprow_merge Merge two sparse rows
sprow_mltadd Sparse row vector multiply-and-add
sprow_set_val Set an entry in a sparse row
sprow_smlt Multiplies a sparse row by a scalar
sprow_sub Subtracts a sparse row from another
sprow_xpd Expand a sparse row
sp_set_val Set the (i,j) entry of a sparse matrix
sp_vm_mlt Compute x^TA for sparse A
sp_zero Zero (but do not remove) all entries of a sparse matrix
svd The SVD of a matrix
sv_mlt Scalar-vector multiply
symmeig Eigenvalues/vectors of a symmetric matrix
tracecatch() Catch and re-raise errors (macro)
trieig Eigenvalues/vectors of a symmetric tridiagonal matrix
Usolve Solve Ux=b where U is upper triangular
UTsolve Solve U^Tx=b where U is upper triangular
v_add Add vectors
v_conv Convolution product of vectors
v_copy Copy vector
v_dump Dump vector data structure to a stream
v_finput Input vector from a stream
v_foutput Output vector to a stream
V_FREE Free (deallocate) a vector (macro)
v_free Free (deallocate) vector (function)
v_free_vars Free a list of vectors
v_get Allocate and initialize a vector
v_get_vars Allocate list of vectors
v_input Input vector from stdin (macro)
v_lincomb Compute sum of a_i x_i for an array of vectors
v_linlist Compute sum of a_i x_i for a list of vectors
v_map Apply function component wise to a vector
v_max Computes max vector entry and index
v_min Computes min vector entry and index
v_mltadd Computes y <- alpha*x+y for vectors x , y
vm_mlt Computes x^TA
vm_mltadd Computes y^T <- y^T+x^TA
v_norm1 Computes ||x||_1 for a vector
v_norm2 Computes ||x||_2 (the Euclidean norm) of a vector
v_norm_inf Computes ||x||_inf for a vector
v_ones Set vector to all 1's
v_output Output vector to stdout (macro)
v_pconv Periodic convolution of two vectors
v_rand Randomize entries of a vector
v_resize Resize a vector
v_resize_vars Resize a list of vectors
v_save Save a vector in MATLAB format
v_slash Computes component wise ratio of vectors
v_sort Sorts vector components
v_star Component wise vector product
v_sub Subtract two vectors
v_sum Sum of components of a vector
v_zero Zero a vector
z_finput Read complex number from file or stream
z_foutput Prints complex number to file or stream
zgivens Compute complex Givens' rotation
zhhtrcols Apply Householder transformation: PA (complex)
zhhtrrows Apply Householder transformation: AP (complex)
zhhtrvec Apply Householder transformation: Px (complex)
zhhvec Compute Householder transformation
zin_prod Complex inner product
z_input Read complex number from stdin
zLAsolve Solve L^*x=b , L complex lower triangular
zLsolve Solve Lx=b , L complex lower triangular
zLUAsolve Solve A^*x=b using complex LU factorization
(A^* - adjoin of A, A is complex)
zLUcondest Complex LU condition estimate
zLUfactor Complex LU factorization
zLUsolve Solve Ax=b using complex LU factorization
zm_add Add complex matrices
zm_adjoint Computes adjoin of complex matrix
zmakeQ Construct Q matrix for complex QR
zmakeR Construct R matrix for complex QR
zmam_mlt Computes A^*B (complex)
zm_dump Dump complex matrix to stream
zm_finput Input complex matrix from stream
ZM_FREE Free (deallocate) complex matrix (macro)
zm_free Free (deallocate) complex matrix (function)
zm_free_vars Free a list of complex matrices
zm_get Allocate complex matrix
zm_get_vars Allocate a list of complex matrices
zm_input Input complex matrix from stdin
zm_inverse Compute inverse of complex matrix
zm_load Load complex matrix in MATLAB format
zmma_mlt Computes AB^* (complex)
zm_mlt Multiply complex matrices
zm_norm1 Complex matrix 1-norm
zm_norm_frob Complex matrix Frobenius norm
zm_norm_inf Complex matrix infinity-norm
zm_rand Randomize complex matrix
zm_resize Resize complex matrix
zm_resize_vars Resize a list of complex matrices
zm_save Save complex matrix in MATLAB format
zm_sub Subtract complex matrices
zmv_mlt Complex matrix-vector multiply
zmv_mltadd Complex matrix-vector multiply and add
zm_zero Zero complex matrix
zQRCPfactor Complex QR factorization with column pivoting
zQRCPsolve Solve Ax = b using complex QR factorization
zQRfactor Complex QR factorization
zQRAsolve Solve A^*x = b using complex QR factorization
zQRsolve Solve Ax = b using complex QR factorization
zrot_cols Complex Givens' rotation of columns
zrot_rows Complex Givens' rotation of rows
z_save Save complex number in MATLAB format
zschur Complex Schur factorization
zset_col Set column of complex matrix
zset_row Set row of complex matrix
zsm_mlt Complex scalar-matrix product
zUAsolve Solve U^*x=b , U complex upper triangular
zUsolve Solve Ux=b , U complex upper triangular
zv_add Add complex vectors
zv_copy Copy complex vector
zv_dump Dump complex vector to a stream
zv_finput Input complex vector from a stream
ZV_FREE Free (deallocate) complex vector (macro)
zv_free Free (deallocate) complex vector (function)
zv_free_vars Free a list of complex vectors
zv_get Allocate complex vector
zv_get_vars Allocate a list of complex vectors
zv_input Input complex vector from a stdin
zv_lincomb Compute sum of a_i x_i for an array of vectors
zv_linlist Compute sum of a_i x_i for a list of vectors
zv_map Apply function to each element of a complex vector
zv_mlt Complex scalar-vector product
zv_mltadd Complex scalar-vector multiply and add
zvm_mlt Computes A^*x (complex)
zvm_mltadd Computes A^*x+y (complex)
zv_norm1 Complex vector 1-norm vnorm1()
zv_norm2 Complex vector 2-norm (Euclidean norm)
zv_norm_inf Complex vector infinity- (or supremum) norm
zv_rand Randomise complex vector
zv_resize Resize complex vector
zv_resize_vars Resize a list of complex vectors
zv_save Save complex vector in MATLAB format
zv_slash Componentwise ratio of complex vectors
zv_star Componentwise product of complex vectors
zv_sub Subtract complex vectors
zv_sum Sum of components of a complex vector
zv_zero Zero complex vector

jacob navia

rouben wrote:

> A [ (size_t) times (size_t) ] matrix is probably much larger
> than any conceivable need in finite elements. If your
> matrix reaches that limit, then you will have so many other
> problems that the matrix storage itself will be the least of
> your worries.
>
> The main reason for the popularity of the finite element
> methods is that the coefficient matrices obtained from it
> are _sparse_. This is because each elements interacts only
> with its neighboring elements and none other. As a results,
> the coefficient matrix consists mostly of zeros. Storing all
> those zeros is a waste --- all you need are the (relatively few)
> nonzero entries.

Yes, that's true. At the moment the only matrix package that I was familiar with was GML, which makes a point
at not implementing a sparse matrix data type. So as it appears that GML stores their matrices as a long
malloc'ed chunk then, if we consider a 32-bit platform and that their elements happened to be doubles, it's
size wouldn't end up being very impressive.


> Look up the UMFPACK package (it's written in standard C) for
> storing and manipulating sparse matrices.

Thanks for the tip, rouben. At the moment I'm not able to test it as the Ubuntu repositories don't offer the
-dev package yet. Still, it appears to be very promising.


Once again thanks for the tip. Kudos!
Rui Maciel

Rui Maciel

On 2009-11-01, Rui Maciel <> wrote:
> rouben wrote:
>
>> Look up the UMFPACK package (it's written in standard C) for storing
>> and manipulating sparse matrices.
>
> Thanks for the tip, rouben. At the moment I'm not able to test it as
> the Ubuntu repositories don't offer the -dev package yet. Still, it
> appears to be very promising.

You can get the UMFPACK source from:

http://www.cise.ufl.edu/research/sparse/umfpack/

It takes some work to set it up but the installation instructions
are quite detailed and helpful.

For a first attempt, compile and install UMFPACK with the "no BLAS"
option (see the installation instructions). This will produce
a working UMFPACK library but will not be optimized for speed.
In the future you may add an optimized BLAS library to take
advantage of the specific features of your hardware.

--
Rouben Rostamian

Rouben Rostamian

Rui Maciel wrote:

> I'm currently looking for a free software matrix algebra library
> written in C (or with C interfaces), preferably if it supports
> indefinitely large matrices. Currently I'm considering the GNU
> scientific library but as it was the only library I found and it
> appears that its matrix data type is only able to provide matrices
> with at most (size_t)-1 elements then I believe I may be missing
> something.
>
> So, what else is there? What libraries do you recommend?

You might find some of these links relevant.

http://en.wikipedia.org/wiki/Basic_L...ra_Subprograms
http://en.wikipedia.org/wiki/Automat...gebra_Software
http://en.wikipedia.org/wiki/LAPACK

Noob

In article <4aecb8db$0$26824$>,
says...
>
> I'm currently looking for a free software matrix algebra library written in C (or with C interfaces),
> preferably if it supports indefinitely large matrices. Currently I'm considering the GNU scientific library
> but as it was the only library I found and it appears that it's matrix data type is only able to provide
> matrices with at most (size_t)-1 elements then I believe I may be missing something.
>
> So, what else is there? What libraries do you recommend?

Here is a standard BLAS:
http://www.netlib.org/blas/

You can get better results with tuned BLAS versions. For instance,
Intel's:
http://software.intel.com/en-us/intel-mkl/

AMD's:
http://developer.amd.com/cpu/Librari...s/default.aspx

Atlas (tunes itself, free of charge):
http://math-atlas.sourceforge.net/

If you really are confounded by (size_t) -1 being too small, almost for
sure you need a sparse package.
See also:
http://www.mathcom.com/corpdir/techinfo.mdir/q207.html
http://www.netlib.org/scalapack/
http://math.nist.gov/spblas/

If (for instance) you have a dense matrix, and it fills up a 64 bit
size_t and then some, then solving the LU system will require many
universe lifetimes on the fastest computer in the world.

Dann Corbit