Re: Complexity explosion in linear solve
- To: mathgroup at smc.vnet.net
- Subject: [mg80264] Re: Complexity explosion in linear solve
- From: carlos at colorado.edu
- Date: Thu, 16 Aug 2007 04:43:31 -0400 (EDT)
- References: <f9s3f8$a8r$1@smc.vnet.net><f9udg9$ajv$1@smc.vnet.net>
On Aug 15, 2:28 am, Jens-Peer Kuska <ku... at informatik.uni-leipzig.de> wrote: > Hi, > > what do you expect ? matrix inversion is an N^3 process > when N is the dimension of the matrix. This mean for numerical > data, that you have to do N^3 operations to get the inverse. > For symbolic data a single operation gives not a single number > it simply add another leaf to your expression because 1+2 is 3 > (leaf count 1) but a+b is a+b (leaf count 2) so every of your N^3 > operations add a new leaf to your expression and you end up (in the > worst case) with N^3 more leafs in your expression .. > > Regards > Jens > > car... at colorado.edu wrote: > > LinearSolve (and cousins Inverse, LUDecomposition etc) works > > efficiently with > > numerical float data. But it is prone to complexity explosion in > > solving symbolic > > systems of moderate size. Here is an example from an actual > > application. > > > Task: solve A x = b for x where A is 16 x 16 and b is 16 x 1: > > > A={{7/12,1/24,-1/2,0,-1/4,-5/24,-1/2,0, > > 0,-1/24,-1/4,0,-1/3,5/24,-1/4,0},{1/24,91/144,-2/3,0, > > 5/24,-3/16,-1/3,0,-1/24,0,-1/3,0,-5/24,-4/9,-2/3, > > 0},{-1/2,-2/3,(-9+(1378*C11*hh)/225)/2,(-4*C13*hh)/15, > > 1/2,-1/3,(-3-2*C11*hh)/2,4*C13*hh,1/4, > > 1/3,(-3-(32*C11*hh)/25)/2,0,-1/4, > > 2/3,(-3-(128*C11*hh)/45)/2,(-64*C13*hh)/15},{0, > > 0,(-4*C13*hh)/15,(224*C22*hh)/5,0,0,-4*C13*hh,16*C23*hh,0,0,0, > > 16*C23*hh,0,0,(64*C13*hh)/15,(64*C23*hh)/5},{-1/4,5/24,1/2,0, > > 7/12,-1/24,1/2,0,-1/3,-5/24,1/4,0,0,1/24,1/4, > > 0},{-5/24,-3/16,-1/3,0,-1/24,91/144,-2/3,0,5/24,-4/9,-2/3, > > 0,1/24,0,-1/3,0},{-1/2,-1/3,(-3-2*C11*hh)/2,-4*C13*hh, > > 1/2,-2/3,(-9+(314*C11*hh)/45)/2,(4*C13*hh)/15,1/4, > > 2/3,(-3-(128*C11*hh)/45)/2,(64*C13*hh)/15,-1/4, > > 1/3,(-3-(32*C11*hh)/15)/2,0},{0,0,4*C13*hh,16*C23*hh,0, > > 0,(4*C13*hh)/15,(832*C22*hh)/15,0, > > 0,(-64*C13*hh)/15,(64*C23*hh)/5,0,0,0,(80*C23*hh)/3},{0,-1/24, > > 1/4,0,-1/3,5/24,1/4,0,7/12,1/24,1/2,0,-1/4,-5/24,1/2, > > 0},{-1/24,0,1/3,0,-5/24,-4/9,2/3,0,1/24,91/144,2/3,0, > > 5/24,-3/16,1/3,0},{-1/4,-1/3,(-3-(32*C11*hh)/25)/2,0, > > 1/4,-2/3,(-3-(128*C11*hh)/45)/2,(-64*C13*hh)/15,1/2, > > 2/3,(-9+(1378*C11*hh)/225)/2,(-4*C13*hh)/15,-1/2, > > 1/3,(-3-2*C11*hh)/2,4*C13*hh},{0,0,0,16*C23*hh,0, > > 0,(64*C13*hh)/15,(64*C23*hh)/5,0,0,(-4*C13*hh)/15,(224*C22*hh)/5, > > 0,0,-4*C13*hh,16*C23*hh},{-1/3,-5/24,-1/4,0,0,1/24,-1/4, > > 0,-1/4,5/24,-1/2,0,7/12,-1/24,-1/2,0},{5/24,-4/9,2/3,0, > > 1/24,0,1/3,0,-5/24,-3/16,1/3,0,-1/24,91/144,2/3, > > 0},{-1/4,-2/3,(-3-(128*C11*hh)/45)/2,(64*C13*hh)/15, > > 1/4,-1/3,(-3-(32*C11*hh)/15)/2,0,1/2, > > 1/3,(-3-2*C11*hh)/2,-4*C13*hh,-1/2, > > 2/3,(-9+(314*C11*hh)/45)/2,(4*C13*hh)/15},{0, > > 0,(-64*C13*hh)/15,(64*C23*hh)/5,0,0,0,(80*C23*hh)/3,0,0, > > 4*C13*hh,16*C23*hh,0,0,(4*C13*hh)/15,(832*C22*hh)/15}} > > > b={0,0,0,0,q*b/2,0,0,0,q*b/2,0,0,0,0,0,0,0} > > > in which all symbolic variables are atoms. Matrix A has rank 13. > > Three BCs: x1=x2=x13=0 are imposed by removing equations 1, 2 and > > 13, > > which gives the reduced 13-system Ahat xhat = bhat. Matrix Ahat is > > symmetric but indefinite, so Cholesky is not an option. Instead > > LinearSolve is used > > > xhat=LinearSolve[Ahat,bhat] > > > Solving takes about 3 mins on a dual G5 running 5.2. LeafCounts of > > the xhat entries reach hundreds of millions: > > > {325197436,235675446,292306655,256982512,146153454,73076907,35324210, > > 18877804,9441784,4745440,2429139,1354800,903023} > > > Obviously Simplify would take a long, long time so I didnt attempt it. > > Another solution method, however, gave this solution in about 10 sec: > > > xhat={q/3, 0, 4*q, 0, q/3, 0, 4*q, 0, q/3, 0, 0, q/3, 0} > > > My question is: is there a way to tell LinearSolve to Simplify as it > > goes > > along? That would preclude or at least alleviate the leafcount > > explosion. Hi Jens, I know *that* - I teach numerical analysis, among other stuff. Also I teach the kids to avoid Cramer's rule for systems of order >3. Then can you explain me why MyLinearSolve[A_,b_]:=Module[{i,j,Ab,detA,detAb,n=Length[A],x}, detA=Simplify[Det[A]]; x=Table[0,{n}]; For [i=1,i<=n,i++, Ab=Transpose[A]; Ab[[i]]=b; detAb=Simplify[Det[Ab]]; x[[i]]=Simplify[detAb/detA] ]; ClearAll[Ab]; Return[x]]; beats LinearSolve for a symbolic system of order 13? Here are my times on a dual G5 Mac: LinearSolve 235 sec unsimplified x MyLinraSolve 6 sec simplified x For the simplification of the LinearSolve solution x I estimate 10^8 years.