Re: Speed of dot product in Mathematica
- To: mathgroup at smc.vnet.net
- Subject: [mg5317] Re: [mg5230] Speed of dot product in Mathematica
- From: carlos at mars.Colorado.EDU (Carlos A. Felippa)
- Date: Wed, 27 Nov 1996 01:47:35 -0500
- Organization: University of Colorado, Boulder
- Sender: owner-wri-mathgroup at wolfram.com
In article <56p225$m3c at dragonfly.wolfram.com> BobHanlon at aol.com writes:
>You can also use
>
>Plus @@ (a b)
>
>If you want the dot product of portions of the factors, use only the
>portions. For example,
>
>Take[a, {100, 1000}].Take[b, {200, 1100}]
>
>
Of course I use Take where possible. The following skyline
factorization module is an example.
SymmSkyMatrixFactor[S_,tol_]:= Module[
{p,a,fail,i,j,k,l,m,n,ii,ij,jj,jk,jmj,d,s,row,v},
row=SymmSkyMatrixRowLengths[S]; s=Max[row];
{p,a}=S; n=Length[p]-1; v=Table[0,{n}]; fail=0;
Do [jj=p[[j+1]]; If [jj<0|row[[j]]==0, Continue[]]; d=a[[jj]];
jmj=Abs[p[[j]]]; jk=jj-jmj;
Do [i=j-jk+k; v[[k]]=0; ii=p[[i+1]];
If [ii<0, Continue[]]; m=Min[ii-Abs[p[[i]]],k]-1;
ij=jmj+k; v[[k]]=a[[ij]];
v[[k]]-=Take[a,{ii-m,ii-1}].Take[v,{k-m,k-1}];
a[[ij]]=v[[k]]*a[[ii]],
{k,1,jk-1}];
d-=Take[a,{jmj+1,jmj+jk-1}].Take[v,{1,jk-1}];
If [Abs[d]<tol*row[[j]], fail=j; a[[jj]]=Infinity; Break[] ];
a[[jj]]=1/d,
{j,1,n}];
Return[{{p,a},fail}]
];
However, the dot product cannot always be done using unit strides
in both arrays, and use of Sum, Do or For exert big speed penalties.
To overcome that problem I once tried the form
but that was rejected as illegal. It seems an unreasonable
restriction for the Take function.