Mr. Roman Pearce

## 1683 Reputation

19 years, 296 days
CECM/SFU
Research Associate

I am a research associate at Simon Fraser University and a member of the Computer Algebra Group at the CECM.

## Sparse Matrix Products in Maple...

Maple

I thought I would share some code for computing sparse matrix products in Maple.  For floating point matrices this is done quickly, but for algebraic datatypes there is a performance problem with the builtin routines, as noted in http://www.mapleprimes.com/questions/205739-How-Do-I-Improve-The-Performance-Of

The code is fairly straightforward in that it uses op(1,A) to extract the dimensions and op(2,A) to extract the non-zero elements of a Matrix or Vector, and then loops over those elements.  I included a sparse map function for cases where you want to map a function (like expand) over non-zero elements only.

`# sparse matrix vector productspmv := proc(A::Matrix,V::Vector)local m,n,Ae,Ve,Vi,R,e;  n, m := op(1,A);  if op(1,V) <> m then error "incompatible dimensions"; end if;  Ae := op(2,A);  Ve := op(2,V);  Vi := map2(op,1,Ve);  R := Vector(n, storage=sparse);  for e in Ae do    n, m := op(1,e);    if member(m, Vi) then R[n] := R[n] + A[n,m]*V[m]; end if;  end do;  return R;end proc:`
`# sparse matrix productspmm := proc(A::Matrix, B::Matrix)local m,n,Ae,Be,Bi,R,l,e,i;  n, m := op(1,A);  i, l := op(1,B);  if i <> m then error "incompatible dimensions"; end if;  Ae := op(2,A);  Be := op(2,B);  R := Matrix(n,l,storage=sparse);  for i from 1 to l do    Bi, Be := selectremove(type, Be, (anything,i)=anything);    Bi := map2(op,[1,1],Bi);    for e in Ae do      n, m := op(1,e);      if member(m, Bi) then R[n,i] := R[n,i] + A[n,m]*B[m,i]; end if;    end do;  end do;  return R;end proc:`
`# sparse mapsmap := proc(f, A::{Matrix,Vector})local B, Ae, e;  if A::Vector then    B := Vector(op(1,A),storage=sparse):  else    B := Matrix(op(1,A),storage=sparse):  end if;  Ae := op(2,A);  for e in Ae do    B[op(1,e)] := f(op(2,e),args[3..nargs]);  end do;  return B;end proc:`

As for how it performs, here is a demo inspired by the original post.

`n := 674;k := 6;`
`A := Matrix(n,n,storage=sparse):for i to n do  for j to k do    A[i,irem(rand(),n)+1] := randpoly(x):  end do:end do:`
`V := Vector(n):for i to k do  V[irem(rand(),n)+1] := randpoly(x):end do:`
`C := CodeTools:-Usage( spmv(A,V) ):  # 7ms, 25x fasterCodeTools:-Usage( A.V ):  # 174 ms`
`B := Matrix(n,n,storage=sparse):for i to n do  for j to k do    B[i,irem(rand(),n)+1] := randpoly(x):  end do:end do:`
`C := CodeTools:-Usage( spmm(A,B) ):  # 2.74 sec, 50x fasterCodeTools:-Usage( A.B ):  # 2.44 min`
`# expand and collect like termsC := CodeTools:-Usage( smap(expand, C) ):`

## Maple 17: High Performance Polynomials...

Maple

As described on the help page ?updates,Maple17,Performance, Maple 17 uses a new data structure for polynomials with integer coefficients. Our goal was to improve the performance and parallel speedup of polynomial algorithms that underpin much of the system and create a platform for large scale polynomial computations. Shown below is the new representation for 9xy3z

## Use solid state disks for swap...

A lot of people want to do very large computations which require a lot of RAM.  But above a certain threshold, the cost of memory explodes.  We had this idea but no excuse to try it.  Buy a good SSD and use it for a swap partition.  I suggest the Intel X25-M.  Make sure you have at least 1/10th of its size in RAM, i.e. 8GB of RAM for an 80GB drive, but of course more is better.  The RAM should act as a massive cache for the SSD, giving you another...

## import branch structure of discussions?...

Is it possible to import the branch structure of discussions from the old MaplePrimes?  In many cases it looks like all the replies are put into a flat list.

## Parallel Sparse Polynomial Multiplicatio...

Maple 14

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.

 1 2 3 4 5 6 7 Last Page 1 of 9
﻿