jonathf / numpoly Goto Github PK

Numpy compatible polynomial representation

Home Page: https://numpoly.readthedocs.io

License: BSD 2-Clause "Simplified" License

Python 100.00%

numpoly's Introduction

Numpoly is a generic library for creating, manipulating and evaluating arrays of polynomials based on numpy.ndarray objects.

Intuitive interface for users experienced with numpy, as the library provides a high level of compatibility with the numpy.ndarray, including fancy indexing, broadcasting, numpy.dtype, vectorized operations to name a few.
Computationally fast evaluations of lots of functionality inherent from numpy.
Vectorized polynomial evaluation.
Support for arbitrary number of dimensions.
Native support for lots of numpy.<name> functions using numpy's compatibility layer (which also exists as numpoly.<name> equivalents).
Support for polynomial division through the operators /, % and divmod.
Extra polynomial specific attributes exposed on the polynomial objects like poly.exponents, poly.coefficients, poly.indeterminants etc.
Polynomial derivation through functions like numpoly.derivative, numpoly.gradient, numpoly.hessian etc.
Decompose polynomial sums into vector of addends using numpoly.decompose.
Variable substitution through numpoly.call.

Installation

Installation should be straight forward:

pip install numpoly

Example Usage

Constructing polynomial is typically done using one of the available constructors:

>>> import numpoly
>>> numpoly.monomial(start=0, stop=3, dimensions=2)
polynomial([1, q0, q0**2, q1, q0*q1, q1**2])

It is also possible to construct your own from symbols together with numpy:

>>> import numpy
>>> q0, q1 = numpoly.variable(2)
>>> numpoly.polynomial([1, q0**2-1, q0*q1, q1**2-1])
polynomial([1, q0**2-1, q0*q1, q1**2-1])

Or in combination with numpy objects using various arithmetics:

>>> q0**numpy.arange(4)-q1**numpy.arange(3, -1, -1)
polynomial([-q1**3+1, -q1**2+q0, q0**2-q1, q0**3-1])

The constructed polynomials can be evaluated as needed:

>>> poly = 3*q0+2*q1+1
>>> poly(q0=q1, q1=[1, 2, 3])
polynomial([3*q1+3, 3*q1+5, 3*q1+7])

Or manipulated using various numpy functions:

>>> numpy.reshape(q0**numpy.arange(4), (2, 2))
polynomial([[1, q0],
            [q0**2, q0**3]])
>>> numpy.sum(numpoly.monomial(13)[::3])
polynomial(q0**12+q0**9+q0**6+q0**3+1)

Installation

Installation should be straight forward from pip:

pip install numpoly

Alternatively, to get the most current experimental version, the code can be installed from Github as follows:

First time around, download the repository:

git clone [email protected]:jonathf/numpoly.git

Every time, move into the repository:
```
cd numpoly/
```
After the first time, you want to update the branch to the most current version of master:
```
git checkout master
git pull
```
Install the latest version of numpoly with:
```
pip install .
```

Development

Installing numpoly for development can be done from the repository root with the command:

pip install -e .[dev]

The deployment of the code is done with Python 3.10 and dependencies are then fixed using:

pip install -r requirements-dev.txt

Testing

To run test:

