SolversProposal - SciPy wiki dump

Proposal of New Solver Classes for SciPy

TODO

agree upon the arguments and their names, see below
implement the linear solvers first, as they are the most used, as a proof of concept
eigenvalue solvers
nonlinear solvers
...

Argument names

Conventions: lowercase letters with words connected using underscores, as it is usual in SciPy.

Examples:

argument names
mtx	matrix (LinearOperator instance)
mtx_a, mtx_b	matrices (LinearOperator instances)
precond	preconditioner (LinearOperator instance)
conf	solver configuration (holding e.g. tol_abs, tol_rel, ... in a single dictionary
tol_abs	absolute tolerance
tol_rel	relatice tolerance
max_iter	max. number of iterations
callback	callback
dtype	result dtype
stats	output argument (a dictionary) holding various statistics a solver might provide (timing, convergence...)

Examples of Usage

"Script" usage

Suppose we have:

   1 class Umfpack( LinearSolver ):
   2     name = 'ls.umfpack'
   3     ...
   4 
   5 class ScipyIterative( LinearSolver ):
   6     name = 'ls.scipy_iterative'
   7     ...

Then:

   1 s = ScipyIterative( tol_abs = 1e-6, mtx = A, method = 'cg' )
   2 x = s( rhs )
   3 x2 = s( rhs2, abs_tol = 1e-8 )
   4 
   5 conf = {'tol_abs' : 1e-5, 'tol_rel' : 1e-5, 'method' : 'bicg'}
   6 s = ScipyIterative( conf = conf )
   7 x = s( rhs, conf = conf, mtx = B )
   8 
   9 ds = Umfpack( mtx = A )
  10 x = ds( rhs )

"Framework" Usage

For example, any linear solver (e.g. umfpack, scipy iterative solvers, ...) can be used like:

   1 # Use a SciPy iterative solver configuration.
   2 iconf = {
   3     'name' : 'ls',
   4     'kind' : 'ls.scipy_iterative',
   5 
   6     'method'   : 'cg',
   7     'max_iter' : 1000,
   8     'tol_abs'  : 1e-12,
   9 }
  10 # Or a direct solver configuration.
  11 dconf = {
  12     'name' : 'ls',
  13     'kind' : 'ls.umfpack',
  14 }
  15 
  16 s = Solver.any_from_conf( conf = dconf, mtx = A ) # Possibly pre-solves by LU.
  17 x1 = s( rhs1 )
  18 x2 = s( rhs2 )
  19 ...
  20 x = s( rhs, mtx = B )

The proper solver class is chosen automatically according to the 'kind' attribute of the configuration.

Abstract classes

Something to start with (ignore the argument names for now...):

   1 class Solver( Struct ):
   2      def __init__( self, conf, **kwargs ):
   3          Struct.__init__( self, conf = conf, **kwargs )
   4      def __call__( self, **kwargs ):
   5          print 'called an abstract Solver instance!'
   6          raise ValueError
   7 
   8 class LinearSolver( Solver ):
   9      def __init__( self, conf, mtx = None, status = None, **kwargs ):
  10          Solver.__init__( self, conf = conf, mtx = mtx, status = status,
  11                           **kwargs )
  12      def __call__( self, rhs, conf = None, mtx = None, status = None ):
  13          print 'called an abstract LinearSolver instance!'
  14          raise ValueError
  15 
  16 class NonlinearSolver( Solver ):
  17      def __init__( self, conf, evaluator = None, linSolver = None,
  18                    status = None, **kwargs ):
  19          Solver.__init__( self, conf = conf, evaluator = evaluator,
  20                           linSolver = linSolver, status = status,
  21                           **kwargs )
  22      def __call__( self, state0, conf = None, evaluator = None,
  23                    linSolver = None, status = None ):
  24          print 'called an abstract NonlinearSolver instance!'
  25          raise ValueError
  26 
  27 class EigenvalueSolver( Solver ):
  28      def __init__( self, conf, mtxA = None, mtxB = None, nEigs = None,
  29                    eigenvectors = None, status = None ):
  30          Solver.__init__( self, conf = conf, mtxA = mtxA, mtxB = mtxB,
  31                           nEigs = nEigs, eigenvectors = eigenvectors,
  32                           status = status )
  33      def __call__( self, mtxA, mtxB = None, nEigs = None,
  34                    eigenvectors = None, status = None, conf = None ):
  35          print 'called an abstract EigenvalueSolver instance!'
  36          raise ValueError

ProposedEnhancements