MathGroup Archive 2007

[Date Index] [Thread Index] [Author Index]

Search the Archive

Complexity explosion in linear solve

  • To: mathgroup at smc.vnet.net
  • Subject: [mg80181] Complexity explosion in linear solve
  • From: carlos at colorado.edu
  • Date: Tue, 14 Aug 2007 07:07:02 -0400 (EDT)

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.



  • Prev by Date: Re: Reading Coordinates from Plot
  • Next by Date: Re: Re: hardware for Mathematica 6.0
  • Previous by thread: Re: problem using Ersek's RootSearch
  • Next by thread: Re: Complexity explosion in linear solve