pytest --doctest-modules numpoly test docs/user_guide/*.rst README.rst

Documentation

To build documentation locally on your system, use make from the doc/ folder:

cd doc/
make html

Run make without argument to get a list of build targets. All targets stores output to the folder doc/.build/html.

Note that the documentation build assumes that pandoc is installed on your system and available in your path.

numpoly's People

Contributors

Stargazers

Watchers

Forkers

fredrikmeyer object-oriented-python utsekaj42 beroda vincentvanlaer

numpoly's Issues

AttributeError with division

Hi,

I'm a brand new user on Github and therefore sorry if it is not how I'm expected to use this topic.

Your package seems very useful and I thank you for implementing all of that. However, I'm finding an error when dividing a polynomial by another. Using "/", "%" and "divmod" all end up in an error :
AttributeError: module 'numpy' has no attribute 'broadcast_shapes'

Do you have any idea what in the division functions could end up in such an error ? I tried with elements of the 'numpoly.baseclass.ndpoly' class of different lengths, and none worked. Note that the same error appears when using the symbol "*".

Thank you,
Archer

Function for gathering all lead coefficients

Get the coefficient in from of the polynomial term with the highest power.

Need to support ordering to deal with multivariate polynomial.

Element-in support: ``x in y`` and ``{numpy,numpoly}.isin(x, y)``

The definition of what those two features means needs defining before implementing.

Min/max/argmin/argmax support

For polynomials, first sort lexicographically by largest exponent: x**3 > x*y > y
For equal order, sort by largest lead coefficient: 4*x**2 > 2*x*y > -8*y**2

This alow for being backward compatible with numpy sorting methods as constant is a lead coefficient.

Addition of array of polynomials with scalar raises ValueError when comparing exponents.

Hello jonathf!

I have been using your modules in my research for a while now and I love them. I have stumbled upon a problem I can't manage to solve, however.

I have an array of polynomials (with 2 elements in the example) and want to add 1 to all of them. Running the following code:

import numpoly
qs = numpoly.variable(2)
qs+1

Raises the following error:

Traceback (most recent call last):
  File "/home/quest/juvdnbro/repos/gpc/src/schrodinger1d/test.py", line 3, in <module>
    qs+1
  File "/home/quest/juvdnbro/.local/lib/python3.10/site-packages/numpoly/baseclass.py", line 242, in __array_ufunc__
    return numpoly.UFUNC_COLLECTION[ufunc](*inputs, **kwargs)
  File "/home/quest/juvdnbro/.local/lib/python3.10/site-packages/numpoly/array_function/add.py", line 60, in add
    return simple_dispatch(
  File "/home/quest/juvdnbro/.local/lib/python3.10/site-packages/numpoly/dispatch.py", line 84, in simple_dispatch
    inputs = numpoly.align_polynomials(*inputs)
  File "/home/quest/juvdnbro/.local/lib/python3.10/site-packages/numpoly/align.py", line 47, in align_polynomials
    polys = align_exponents(*polys)
  File "/home/quest/juvdnbro/.local/lib/python3.10/site-packages/numpoly/align.py", line 211, in align_exponents
    if numpy.all(poly.exponents == global_exponents):
ValueError: operands could not be broadcast together with shapes (2,2) (3,2)

It seems that the scalar 1 gets converted into polynomial([1,1]), but the exponents don't have the correct shape. Is there a better way to do this, or would it be possible to allow for this kind of arithmetic?

Thanks in advance!
Jul

Deprecation pkg_resources

This is a duplicate of jonathf/chaospy#408 for this package.
I can also propose a PR for this one if it is fine for you.

Infinity loop in poly_divmod when high relative error is introduced

Issue in Chaospy: jonathf/chaospy#257

six is not in install_requires but numpoly/baseclass.py imports six

Updating the PyPi package

PR #111 removed a bug with the newest version of numpy, but unfortunately has not yet found its way into the official repo.

Would you be willing to update the PyPi package?

Thanks!

Problem with flatiter

I believe the flatiter might not be handled correctly.

import numpoly
X = numpoly.variable(6).reshape((3,2))
print('using flatten:')
for p in X.flatten():
    print(p)
print('using flat:')
for p in X.flat:
    print(p)

output:

using flatten:
q0
q1
q2
q3
q4
q5
using flat:
(1, 0, 0, 0, 0, 0)
(0, 1, 0, 0, 0, 0)
(0, 0, 1, 0, 0, 0)
(0, 0, 0, 1, 0, 0)
(0, 0, 0, 0, 1, 0)
(0, 0, 0, 0, 0, 1)

``prod`` implementation

Current implementation uses multiply repeatedly along axis. This can likely be
done much more efficiently with dedicated code.

Somewhat the same problem as issue #46, but a lot more book keeping.

test failed

test/test_align.py . [ 0%]
test/test_array_function.py ............................................ [ 35%]
....................................................................... [ 90%]
test/test_baseclass.py F [ 91%]
test/test_construct.py .. [ 92%]
test/test_poly_function.py ......... [100%]
=================================== FAILURES ===================================
_________________________________ test_scalars _________________________________
def test_scalars():
assert XY.shape == (2,)
assert XY.size == 2
assert X.shape == ()
assert X.size == 1
assert EMPTY.shape in [(), (0,)] # different behavior in py2/3
assert EMPTY.size == 0

    assert numpy.all(numpy.array(XY.coefficients) == [[1, 0], [0, 1]])
    assert X.coefficients == [1]
    assert EMPTY.coefficients == []

    assert numpy.all(XY.exponents == [[1, 0], [0, 1]])
    assert XY.exponents.shape == (2, 2)
    assert X.exponents == 1
    assert X.exponents.shape == (1, 1)
    assert numpy.all(EMPTY.exponents == 0)
    assert EMPTY.exponents.shape == (1, 1)

    assert numpy.all(XY.indeterminants == XY)
    assert X.indeterminants == X

    assert numpy.all(XY.values == numpy.array(
        [(1, 0), (0, 1)], dtype=[("<;", "<i8",), (";<", "<i8",)]))
    assert X.values == numpy.array((1,), dtype=[("<", "<i8",)])

  assert EMPTY.values.dtype == numpy.dtype([(";", "<i8")])

E AssertionError: assert dtype([(';', '<i4')]) == dtype([(';', '<i8')])
E + where dtype([(';', '<i4')]) = array([], dtype=[(';', '<i4')]).dtype
E + where array([], dtype=[(';', '<i4')]) = EMPTY.values
E + and dtype([(';', '<i8')]) = <class 'numpy.dtype'>([(';', '<i8')])
E + where <class 'numpy.dtype'> = numpy.dtype
test/test_baseclass.py:33: AssertionError
=============================== warnings summary ===============================
test/test_array_function.py::test_numpy_isfinite[func_interface0]
test/test_array_function.py::test_numpy_isfinite[func_interface1]
/builddir/build/BUILD/numpoly-0.1.16/test/test_array_function.py:236: RuntimeWarning: invalid value encountered in log
poly = numpoly.polynomial([numpy.log(-1.), X, numpy.log(0)])
test/test_array_function.py::test_numpy_isfinite[func_interface0]
test/test_array_function.py::test_numpy_isfinite[func_interface1]
/builddir/build/BUILD/numpoly-0.1.16/test/test_array_function.py:236: RuntimeWarning: divide by zero encountered in log
poly = numpoly.polynomial([numpy.log(-1.), X, numpy.log(0)])
test/test_array_function.py::test_numpy_nonzero[interface0]
test/test_array_function.py::test_numpy_nonzero[interface1]
/usr/lib/python3.8/site-packages/numpy/core/fromnumeric.py:61: DeprecationWarning: Calling nonzero on 0d arrays is deprecated, as it behaves surprisingly. Use atleast_1d(cond).nonzero() if the old behavior was intended. If the context of this warning is of the form arr[nonzero(cond)], just use arr[cond].
return bound(*args, **kwds)
test/test_array_function.py::test_numpy_nonzero[interface2]
/builddir/build/BUILD/numpoly-0.1.16/conftest.py:17: DeprecationWarning: Calling nonzero on 0d arrays is deprecated, as it behaves surprisingly. Use atleast_1d(arr).nonzero() if the old behavior was intended.
return func(*args, **kwargs)
test/test_poly_function.py::test_numpoly_cross_truncate
/builddir/build/BUILD/numpoly-0.1.16/numpoly/poly_function/monomial/cross_truncation.py:49: RuntimeWarning: divide by zero encountered in true_divide
out = numpy.sum((indices/bound)norm, axis=-1)(1./norm) <= 1
test/test_poly_function.py::test_numpoly_cross_truncate
/builddir/build/BUILD/numpoly-0.1.16/numpoly/poly_function/monomial/cross_truncation.py:49: RuntimeWarning: invalid value encountered in true_divide
out = numpy.sum((indices/bound)norm, axis=-1)(1./norm) <= 1
test/test_poly_function.py::test_numpoly_cross_truncate
/builddir/build/BUILD/numpoly-0.1.16/numpoly/poly_function/monomial/cross_truncation.py:49: RuntimeWarning: invalid value encountered in less_equal
out = numpy.sum((indices/bound)norm, axis=-1)(1./norm) <= 1
test/test_poly_function.py::test_numpoly_cross_truncate
/builddir/build/BUILD/numpoly-0.1.16/numpoly/poly_function/monomial/cross_truncation.py:47: RuntimeWarning: divide by zero encountered in true_divide
out = numpy.max(indices/bound, axis=-1) <= 1
test/test_poly_function.py::test_numpoly_cross_truncate
/builddir/build/BUILD/numpoly-0.1.16/numpoly/poly_function/monomial/cross_truncation.py:47: RuntimeWarning: invalid value encountered in true_divide
out = numpy.max(indices/bound, axis=-1) <= 1
test/test_poly_function.py::test_numpoly_cross_truncate
/builddir/build/BUILD/numpoly-0.1.16/numpoly/poly_function/monomial/cross_truncation.py:47: RuntimeWarning: invalid value encountered in less_equal
out = numpy.max(indices/bound, axis=-1) <= 1
-- Docs: https://docs.pytest.org/en/latest/warnings.html
============== 1 failed, 127 passed, 13 warnings in 15.83 seconds ==============

Add support save/load to disk

Pickle should work out of the box. Preferably using the numpy.save
interface.

Please update the version on the PyPI website

here

Diagonalize parametric matrix

Hello,

is it possible to get a parametric eigendecomposition of a matrix whose elements are a mix of known numerical values and parametric variables, expressed as ndpoly entries?

Thank you.

``power`` implementation

Current implementation is a little hacky: repeat multiplying against itself n
times. For-loop over n, if array.

Create a method that does not rely on repeated calling multiply, and instead
allocates a single chunk of memory and fill inn results there.
Somehow avoid using element-by-element for-loop over exponents.

Make proxy array function

The features suggested in #54 depend on sorting the polynomial. Instead of sorting in place, create a function that substitute polynomial for integers that whose order correspond to the sorting in question.

Easies way to solve #54 for now.

Problem with cross_truncate

I found a problem with cross_truncation. A concrete case is:
cross_truncate(numpy.array([[5,1,1,1,1]]), 10, 1.0)
We are expecting the answer to be true (the sum is 9 which is below 9), but the answer is False. I think this is due to numerical errors, and maybe the definition of the norm should be reformulated (to avoid divisions)?

``cumprod`` implementation

Same as #47, but quite a bit harder as even more bookeeping is needed. Might not even be possible without a element-by-element for-loop.

Preserving indeterminants with numpoly.multiply etc

Is there a way to preserve the indeterminants of a set of polynomials when one applies multiply, sum, etc? For example, I create a set of 3-dimensional polynomials which may contain indeterminants q0, q1, q2 but sometimes only q0 (like when you create a basis of graded lexicographic monomials in d-dimensions). After applying align_indeterminants() I now have a set where all polynomials return ('q0, 'q1', 'q2') as "names" even though a polynomial may only be 'q0^2' etc. However, when I use numpoly.multiply(q0^2, q1) I get a new polynomial that has as names ('q0', 'q1') only. I would like to preserve the structure I had before where the new polynomial still returns ('q0, 'q1', 'q2') as "names" even though indeterminate 'q2' does not appear in the polynomial. This is useful because I have a dictionary that has as keys the exponents in the form [., ., .].

function for gathering all lead exponents

Same as #53, but returns the actual exponent instead of the coefficient.

jonathf / numpoly Goto Github PK

numpoly's Introduction

Installation

Example Usage

Installation

Development

Testing

Documentation

numpoly's People

Contributors

Stargazers

Watchers

Forkers

numpoly's Issues

Recommend Projects

Recommend Topics

Recommend Org