Items tagged with sparse sparse Tagged Items Feed

Hi everybody, I have a rather large code in Maple for my research which needs fast eigensolving for large sparse matrices, thus, I am using the Matlab link to parse my matrices and get them eigensolved in Matlab. The problem is that, I am unable to directly parse a Maple sparse matrix to a Matlab sparse matrix which substantially slows down the calculation.

This is the best I can do so far:

 

> H:=Matrix(5000,5000,shape=symmetric,storage=sparse,datatype=float[8]):

I am dissapointed. I recently upgraded to Maple 15 and I bought a license for the NAG libraries thinking that I could use NAG mark 8 chapter F12 to solve my Eigenvalue problem with large sparse matrix, in reality, there seems to be no implementation of these NAG routines in Maple. is there a way to get around this problem? I know I can also connect to MATLAB to take advantage of the integrated sparse handling but I want to avoid buying another piece of software(you know we scientist aren't millionaries).

I'm trying to compute the determinant of a 19x19 matrix, each with entries consisting of multivariate polynomials (polys in 6 vars, of degree at most 8).  Specifically I want to factor the determinant.  Running this on Maple 15 quickly gives an "object too large" error in the determinant function, after using only 700 MB of RAM in about 2 minutes.  My machine has 3GB of RAM, so I am confused as to why this object is "too large". 

 

I...

Maple NAG connector

April 07 2011 by Alger 359 Maple 14

Hi,

How to load Maple-NAG connector to execut this example ?

Is it an additional program to install ?

n := 4:
nnz := 11:
la := 22:
lfill := 1:
dtol := 0:
pstrat := "Nag_SparseNsym_CompletePiv":
milu := "Nag_SparseNsym_UnModFact":
a := Vector[row](la,[1.0,1.0,-1.0,2.0,2.0,3.0,-2.0,1.0,-2.0,1.0,1.0],datatype=float[8]):
irow := Vector[row](22,[1,1,2,2,2,3,3,4,4,4,4],datatype=integer[kernelopts('wordsize')/8]):
icol := Vector[row...

I am using finite element analysis and have 47008 nodes, so require a 47008 by 47008 matrix.  I wish to create the matrix using a procedure and using if statements to fill in the matrix with the equations (mostly repeating and every element in each row is 0 excpet for 4 or 5).  Anyway, I have tested my solution method on a much smalller matrix and know it works.  My question is if there is a way to create such a large matrix without getting a memory error, and...

Our previous article described the design of fast algorithms for multiplying and dividing sparse polynomials. We have integrated these algorithms into the expand and divide commands of Maple 14. In this post I want to talk a bit about what you might see when you try Maple 14. Keep in mind that the product isn't released yet and I don't work for Maplesoft, so general disclaimers apply. Nevertheless, one of the first things you may notice is this.

task manager with maple 14

In our previous article we described a packed representation for sparse polynomials is designed for scalability and high performance. The expand and divide commands in Maple 14 use this representation internally to multiply and divide polynomials with integer coefficients, converting to and from Maple's generic data structure described here. In this post I want to show you how these algorithms work and why they are fast. It's a critical stepping stone for our next topic, which is parallelization.

sdmp multiplication

Our first article introduced Maple's polynomial data structure and explained how Maple spends a lot of time working with monomials. To multiply polynomials having n and m terms, Maple must construct, simplify, hash, and sort all nm pairwise products to determine what monomials are equal. This work is performed even if the result has far fewer than nm terms, making it a rather inefficient way to multiply large multivariate polynomials. This article describes a new data structure for multivariate polynomials that is being added to Maple for a future release.

sdmp packed arrays

9xyz  -  4yz  -  6xyz  -  8x  -  5

This is the first in a series of short informal articles about our efforts to speed up polynomial arithmetic in Maple. We begin with an example of how polynomials are represented in Maple right now.

Maple sum structure

9xyz  -  4yz  -  6xyz  -  8x  -  5

When you enter a polynomial in Maple, it creates a generic data structure like the one above. In Maple's representation this polynomial is a sum of terms that is 11 words of memory long where each word is either 32 or 64 bits. For each term it stores a pointer to a monomial followed by a coefficient.

I was reminded of this by another thread.

It is faster to add in-place a large size storage=sparse float[8] Matrix into a new empty storage=rectangular float[8] Matrix than it is to convert it that way using the Matrix() or rtable() constructors.

Here's an example. First I'll do it with in-place Matrix addition. And then after that with a call to Matrix(). I measure the time to execute as well as the increase in bytes-allocated and bytes-used.

> with(LinearAlgebra):

> N := 500:
> A := RandomMatrix(N,'density'=0.1,
>                   'outputoptions'=['storage'='sparse',
>                                    'datatype'=float[8]]):

> st,ba,bu := time(),kernelopts(bytesalloc),kernelopts(bytesused):

> B := Matrix(N,'datatype'=float[8]):
> MatrixAdd(B,A,'inplace'=true):

> time()-st,kernelopts(bytesalloc)-ba,kernelopts(bytesused)-bu;
                            0.022, 2489912, 357907

"I've seen this element before..." Often we are faced with the problem of building up sets incrementally, by removing pieces one at a time from a larger whole. The bottlenecks in this case are usually: 1) adding a small set X to a large set S (copies S and X, making this ~O(|S|+|X|)) 2) removing elements of the large set S from the small set X (binary search: |X|*log(|S|)) A classic example of this is a breadth-first-search. We start at one vertex of a graph and in each iteration we add the set of new neighbors X to the set of vertices S that have already been found. We can make this more useful by making the program return the sets of new neighbors found in each iteration, that is, the sets of vertices that are distance 1, 2, 3, etc. from the initial vertex.

When working with large sparse linear systems you often want to look at their non-zero structure, however Maple's existing tools are all designed for dense matrices. I wrote a little tool to produce images like this in reasonable time. You can download the code here, and the rest of this post is a quick tutorial on how to use the included command. Maple 11 is required.

What is the largest linear system that Maple can solve? You might be surprised to find out. In this article we present strategies for solving sparse linear systems over the rationals. An example implementation is provided, but first we present a bit of background. Sparse linear systems arise naturally from problems in mathematics, science, and engineering. Typically many quantities are related, but because of an underlying structure only a small subset of the elements appear in most equations. Consider networks, finite element models, structural analysis problems, and linear programming problems.

This tip comes care of Dr. Michael Monagan at Simon Fraser University. Represent your sparse matrix as a list of rows, and represent each row as a linear equation in an indexed name. For example:

A := [[1,0,3],[2,0,0],[0,4,5]];

S := [ 1*x[1] + 3*x[3], 2*x[1], 4*x[2]+5*x[3] ];

To compute the product of the matrix A with a Vector X, assign x[i] := V[i] and evaluate. This can be done inside of a procedure because x is a table.

V := [7,8,9]: for i to 3 do x[i...

This is a quick programming exercise to correct the following problem in Maple:

for n from 8 to 12 do 
  A := Matrix(2^n, 2^n, storage='sparse'):  # zero matrix
  print( 2^n, time(LinearAlgebra:-Transpose(A)) );
end do:

The problem is that the LinearAlgebra:-Transpose command is not sparse. That is, the time it takes is proportional to the overall size of the matrix and not to the number of non-zero entries - even when your matrix uses sparse storage. In this post we will look at what is required to program a new Transpose command which can handle much larger matrices.

Page 1 of 1