pyneqsys package

Package for numerically solving symbolically defined systems of non-linear equations.

pyneqsys ties computer algebra systems like:

  • SymPy
  • symengine
  • SymCXX

and numerical solvers such as:

  • MINPACK in SciPy
  • NLEQ2 in pynleq2
  • KINSOL in pykinsol

together. In addition pyneqsys provides helper classes for handling e.g. conditional equations in a system.

Submodules

pyneqsys.core module

class pyneqsys.core.ChainedNeqSys(neqsystems, **kwargs)[source]

Bases: pyneqsys.core._NeqSysBase

Chain multiple formulations of non-linear systems for using the result of one as starting guess for the other

Examples

>>> neqsys_lin = NeqSys(1, 1, lambda x, p: [x[0]**2 - p[0]])
>>> from math import log, exp
>>> neqsys_log = NeqSys(1, 1, lambda x, p: [2*x[0] - log(p[0])],
...    pre_processors=[lambda x, p: ([log(x[0]+1e-60)], p)],
...    post_processors=[lambda x, p: ([exp(x[0])], p)])
>>> chained = ChainedNeqSys([neqsys_log, neqsys_lin])
>>> x, info = chained.solve([1, 1], [4])
>>> assert info['success']
>>> print(x)
[ 2.]
>>> print(info['intermediate_info'][0]['nfev'],
...       info['intermediate_info'][1]['nfev'])  
4 3

Methods

plot_series(xres, varied_data, varied_idx, …) Plots the results from solve_series().
plot_series_residuals(xres, varied_data, …) Analogous to plot_series() but will plot residuals.
plot_series_residuals_internal(varied_data, …) Analogous to plot_series() but for internal residuals from last run.
post_process(x, params) Used internally for transformation of variables.
pre_process(x, params[, conds]) Used internally for transformation of variables.
rms(x[, params]) Returns root mean square value of f(x, params)
solve(x0[, params, internal_x0, solver]) Solve with user specified solver choice.
solve_and_plot_series(x0, params, …[, …]) Solve and plot for a series of a varied parameter.
solve_series(x0, params, varied_data, varied_idx) Solve system for a set of parameters in which one is varied
post_process(x, params)[source]

Used internally for transformation of variables.

pre_process(x, params, conds=None)[source]

Used internally for transformation of variables.

solve(x0, params=(), internal_x0=None, solver=None, **kwargs)[source]

Solve with user specified solver choice.

Parameters:

x0: 1D array of floats

Guess (subject to self.post_processors)

params: 1D array_like of floats

Parameters (subject to self.post_processors)

internal_x0: 1D array of floats

When given it overrides (processed) x0. internal_x0 is not subject to self.post_processors.

solver: str or callable or None or iterable of such

if str: uses _solve_``solver``(*args, **kwargs). if None: chooses from PYNEQSYS_SOLVER environment variable. if iterable: chain solving.

attached_solver: callable factory

Invokes: solver = attached_solver(self).

Returns:

array:

solution vector (post-processed by self.post_processors)

dict:

info dictionary containing ‘success’, ‘nfev’, ‘njev’ etc.

Examples

>>> neqsys = NeqSys(2, 2, lambda x, p: [
...     (x[0] - x[1])**p[0]/2 + x[0] - 1,
...     (x[1] - x[0])**p[0]/2 + x[1]
... ])
>>> x, sol = neqsys.solve([1, 0], [3], solver=(None, 'mpmath'))
>>> assert sol['success']
>>> print(x)
[0.841163901914009663684741869855]
[0.158836098085990336315258130144]
class pyneqsys.core.ConditionalNeqSys(condition_cb_pairs, neqsys_factory, **kwargs)[source]

Bases: pyneqsys.core._NeqSysBase

Collect multiple systems of non-linear equations with different conditionals.

If a problem in a fixed number of variables is described by different systems of equations for different domains, then this class may be used to describe that set of systems.

