diff liboctave/UMFPACK/UMFPACK/OCTAVE/umfpack_solve.m @ 5164:57077d0ddc8e

[project @ 2005-02-25 19:55:24 by jwe]
author jwe
date Fri, 25 Feb 2005 19:55:28 +0000
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/liboctave/UMFPACK/UMFPACK/OCTAVE/umfpack_solve.m	Fri Feb 25 19:55:28 2005 +0000
@@ -0,0 +1,97 @@
+function x = umfpack_solve (arg1, op, arg2, Control)
+% UMFPACK_SOLVE
+%
+% x = umfpack_solve (A, '\', b, Control)
+% x = umfpack_solve (b, '/', A, Control)
+%
+% Computes x = A\b, or b/A, where A is square.  Uses UMFPACK if A is sparse.
+% The Control argument is optional.
+%
+% See also umfpack, umfpack_make, umfpack_details, umfpack_report,
+% and umfpack_simple.
+
+% UMFPACK Version 4.3 (Jan. 16, 2004), Copyright (c) 2004 by Timothy A.
+% Davis.  All Rights Reserved.  Type umfpack_details for License.
+
+%-------------------------------------------------------------------------------
+% check inputs and get default control parameters
+%-------------------------------------------------------------------------------
+
+if (op == '\\')
+    A = arg1 ;
+    b = arg2 ;
+elseif (op == '/')
+    A = arg2 ;
+    b = arg1 ;
+else
+    help umfack_solve
+    error ('umfpack_solve:  unrecognized operator') ;
+end
+
+[m n] = size (A) ;
+if (m ~= n)
+    help umfpack_solve
+    error ('umfpack_solve:  A must be square') ;
+end
+
+[m1 n1] = size (b) ;
+if ((op == '\\' & n ~= m1) | (op == '/' & n1 ~= m))
+    help umfpack_solve
+    error ('umfpack_solve:  b has the wrong dimensions') ;
+end
+
+if (nargin < 4)
+    Control = umfpack ;
+end
+
+%-------------------------------------------------------------------------------
+% solve the system
+%-------------------------------------------------------------------------------
+
+if (op == '\\')
+
+    if (~issparse (A))
+
+	% A is not sparse, so just use MATLAB
+	x = A\b ;
+
+    elseif (n1 == 1 & ~issparse (b))
+
+	% the UMFPACK '\' requires b to be a dense column vector
+	x = umfpack (A, '\\', b, Control) ;
+
+    else
+
+	% factorize with UMFPACK and do the forward/back solves in MATLAB 
+	[L, U, P, Q, R] = umfpack (A, Control) ;
+	keyboard
+	x = Q * (U \ (L \ (P * (R \ b)))) ;
+
+    end
+
+else
+
+    if (~issparse (A))
+
+	% A is not sparse, so just use MATLAB
+	x = b/A ;
+
+    elseif (m1 == 1 & ~issparse (b))
+
+	% the UMFPACK '/' requires b to be a dense column vector
+	x = umfpack (b, '/', A, Control) ;
+
+    else
+
+	% factorize with UMFPACK and do the forward/back solves in MATLAB
+	% this mimics the behavior of x = b/A, except for the row scaling
+	[L, U, P, Q, R] = umfpack (A.', Control) ;
+	x = (Q * (U \ (L \ (P * (R \ (b.')))))).' ;
+
+	% an alternative method:
+	% [L, U, P, Q, r] = umfpack (A, Control) ;
+	% x = (R \ (P' * (L.' \ (U.' \ (Q' * b.'))))).' ;
+
+    end
+
+end