[Scipy-svn] r2152 - in trunk/Lib/linsolve: . umfpack

scipy-svn at scipy.org scipy-svn at scipy.org
Wed Aug 9 04:47:11 CDT 2006


Author: rc
Date: 2006-08-09 04:47:04 -0500 (Wed, 09 Aug 2006)
New Revision: 2152

Modified:
   trunk/Lib/linsolve/__init__.py
   trunk/Lib/linsolve/info.py
   trunk/Lib/linsolve/linsolve.py
   trunk/Lib/linsolve/umfpack/__init__.py
   trunk/Lib/linsolve/umfpack/info.py
   trunk/Lib/linsolve/umfpack/umfpack.py
Log:
updated docs

Modified: trunk/Lib/linsolve/__init__.py
===================================================================
--- trunk/Lib/linsolve/__init__.py	2006-08-05 20:24:43 UTC (rev 2151)
+++ trunk/Lib/linsolve/__init__.py	2006-08-09 09:47:04 UTC (rev 2152)
@@ -2,6 +2,10 @@
 
 from info import __doc__
 
+import umfpack
+__doc__ = '\n\n'.join( (__doc__,  umfpack.__doc__) )
+del umfpack
+
 from linsolve import *
 
 __all__ = filter(lambda s:not s.startswith('_'),dir())

Modified: trunk/Lib/linsolve/info.py
===================================================================
--- trunk/Lib/linsolve/info.py	2006-08-05 20:24:43 UTC (rev 2151)
+++ trunk/Lib/linsolve/info.py	2006-08-09 09:47:04 UTC (rev 2152)
@@ -2,5 +2,50 @@
 Linear Solvers
 ==============
 
+The default solver is SuperLU (included in the scipy distribution), which can
+solve real or complex linear systems in both single and double precisions.  It
+is automatically replaced by UMFPACK, if available. Note that UMFPACK works in
+double precision only, so switch it off by
+>>> use_solver( use = {'useUmfpack': False} )
+to solve in the single precision.
+
+Example session:
+
+>>> from scipy.sparse import csc_matrix
+>>> from numpy import array
+>>> from scipy.linsolve import spdiags, spsolve, use_solver
+>>>
+>>> print "Inverting a sparse linear system:"
+>>> print "The sparse matrix (constructed from diagonals):"
+>>> a = spdiags([[1, 2, 3, 4, 5], [6, 5, 8, 9, 10]], [0, 1], 5, 5)
+>>> b = array([1, 2, 3, 4, 5])
+>>> print "Solve: single precision complex:"
+>>> use_solver( use = {'useUmfpack' : False} )
+>>> a = a.astype('F')
+>>> x = spsolve(a, b)
+>>> print x
+>>> print "Error: ", a*x-b
+>>>
+>>> print "Solve: double precision complex:"
+>>> use_solver( use = {'useUmfpack' : True} )
+>>> a = a.astype('D')
+>>> x = spsolve(a, b)
+>>> print x
+>>> print "Error: ", a*x-b
+>>>
+>>> print "Solve: double precision:"
+>>> a = a.astype('d')
+>>> x = spsolve(a, b)
+>>> print x
+>>> print "Error: ", a*x-b
+>>>
+>>> print "Solve: single precision:"
+>>> use_solver( use = {'useUmfpack' : False} )
+>>> a = a.astype('f')
+>>> x = spsolve(a, b.astype('f'))
+>>> print x
+>>> print "Error: ", a*x-b
+
+
 """
 postpone_import = 1

Modified: trunk/Lib/linsolve/linsolve.py
===================================================================
--- trunk/Lib/linsolve/linsolve.py	2006-08-05 20:24:43 UTC (rev 2151)
+++ trunk/Lib/linsolve/linsolve.py	2006-08-09 09:47:04 UTC (rev 2152)
@@ -90,22 +90,23 @@
                  diag_pivot_thresh, drop_tol, relax, panel_size)
 
 def _testme():
-    from scipy.sparse import csc_matrix, dok_matrix
-    from numpy import transpose, array, arange
-    
+    from scipy.sparse import csc_matrix
+    from numpy import array
+    from scipy.linsolve import spdiags, spsolve, use_solver
+
     print "Inverting a sparse linear system:"
     print "The sparse matrix (constructed from diagonals):"
     a = spdiags([[1, 2, 3, 4, 5], [6, 5, 8, 9, 10]], [0, 1], 5, 5)
     b = array([1, 2, 3, 4, 5])
     print "Solve: single precision complex:"
-    globals()['useUmfpack'] = False
+    use_solver( use = {'useUmfpack' : False} )
     a = a.astype('F')
     x = spsolve(a, b)
     print x
     print "Error: ", a*x-b
 
     print "Solve: double precision complex:"
-    globals()['useUmfpack'] = True
+    use_solver( use = {'useUmfpack' : True} )
     a = a.astype('D')
     x = spsolve(a, b)
     print x
@@ -118,7 +119,7 @@
     print "Error: ", a*x-b
 
     print "Solve: single precision:"
-    globals()['useUmfpack'] = False
+    use_solver( use = {'useUmfpack' : False} )
     a = a.astype('f')
     x = spsolve(a, b.astype('f'))
     print x

Modified: trunk/Lib/linsolve/umfpack/__init__.py
===================================================================
--- trunk/Lib/linsolve/umfpack/__init__.py	2006-08-05 20:24:43 UTC (rev 2151)
+++ trunk/Lib/linsolve/umfpack/__init__.py	2006-08-09 09:47:04 UTC (rev 2152)
@@ -1,5 +1,7 @@
+from info import __doc__
+
 from umfpack import *
 
+__all__ = filter(lambda s:not s.startswith('_'),dir())
 from numpy.testing import ScipyTest
 test = ScipyTest().test
-

Modified: trunk/Lib/linsolve/umfpack/info.py
===================================================================
--- trunk/Lib/linsolve/umfpack/info.py	2006-08-05 20:24:43 UTC (rev 2151)
+++ trunk/Lib/linsolve/umfpack/info.py	2006-08-09 09:47:04 UTC (rev 2152)
@@ -1,5 +1,119 @@
 """