The user provides a set of conditions which governs what system of equations to apply. The set of conditions then represent a vector of booleans which is passed to a user provided factory function of NeqSys instances. The conditions may be asymmetrical (each condition consits of two callbacks, one for evaluating when the condition was previously False, and one when it was previously True. The motivation for this asymmetry is that a user may want to introduce a tolerance for numerical noise in the solution (and avoid possibly endless loops).

If fastcache is available an LRU cache will be used for neqsys_factory, it is therefore important that the factory function is idempotent.

Parameters:

condition_cb_pairs : list of (callback, callback) tuples

Callbacks should have the signature: f(x, p) -> bool.

neqsys_factory : callback

Should have the signature f(conds) -> NeqSys instance, where conds is a list of bools.

names : list of strings

Examples

>>> from math import sin, pi
>>> f_a = lambda x, p: [sin(p[0]*x[0])]  # when x <= 0
>>> f_b = lambda x, p: [x[0]*(p[1]-x[0])]  # when x >= 0
>>> factory = lambda conds: NeqSys(1, 1, f_b) if conds[0] else NeqSys(
...     1, 1, f_a)
>>> cneqsys = ConditionalNeqSys([(lambda x, p: x[0] > 0,  # no 0-switch
...                               lambda x, p: x[0] >= 0)],  # no 0-switch
...                             factory)
>>> x, sol = cneqsys.solve([0], [pi, 3])
>>> assert sol['success']
>>> print(x)
[ 0.]
>>> x, sol = cneqsys.solve([-1.4], [pi, 3])
>>> assert sol['success']
>>> print(x)
[-1.]
>>> x, sol = cneqsys.solve([2], [pi, 3])
>>> assert sol['success']
>>> print(x)
[ 3.]
>>> x, sol = cneqsys.solve([7], [pi, 3])
>>> assert sol['success']
>>> print(x)
[ 3.]

Methods

f_cb(x, params[, conds])
get_conds(x, params[, prev_conds])
plot_series(xres, varied_data, varied_idx, …) Plots the results from solve_series().
plot_series_residuals(xres, varied_data, …) Analogous to plot_series() but will plot residuals.
plot_series_residuals_internal(varied_data, …) Analogous to plot_series() but for internal residuals from last run.
post_process(x, params[, conds]) Used internally for transformation of variables.
pre_process(x, params[, conds]) Used internally for transformation of variables.
rms(x[, params]) Returns root mean square value of f(x, params)
solve(x0[, params, internal_x0, solver, …]) Solve the problem (systems of equations)
solve_and_plot_series(x0, params, …[, …]) Solve and plot for a series of a varied parameter.
solve_series(x0, params, varied_data, varied_idx) Solve system for a set of parameters in which one is varied
f_cb(x, params, conds=None)[source]
get_conds(x, params, prev_conds=None)[source]
post_process(x, params, conds=None)[source]

Used internally for transformation of variables.

pre_process(x, params, conds=None)[source]

Used internally for transformation of variables.

solve(x0, params=(), internal_x0=None, solver=None, conditional_maxiter=20, initial_conditions=None, **kwargs)[source]

Solve the problem (systems of equations)

Parameters:

x0 : array

Guess.

params : array

internal_x0 : array

solver : str or callable or iterable of such.

conditional_maxiter : int

Maximum number of switches between conditions.

initial_conditions : iterable of bools

Corresponding conditions to x0

**kwargs :

Keyword arguments passed on to NeqSys.solve().

class pyneqsys.core.NeqSys(nf, nx=None, f=None, jac=None, band=None, pre_processors=None, post_processors=None, internal_x0_cb=None, **kwargs)[source]

Bases: pyneqsys.core._NeqSysBase

Represents a system of non-linear equations.

This class provides a unified interface to:

  • scipy.optimize.root
  • NLEQ2
  • KINSOL
  • mpmath
  • levmar
Parameters:

nf : int

Number of functions.

nx : int

Number of independent variables.

f : callback

Function to solve for. Signature f(x) -> y where len(x) == nx and len(y) == nf.

jac : callback or None (default)

Jacobian matrix (dfdy).

band : tuple (default: None)

Number of sub- and super-diagonals in jacobian.

names : iterable of str (default: None)

Names of variables, used for plotting and for referencing by name.

param_names : iterable of strings (default: None)

Names of the parameters, used for referencing parameters by name.

x_by_name : bool, default: False

Will values for x be referred to by name (in dictionaries) instead of by index (in arrays)?

par_by_name : bool, default: False

Will values for parameters be referred to by name (in dictionaries) instead of by index (in arrays)?

latex_names : iterable of str, optional

Names of variables in LaTeX format.

latex_param_names : iterable of str, optional

Names of parameters in LaTeX format.

pre_processors : iterable of callables (optional)

(Forward) transformation of user-input to solve() signature: f(x1[:], params1[:]) -> x2[:], params2[:]. Insert at beginning.

post_processors : iterable of callables (optional)

(Backward) transformation of result from solve() signature: f(x2[:], params2[:]) -> x1[:], params1[:]. Insert at end.

internal_x0_cb : callback (optional)

callback with signature f(x[:], p[:]) -> x0[:] if not specified, x from self.pre_processors will be used.

See also

pyneqsys.symbolic.SymbolicSys
use a CAS (SymPy by default) to derive the jacobian.

Examples

>>> neqsys = NeqSys(2, 2, lambda x, p: [(x[0] - x[1])**p[0]/2 + x[0] - 1,
...                                     (x[1] - x[0])**p[0]/2 + x[1]])
>>> x, sol = neqsys.solve([1, 0], [3])
>>> assert sol['success']
>>> print(x)
[ 0.8411639  0.1588361]

Methods

plot_series(xres, varied_data, varied_idx, …) Plots the results from solve_series().
plot_series_residuals(xres, varied_data, …) Analogous to plot_series() but will plot residuals.
plot_series_residuals_internal(varied_data, …) Analogous to plot_series() but for internal residuals from last run.
post_process(xout, params_out) Used internally for transformation of variables.
pre_process(x0[, params]) Used internally for transformation of variables.
rms(x[, params]) Returns root mean square value of f(x, params)
solve(x0[, params, internal_x0, solver, …]) Solve with user specified solver choice.
solve_and_plot_series(x0, params, …[, …]) Solve and plot for a series of a varied parameter.
solve_series(x0, params, varied_data, varied_idx) Solve system for a set of parameters in which one is varied
post_process(xout, params_out)[source]

Used internally for transformation of variables.

pre_process(x0, params=())[source]

Used internally for transformation of variables.

solve(x0, params=(), internal_x0=None, solver=None, attached_solver=None, **kwargs)[source]

Solve with user specified solver choice.

Parameters:

x0: 1D array of floats

Guess (subject to self.post_processors)

params: 1D array_like of floats

Parameters (subject to self.post_processors)

internal_x0: 1D array of floats

When given it overrides (processed) x0. internal_x0 is not subject to self.post_processors.

solver: str or callable or None or iterable of such

if str: uses _solve_``solver``(*args, **kwargs). if None: chooses from PYNEQSYS_SOLVER environment variable. if iterable: chain solving.

attached_solver: callable factory

Invokes: solver = attached_solver(self).

Returns:

array:

solution vector (post-processed by self.post_processors)

dict:

info dictionary containing ‘success’, ‘nfev’, ‘njev’ etc.

Examples

>>> neqsys = NeqSys(2, 2, lambda x, p: [
...     (x[0] - x[1])**p[0]/2 + x[0] - 1,
...     (x[1] - x[0])**p[0]/2 + x[1]
... ])
>>> x, sol = neqsys.solve([1, 0], [3], solver=(None, 'mpmath'))
>>> assert sol['success']
>>> print(x)
[0.841163901914009663684741869855]
[0.158836098085990336315258130144]

pyneqsys.plotting module

pyneqsys.plotting.mpl_outside_legend(ax, **kwargs)[source]

Places a legend box outside a matplotlib Axes instance.

pyneqsys.plotting.plot_series(xres, varied_data, indices=None, info=None, fail_vline=None, plot_kwargs_cb=None, ls=('-', '--', ':', '-.'), c=('k', 'r', 'g', 'b', 'c', 'm', 'y'), labels=None, ax=None, names=None, latex_names=None)[source]

Plot the values of the solution vector vs the varied parameter.

Parameters:

xres : array

Solution vector of shape (varied_data.size, x0.size).

varied_data : array

Numerical values of the varied parameter.

indices : iterable of integers, optional

Indices of variables to be plotted. default: all

fail_vline : bool

Show vertical lines where the solver failed.

plot_kwargs_cb : callable

Takes the index as single argument, returns a dict passed to the plotting function

ls : iterable of str

Linestyles.

c : iterable of str

Colors.

labels : iterable of str

ax : matplotlib Axes instance

names : iterable of str

latex_names : iterable of str

pyneqsys.solvers module

Currently this module contains a few solvers for demonstration purposes and it is not in the scope of pyneqsys to provide “production” solvers, but instead be a package for using well established solvers for solving systems of equations defined in a uniform (optionally symbolic) way.

Nevertheless, the solvers provided here may serve as a starting point for writing new types of solvers for systems of non-linear equations (with an emphasis on a concise API rather than optimal performance).

Do not rely on any of the classes and functions in this module (pyneqsys.solvers) since they are all to be regarded as provisional, i.e. they may be renamed, change behavior and/or signature without any prior notice.

class pyneqsys.solvers.AutoDampedGradientDescentSolver(tgt_oscill=0.1, start_damp=0.1, nhistory=4, tgt_pow=0.3)[source]

Bases: pyneqsys.solvers.GradientDescentSolver

Methods

__call__(inst) Call self as a function.
alloc()
cb_factory()
damping(iter_idx, mx_iter)
f(x)
j(x)
line_search(x, dx[, mxiter, alpha]) Goldstein-Armijo linesearch (backtracking)
step(x, iter_idx, maxiter)
damping(iter_idx, mx_iter)[source]
class pyneqsys.solvers.DampedGradientDescentSolver(base_damp=0.5, exp_damp=0.5)[source]

Bases: pyneqsys.solvers.GradientDescentSolver

Methods

__call__(inst) Call self as a function.
alloc()
cb_factory()
damping(iter_idx, mx_iter)
f(x)
j(x)
line_search(x, dx[, mxiter, alpha]) Goldstein-Armijo linesearch (backtracking)
step(x, iter_idx, maxiter)
damping(iter_idx, mx_iter)[source]
class pyneqsys.solvers.GradientDescentSolver[source]

Bases: pyneqsys.solvers.SolverBase

Example of a custom solver

Parameters:instance: NeqSys instance (passed by NeqSys.solve)

Notes

Note that this is an inefficient solver for demonstration purposes.

Attributes

damping: callback with signature f(iter_idx, iter_max) -> float default: lambda idx, mx: exp(-idx/mx)/2

Methods

__call__(inst) Call self as a function.
alloc()
cb_factory()
damping(iter_idx, mx_iter)
f(x)
j(x)
line_search(x, dx[, mxiter, alpha]) Goldstein-Armijo linesearch (backtracking)
step(x, iter_idx, maxiter)
damping(iter_idx, mx_iter)[source]
step(x, iter_idx, maxiter)[source]
class pyneqsys.solvers.LineSearchingGradientDescentSolver[source]

Bases: pyneqsys.solvers.SolverBase

Methods

__call__(inst) Call self as a function.
alloc()
cb_factory()
f(x)
j(x)
line_search(x, dx[, mxiter, alpha]) Goldstein-Armijo linesearch (backtracking)
step(x, iter_idx, maxiter)
step(x, iter_idx, maxiter)[source]
class pyneqsys.solvers.PolakRibiereConjugateGradientSolver(reset_freq=10)[source]

Bases: pyneqsys.solvers.SolverBase

Methods

__call__(inst) Call self as a function.
alloc()
cb_factory()
f(x)
j(x)
line_search(x, dx[, mxiter, alpha]) Goldstein-Armijo linesearch (backtracking)
step(x, iter_idx, maxiter)
alloc()[source]
step(x, iter_idx, maxiter)[source]
class pyneqsys.solvers.SolverBase[source]

Bases: object

Methods

__call__(inst) Call self as a function.
alloc()
cb_factory()
f(x)
j(x)
line_search(x, dx[, mxiter, alpha]) Goldstein-Armijo linesearch (backtracking)
step(x)
alloc()[source]
cb_factory()[source]
f(x)[source]
j(x)[source]

Goldstein-Armijo linesearch (backtracking)

step(x)[source]
pyneqsys.solvers.rms(x)[source]

pyneqsys.symbolic module

class pyneqsys.symbolic.SymbolicSys(x, exprs, params=(), jac=True, backend=None, **kwargs)[source]

Bases: pyneqsys.core.NeqSys

Symbolically defined system of non-linear equations.

This object is analogous to pyneqsys.NeqSys but instead of providing a callable, the user provides symbolic expressions.

Parameters:

x : iterable of Symbols

exprs : iterable of expressions for f

params : iterable of Symbols (optional)

list of symbols appearing in exprs which are parameters

jac : ImmutableMatrix or bool

If True:
  • Calculate Jacobian from exprs.
If False:
  • Do not compute Jacobian (numeric approximation).
If ImmutableMatrix:
  • User provided expressions for the Jacobian.

backend : str or sym.Backend

See documentation of sym.Backend.

module : str

module keyword argument passed to backend.Lambdify.

**kwargs:

Notes

When using SymPy as the backend, a limited number of unknowns is supported. The reason is that (currently) sympy.lambdify has an upper limit on number of arguments.

Examples

>>> import sympy as sp
>>> e = sp.exp
>>> x = x0, x1 = sp.symbols('x:2')
>>> params = a, b = sp.symbols('a b')
>>> neqsys = SymbolicSys(x, [a*(1 - x0), b*(x1 - x0**2)], params)
>>> xout, sol = neqsys.solve('scipy', [-10, -5], [1, 10])
>>> print(xout)
[ 1.  1.]
>>> print(neqsys.get_jac()[0, 0])
-a

Methods

from_callback(cb[, nx, nparams]) Generate a SymbolicSys instance from a callback.
get_jac() Return the jacobian of the expressions
plot_series(xres, varied_data, varied_idx, …) Plots the results from solve_series().
plot_series_residuals(xres, varied_data, …) Analogous to plot_series() but will plot residuals.
plot_series_residuals_internal(varied_data, …) Analogous to plot_series() but for internal residuals from last run.
post_process(xout, params_out) Used internally for transformation of variables.
pre_process(x0[, params]) Used internally for transformation of variables.
rms(x[, params]) Returns root mean square value of f(x, params)
solve(x0[, params, internal_x0, solver, …]) Solve with user specified solver choice.
solve_and_plot_series(x0, params, …[, …]) Solve and plot for a series of a varied parameter.
solve_series(x0, params, varied_data, varied_idx) Solve system for a set of parameters in which one is varied
classmethod from_callback(cb, nx=None, nparams=None, **kwargs)[source]

Generate a SymbolicSys instance from a callback.

Parameters:

cb : callable

Should have the signature cb(x, p, backend) -> list of exprs.

nx : int

Number of unknowns, when not given it is deduced from kwargs['names'].

nparams : int

Number of parameters, when not given it is deduced from kwargs['param_names'].

**kwargs :

Keyword arguments passed on to SymbolicSys. See also pyneqsys.NeqSys.

Examples

>>> symbolicsys = SymbolicSys.from_callback(lambda x, p, be: [
...     x[0]*x[1] - p[0],
...     be.exp(-x[0]) + be.exp(-x[1]) - p[0]**-2
... ], 2, 1)
...
get_jac()[source]

Return the jacobian of the expressions

class pyneqsys.symbolic.TransformedSys(x, exprs, transf, params=(), post_adj=None, **kwargs)[source]

Bases: pyneqsys.symbolic.SymbolicSys

A system which transforms the equations and variables internally

Can be used to reformulate a problem in a numerically more stable form.

Parameters:

x : iterable of variables

exprs : iterable of expressions

Expressions to find root for (untransformed).

transf : iterable of pairs of expressions

Forward, backward transformed instances of x.

params : iterable of symbols

post_adj : callable (default: None)

To tweak expression after transformation.

**kwargs :

Keyword arguments passed onto SymbolicSys.

Methods

from_callback(cb, transf_cbs, nx[, nparams, …]) Generate a TransformedSys instance from a callback
get_jac() Return the jacobian of the expressions
plot_series(xres, varied_data, varied_idx, …) Plots the results from solve_series().
plot_series_residuals(xres, varied_data, …) Analogous to plot_series() but will plot residuals.
plot_series_residuals_internal(varied_data, …) Analogous to plot_series() but for internal residuals from last run.
post_process(xout, params_out) Used internally for transformation of variables.
pre_process(x0[, params]) Used internally for transformation of variables.
rms(x[, params]) Returns root mean square value of f(x, params)
solve(x0[, params, internal_x0, solver, …]) Solve with user specified solver choice.
solve_and_plot_series(x0, params, …[, …]) Solve and plot for a series of a varied parameter.
solve_series(x0, params, varied_data, varied_idx) Solve system for a set of parameters in which one is varied
classmethod from_callback(cb, transf_cbs, nx, nparams=0, pre_adj=None, **kwargs)[source]

Generate a TransformedSys instance from a callback

Parameters:

cb : callable

Should have the signature cb(x, p, backend) -> list of exprs. The callback cb should return untransformed expressions.

transf_cbs : pair or iterable of pairs of callables

Callables for forward- and backward-transformations. Each callable should take a single parameter (expression) and return a single expression.

nx : int

Number of unkowns.

nparams : int

Number of parameters.

pre_adj : callable, optional

To tweak expression prior to transformation. Takes a sinlge argument (expression) and return a single argument rewritten expression.

**kwargs :

Keyword arguments passed on to TransformedSys. See also SymbolicSys and pyneqsys.NeqSys.

Examples

>>> import sympy as sp
>>> transformed = TransformedSys.from_callback(lambda x, p, be: [
...     x[0]*x[1] - p[0],
...     be.exp(-x[0]) + be.exp(-x[1]) - p[0]**-2
... ], (sp.log, sp.exp), 2, 1)
...
pyneqsys.symbolic.linear_exprs(A, x, b=None, rref=False, Matrix=None)[source]

Returns Ax - b

Parameters:

A : matrix_like of numbers

Of shape (len(b), len(x)).

x : iterable of symbols

b : array_like of numbers (default: None)

When None, assume zeros of length len(x).

Matrix : class

When rref == True: A matrix class which supports slicing, and methods __mul__ and rref. Defaults to sympy.Matrix.

rref : bool

Calculate the reduced row echelon form of (A | -b).

Returns:

A list of the elements in the resulting column vector.

pyneqsys.symbolic.linear_rref(A, b, Matrix=None, S=None)[source]

Transform a linear system to reduced row-echelon form

Transforms both the matrix and right-hand side of a linear system of equations to reduced row echelon form

Parameters:

A : Matrix-like

Iterable of rows.

b : iterable

Returns:

A’, b’ - transformed versions