If you use Python and SciPy, the scipy.sparse allows you to use single-precision floats. And you can wrap Python stuff into Matlab, if you choose to go that way.

However, even if you use single precision, you may not gain a factor two in memory. On my computer (64 bit Matlab on Linux), an NxN sparse matrix A in Matlab has a storage cost of 8*(2*nnz(A) + N + 1) bytes, since you store one int64 and one double for each numeric value. For a factor of two savings, you would have to make sure that the integer size used is int32, and that the maximum value of that is enough for you (which it probably is).

My recommendation for you would be to just buy more memory for your machine, or find another machine to run on, rather than trying to reimplement stuff for the sake of saving some space. Memory is cheap nowadays, while your time is limited.

Alternatively, you could implement a routine for computing the matrix-vector product of your system on the fly. That way you never have to worry about memory again, at least not until your geometry gets huge. The performance hit might be an issue though, depending on how costly your matrix assembly is.

Dear Kevin,

It seems the matlab mex library has actually an sparse implementation, for general numeric matrices, this function is called

mxCreateSparseNumericMatrix

which I can't find it in the documentation but if you do (Assuming you are using linux)

objdump -T libmx.so | c++filt > /tmp/mexLib

(where libmx should be available in /matlabroot/bin/glnxa64) you can see that the libmx actually provides this function, so you might be able to transfer the arrays containing the values with single precision to a mex function then do your calculations using C or Fortran and then transfer back the results, of course you will probabily need to spend a few night awake but it's doable.

Meanwhile, just out of curiosity, why you are running out of memory? assuming that you grid is connected to 6 neighboring grid you sparse presentation will only have 7*N^3 (approx) double elements and 56*N^3 bytes (fjust for doubles) and you can easily accomodate maybe ten million grid points on a few gigabytes of RAM,

If you really need to solve a larger problem you should think of parallelizing your code and run it on a cluster rather than a single workstation, since even if you can somehow manage to fit say a 100 million grid points in memory (by buying memory of using float instead of double and optimising your code etc), it is almost certainly not feasible to run such huge calculations on a single CPU

Goodluck

Sina

I suggest you using sparse2 from SuiteSparse and/or suite_sparse from MILAMIN.

These are functions that improve significantly the performance of the standard sparse Matlab function.

There is not additional work to do since these packages come with precompiled mex-functions. 

Kind regards!

http://faculty.cse.tamu.edu/davis/suitesparse.html

http://milamin.org/

As you mentioned, MATLAB sparse solvers do not appear to allow single precision. SuiteSparse mentioned above also appears to be limited to either complex or double precision numbers. Several direct solvers have the ability to run in single precision (MUMPS, PARADISO, SUPERLU, etc) The problem is direct solvers generally require a lot more memory than iterative/krylov solvers. Some codes (WSMP, TAUCS, etc) get around this by providing out-of-core variants. There has been recent work on both direct and iterative mixed precision solvers (i.e., most of operations in single precision, but residual updates in double precision). The work seemed to at least initially be motivated by GPU computing. Older GPU's performed single precision operations significantly faster (i.e., 6X I believe) than double precision operations. For more recent GPU's the speed-up is only 2X, which is about the same as for CPU'. Some iterative solver packages such as PETSc could have single precision capabilities. A good starting point for available software is

http://www.netlib.org/utk/people/JackDongarra/la-sw.html

If memory is such a bottle-neck, then you may want to look at the Jacobian-Free Newton Krylov work of Prof David Keyes OR Prof Rainald Lohner. Lohner was a co-author of  "A Fast, Matrix-free Implicit Method for Compressible Flows on Unstructured Grids".  The basic idea is that Krylov methods do not really need the Jacobian matrix itself, but merely products of the Jacobian matrix times a vector which can be approximated with scaled changes to the residual. Prof T C Kelley actually provides simple MATLAB code for Newton-Krylov methods  with his book on "Solving Nonlinear Equations with Newton Methods"

http://www4.ncsu.edu/~ctk/newtony.html

The disadvantage of the method is you typically require more GMRES/Krylov iterations per solve. The advantage is that one does not have to form and/or store the Jacobian AND the approach eliminates the need for an outer iteration for nonlinear equations.

The undocumented mxCreateSparseNumericArray is probably not a good idea. See this StackOverflow question and answer:

http://stackoverflow.com/questions/20552421/how-to-create-a-single-float-sparse-matrix-in-mex-files

There is SuiteSparse, but it'd probably be better to just get more memory, assuming you're using 64-bit Matlab. If your Matlab install is 32-bit, then that's probably the root of your problem.