-UMFPACK v4.4 wrappers.
+Interface to the UMFPACK library.
+=================================
+
+Routines for symbolic and numeric LU factorization of sparse
+matrices and for solving systems of linear equations with sparse matrices.
+
+Tested with UMFPACK V4.4 (Jan. 28, 2005), V5.0 (May 5, 2006)
+Copyright (c) 2005 by Timothy A. Davis.  All Rights Reserved.
+UMFPACK homepage: http://www.cise.ufl.edu/research/sparse/umfpack
+
+Contains: UmfpackContext class
+
+Use 'print UmfpackContext().funs' to see all UMFPACK library functions the
+module exposes, if you need something not covered by the examples below.
+
+Installation:
+=============
+Example site.cfg entry:
+
+UMFPACK v4.4 in <dir>:
+
+[amd]
+library_dirs = <dir>/UMFPACK/AMD/Lib
+include_dirs = <dir>/UMFPACK/AMD/Include
+amd_libs = amd
+
+[umfpack]
+library_dirs = <dir>/UMFPACK/UMFPACK/Lib
+include_dirs = <dir>/UMFPACK/UMFPACK/Include
+umfpack_libs = umfpack
+
+UMFPACK v5.0 (as part of UFsparse package) in <dir>:
+
+[amd]
+library_dirs = <dir>/UFsparse/AMD/Lib
+include_dirs = <dir>/UFsparse/AMD/Include, <dir>/UFsparse/UFconfig
+amd_libs = amd
+
+[umfpack]
+library_dirs = <dir>/UFsparse/UMFPACK/Lib
+include_dirs = <dir>/UFsparse/UMFPACK/Include, <dir>/UFsparse/UFconfig
+umfpack_libs = umfpack
+
+Examples:
+=========
+
+Assuming this module imported as um (import scipy.linsolve.umfpack as um)
+
+Sparse matrix in CSR or CSC format: mtx
+Righthand-side: rhs
+Solution: sol
+
+# Contruct the solver.
+umfpack = um.UmfpackContext() # Use default 'di' family of UMFPACK routines.
+
+# One-shot solution.
+sol = umfpack( um.UMFPACK_A, mtx, rhs, autoTranspose = True )
+# same as:
+sol = umfpack.linsolve( um.UMFPACK_A, mtx, rhs, autoTranspose = True )
+
+-or-
+
+# Make LU decomposition.
+umfpack.numeric( mtx )
+...
+# Use already LU-decomposed matrix.
+sol1 = umfpack( um.UMFPACK_A, mtx, rhs1, autoTranspose = True )
+sol2 = umfpack( um.UMFPACK_A, mtx, rhs2, autoTranspose = True )
+# same as:
+sol1 = umfpack.solve( um.UMFPACK_A, mtx, rhs1, autoTranspose = True )
+sol2 = umfpack.solve( um.UMFPACK_A, mtx, rhs2, autoTranspose = True )
+
+-or-
+
+# Make symbolic decomposition.
+umfpack.symbolic( mtx0 )
+# Print statistics.
+umfpack.report_symbolic()
+
+...
+
+# Make LU decomposition of mtx1 which has same structure as mtx0.
+umfpack.numeric( mtx1 )
+# Print statistics.
+umfpack.report_numeric()
+
+# Use already LU-decomposed matrix.
+sol1 = umfpack( um.UMFPACK_A, mtx1, rhs1, autoTranspose = True )
+
+...
+
+# Make LU decomposition of mtx2 which has same structure as mtx0.
+umfpack.numeric( mtx2 )
+sol2 = umfpack.solve( um.UMFPACK_A, mtx2, rhs2, autoTranspose = True )
+
+# Print all statistics.
+umfpack.report_info()
+
+Setting control parameters:
+===========================
+Assuming this module imported as um:
+
+List of control parameter names is accessible as 'um.umfControls' - their
+meaning and possible values are described in the UMFPACK documentation.
+To each name corresponds an attribute of the 'um' module, such as,
+for example 'um.UMFPACK_PRL' (controlling the verbosity of umfpack report
+functions).  These attributes are in fact indices into the control array
+- to set the corresponding control array value, just do the following:
+
+umfpack = um.UmfpackContext()
+umfpack.control[um.UMFPACK_PRL] = 4 # Let's be more verbose.
+
+--
+Author: Robert Cimrman
 """
 
+postpone_import = 1
 global_symbols = ['UmfpackContext']

Modified: trunk/Lib/linsolve/umfpack/umfpack.py
===================================================================
--- trunk/Lib/linsolve/umfpack/umfpack.py	2006-08-05 20:24:43 UTC (rev 2151)
+++ trunk/Lib/linsolve/umfpack/umfpack.py	2006-08-09 09:47:04 UTC (rev 2152)
@@ -1,3 +1,11 @@
+"""
+Interface to the UMFPACK library.
+
+--
+Author: Robert Cimrman
+"""
+
+
 #from base import Struct, pause
 import numpy as nm
 import scipy.sparse as sp
@@ -7,93 +15,6 @@
 except:
     _um = None
 
-__doc__ = """
-Interface to the UMFPACK library.
-Routines for symbolic and numeric LU factorization of sparse
-matrices and for solving systems of linear equations with sparse matrices.
-
-Tested with UMFPACK V4.4 (Jan. 28, 2005),
-Copyright (c) 2005 by Timothy A. Davis.  All Rights Reserved.
-UMFPACK homepage: http://www.cise.ufl.edu/research/sparse/umfpack
-
-Contains: UmfpackContext class
-
-Use 'print UmfpackContext().funs' to see all UMFPACK library functions the
-module exposes, if you need something not covered by the examples below.
-
-Examples:
-=========
-
-Assuming this module imported as um
-
-Sparse matrix in CSR or CSC format: mtx
-Righthand-side: rhs
-Solution: sol
-
-# Contruct the solver.
-umfpack = um.UmfpackContext() # Use default 'di' family of UMFPACK routines.
-
-# One-shot solution.
-sol = umfpack( um.UMFPACK_A, mtx, rhs, autoTranspose = True )
-# same as:
-sol = umfpack.linsolve( um.UMFPACK_A, mtx, rhs, autoTranspose = True )
-
--or-
-
-# Make LU decomposition.
-umfpack.numeric( mtx )
-...
-# Use already LU-decomposed matrix.
-sol1 = umfpack( um.UMFPACK_A, mtx, rhs1, autoTranspose = True )
-sol2 = umfpack( um.UMFPACK_A, mtx, rhs2, autoTranspose = True )
-# same as:
-sol1 = umfpack.solve( um.UMFPACK_A, mtx, rhs1, autoTranspose = True )
-sol2 = umfpack.solve( um.UMFPACK_A, mtx, rhs2, autoTranspose = True )
-
--or-
-
-# Make symbolic decomposition.
-umfpack.symbolic( mtx0 )
-# Print statistics.
-umfpack.report_symbolic()
-
-...
-
-# Make LU decomposition of mtx1 which has same structure as mtx0.
-umfpack.numeric( mtx1 )
-# Print statistics.
-umfpack.report_numeric()
-
-# Use already LU-decomposed matrix.
-sol1 = umfpack( um.UMFPACK_A, mtx1, rhs1, autoTranspose = True )
-
-...
-
-# Make LU decomposition of mtx2 which has same structure as mtx0.
-umfpack.numeric( mtx2 )
-sol2 = umfpack.solve( um.UMFPACK_A, mtx2, rhs2, autoTranspose = True )
-
-# Print all statistics.
-umfpack.report_info()
-
-Setting control parameters:
-===========================
-Assuming this module imported as um:
-
-List of control parameter names is accessible as 'um.umfControls' - their
-meaning and possible values are described in the UMFPACK documentation.
-To each name corresponds an attribute of the 'um' module, such as,
-for example 'um.UMFPACK_PRL' (controlling the verbosity of umfpack report
-functions).  These attributes are in fact indices into the control array
-- to set the corresponding control array value, just do the following:
-
-umfpack = um.UmfpackContext()
-umfpack.control[um.UMFPACK_PRL] = 4 # Let's be more verbose.
-
---
-Author: Robert Cimrman
-"""
-
 ##
 # 30.11.2005, c
 def updateDictWithVars( adict, module, pattern, group = None ):



More information about the Scipy-svn mailing list