1

Solving linear systems: Which technique is fastest?

 2 years ago
source link: https://blogs.sas.com/content/iml/2011/08/17/solving-linear-systems-which-technique-is-fastest.html
Go to the source link to view the article. You can view the picture content, updated content and better typesetting reading experience. If the link is broken, please click the button below to view the snapshot at that time.

Solving linear systems: Which technique is fastest?

16

I've previously described ways to solve systems of linear equations, A*b = c. While discussing the relative merits of the solving a system for a particular right hand side versus solving for the inverse matrix, I made the assertion that it is faster to solve a particular system than it is to compute an inverse and use the inverse to solve the system.

In particular, in terms of the SAS/IML SOLVE function and INV function, I asserted that it is faster to run b = solve(A,c); than Ainv = inv(A); b = Ainv * c;.

A colleague asked a good question: "How much faster?"

The following SAS/IML program answers this question. The program generates a random n x n matrix for a range of values for n. For each matrix, the program times how long it takes to solve the linear system. The program is adapted from Chapter 15 of Wicklin (2010), Statistical Programming with SAS/IML Software. 

proc iml;
/* this program computes the solution to a linear system in two 
   different ways and compares the performance of each method */
size = T(do(100, 1000, 100)); /* 100, 200, ... 1000 */
results = j(nrow(size), 2);   /* allocate room for results */
do i = 1 to nrow(size);
   n = size[i];
   A = rannor(j(n,n,1));           /* n x n matrix */
   b = rannor(j(n,1,1));           /* n x 1 vector */
 
   /* use the INV function to solve a linear system Ax=b */
   t0 = time();                    /* begin timing INV */
      AInv = inv(A);               /* compute inverse of A */
      x = AInv*b;                  /* solve linear equation A*x=b */
   results[i,1] = time() - t0;     /* end timing */
 
   /* use the SOLVE function to solve the same linear system */
   t0 = time();                    /* begin timing SOLVE */
      x = solve(A,b);              /* solve linear equation directly */
   results[i,2] = time() - t0;     /* end timing */
end;

You can save the times to a SAS data set and use the SGPLOT procedure to compare the performance of the two methods: 

/* write results to a data set */
y = size || results;
create Performance from y[colname={"Size" "INV" "SOLVE"}];
append from y;
quit;
 
title "Time Required to Solve a Linear System: INV versus SOLVE";
title2 "From Wicklin (2010), Statistical Programming with SAS/IML Software";
proc sgplot data=Performance;
  series x=Size y=INV / curvelabel;
  series x=Size y=SOLVE /curvelabel;
  yaxis grid label="Time (s)";
  xaxis grid label="Size of Matrix";
run;

For a 1000 x 1000 matrix, it takes about 0.8 seconds to solve the system by computing the matrix inverse, whereas it takes 0.2 seconds to solve the system directly. That ratio is fairly typical: it takes about four times longer to solve a linear system with INV as with SOLVE.

This result is not unique to SAS/IML software. Although the algorithm that is used to compute each solution will affect the shape of the curves, solving a linear system directly should be faster than a solution that involves computing a matrix inverse.


report

Recommend

About Joyk


Aggregate valuable and interesting links.
Joyk means Joy of geeK