/usr/share/pyshared/numpy/polynomial/polytemplate.py

"""
Template for the Chebyshev and Polynomial classes.

This module houses a Python string module Template object (see, e.g.,
http://docs.python.org/library/string.html#template-strings) used by
the `polynomial` and `chebyshev` modules to implement their respective
`Polynomial` and `Chebyshev` classes.  It provides a mechanism for easily
creating additional specific polynomial classes (e.g., Legendre, Jacobi,
etc.) in the future, such that all these classes will have a common API.

"""
import string
import sys

if sys.version_info[0] >= 3:
    rel_import = "from . import"
else:
    rel_import = "import"

polytemplate = string.Template('''
from __future__ import division
REL_IMPORT polyutils as pu
import numpy as np

class $name(pu.PolyBase) :
    """A $name series class.

    $name instances provide the standard Python numerical methods '+',
    '-', '*', '//', '%', 'divmod', '**', and '()' as well as the listed
    methods.

    Parameters
    ----------
    coef : array_like
        $name coefficients, in increasing order.  For example,
        ``(1, 2, 3)`` implies ``P_0 + 2P_1 + 3P_2`` where the
        ``P_i`` are a graded polynomial basis.
    domain : (2,) array_like, optional
        Domain to use. The interval ``[domain[0], domain[1]]`` is mapped to
        the interval ``[window[0], window[1]]`` by shifting and scaling.
        The default value is $domain.
    window : (2,) array_like, optional
        Window, see ``domain`` for its use. The default value is $domain.
        .. versionadded:: 1.6.0

    Attributes
    ----------
    coef : (N,) array
        $name coefficients, from low to high.
    domain : (2,) array
        Domain that is mapped to ``window``.
    window : (2,) array
        Window that ``domain`` is mapped to.

    Class Attributes
    ----------------
    maxpower : int
        Maximum power allowed, i.e., the largest number ``n`` such that
        ``p(x)**n`` is allowed. This is to limit runaway polynomial size.
    domain : (2,) ndarray
        Default domain of the class.
    window : (2,) ndarray
        Default window of the class.

    Notes
    -----
    It is important to specify the domain in many cases, for instance in
    fitting data, because many of the important properties of the
    polynomial basis only hold in a specified interval and consequently
    the data must be mapped into that interval in order to benefit.

    Examples
    --------

    """
    # Limit runaway size. T_n^m has degree n*2^m
    maxpower = 16
    # Default domain
    domain = np.array($domain)
    # Default window
    window = np.array($domain)
    # Don't let participate in array operations. Value doesn't matter.
    __array_priority__ = 0

    def has_samecoef(self, other):
        """Check if coefficients match.

        Parameters
        ----------
        other : class instance
            The other class must have the ``coef`` attribute.

        Returns
        -------
        bool : boolean
            True if the coefficients are the same, False otherwise.

        Notes
        -----
        .. versionadded:: 1.6.0

        """
        if len(self.coef) != len(other.coef):
            return False
        elif not np.all(self.coef == other.coef):
            return False
        else:
            return True

    def has_samedomain(self, other):
        """Check if domains match.

        Parameters
        ----------
        other : class instance
            The other class must have the ``domain`` attribute.

        Returns
        -------
        bool : boolean
            True if the domains are the same, False otherwise.

        Notes
        -----
        .. versionadded:: 1.6.0

        """
        return np.all(self.domain == other.domain)

    def has_samewindow(self, other):
        """Check if windows match.

        Parameters
        ----------
        other : class instance
            The other class must have the ``window`` attribute.

        Returns
        -------
        bool : boolean
            True if the windows are the same, False otherwise.

        Notes
        -----
        .. versionadded:: 1.6.0

        """
        return np.all(self.window == other.window)

    def has_samewindow(self, other):
        """Check if windows match.

        Parameters
        ----------
        other : class instance
            The other class must have the ``window`` attribute.

        Returns
        -------
        bool : boolean
            True if the windows are the same, False otherwise.

        """
        return np.all(self.window == other.window)

    def __init__(self, coef, domain=$domain, window=$domain) :
        [coef, dom, win] = pu.as_series([coef, domain, window], trim=False)
        if len(dom) != 2 :
            raise ValueError("Domain has wrong number of elements.")
        if len(win) != 2 :
            raise ValueError("Window has wrong number of elements.")
        self.coef = coef
        self.domain = dom
        self.window = win

    def __repr__(self):
        format = "%s(%s, %s, %s)"
        coef = repr(self.coef)[6:-1]
        domain = repr(self.domain)[6:-1]
        window = repr(self.window)[6:-1]
        return format % ('$name', coef, domain, window)

    def __str__(self) :
        format = "%s(%s)"
        coef = str(self.coef)
        return format % ('$nick', coef)

    # Pickle and copy

    def __getstate__(self) :
        ret = self.__dict__.copy()
        ret['coef'] = self.coef.copy()
        ret['domain'] = self.domain.copy()
        ret['window'] = self.window.copy()
        return ret

    def __setstate__(self, dict) :
        self.__dict__ = dict

    # Call

    def __call__(self, arg) :
        off, scl = pu.mapparms(self.domain, self.window)
        arg = off + scl*arg
        return ${nick}val(arg, self.coef)

    def __iter__(self) :
        return iter(self.coef)

    def __len__(self) :
        return len(self.coef)

    # Numeric properties.

    def __neg__(self) :
        return self.__class__(-self.coef, self.domain, self.window)

    def __pos__(self) :
        return self

    def __add__(self, other) :
        """Returns sum"""
        if isinstance(other, self.__class__) :
            if self.has_samedomain(other) and self.has_samewindow(other):
                coef = ${nick}add(self.coef, other.coef)
            else :
                raise PolyDomainError()
        else :
            try :
                coef = ${nick}add(self.coef, other)
            except :
                return NotImplemented
        return self.__class__(coef, self.domain, self.window)

    def __sub__(self, other) :
        """Returns difference"""
        if isinstance(other, self.__class__) :
            if self.has_samedomain(other) and self.has_samewindow(other):
                coef = ${nick}sub(self.coef, other.coef)
            else :
                raise PolyDomainError()
        else :
            try :
                coef = ${nick}sub(self.coef, other)
            except :
                return NotImplemented
        return self.__class__(coef, self.domain, self.window)

    def __mul__(self, other) :
        """Returns product"""
        if isinstance(other, self.__class__) :
            if self.has_samedomain(other) and self.has_samewindow(other):
                coef = ${nick}mul(self.coef, other.coef)
            else :
                raise PolyDomainError()
        else :
            try :
                coef = ${nick}mul(self.coef, other)
            except :
                return NotImplemented
        return self.__class__(coef, self.domain, self.window)

    def __div__(self, other):
        # set to __floordiv__ /.
        return self.__floordiv__(other)

    def __truediv__(self, other) :
        # there is no true divide if the rhs is not a scalar, although it
        # could return the first n elements of an infinite series.
        # It is hard to see where n would come from, though.
        if isinstance(other, self.__class__) :
            if len(other.coef) == 1 :
                coef = div(self.coef, other.coef)
            else :
                return NotImplemented
        elif np.isscalar(other) :
            # this might be overly restrictive
            coef = self.coef/other
        else :
            return NotImplemented
        return self.__class__(coef, self.domain, self.window)

    def __floordiv__(self, other) :
        """Returns the quotient."""
        if isinstance(other, self.__class__) :
            if np.all(self.domain == other.domain) :
                quo, rem = ${nick}div(self.coef, other.coef)
            else :
                raise PolyDomainError()
        else :
            try :
                quo, rem = ${nick}div(self.coef, other)
            except :
                return NotImplemented
        return self.__class__(quo, self.domain, self.window)

    def __mod__(self, other) :
        """Returns the remainder."""
        if isinstance(other, self.__class__) :
            if self.has_samedomain(other) and self.has_samewindow(other):
                quo, rem = ${nick}div(self.coef, other.coef)
            else :
                raise PolyDomainError()
        else :
            try :
                quo, rem = ${nick}div(self.coef, other)
            except :
                return NotImplemented
        return self.__class__(rem, self.domain, self.window)

    def __divmod__(self, other) :
        """Returns quo, remainder"""
        if isinstance(other, self.__class__) :
            if self.has_samedomain(other) and self.has_samewindow(other):
                quo, rem = ${nick}div(self.coef, other.coef)
            else :
                raise PolyDomainError()
        else :
            try :
                quo, rem = ${nick}div(self.coef, other)
            except :
                return NotImplemented
        quo = self.__class__(quo, self.domain, self.window)
        rem = self.__class__(rem, self.domain, self.window)
        return quo, rem

    def __pow__(self, other) :
        try :
            coef = ${nick}pow(self.coef, other, maxpower = self.maxpower)
        except :
            raise
        return self.__class__(coef, self.domain, self.window)

    def __radd__(self, other) :
        try :
            coef = ${nick}add(other, self.coef)
        except :
            return NotImplemented
        return self.__class__(coef, self.domain, self.window)

    def __rsub__(self, other):
        try :
            coef = ${nick}sub(other, self.coef)
        except :
            return NotImplemented
        return self.__class__(coef, self.domain, self.window)

    def __rmul__(self, other) :
        try :
            coef = ${nick}mul(other, self.coef)
        except :
            return NotImplemented
        return self.__class__(coef, self.domain, self.window)

    def __rdiv__(self, other):
        # set to __floordiv__ /.
        return self.__rfloordiv__(other)

    def __rtruediv__(self, other) :
        # there is no true divide if the rhs is not a scalar, although it
        # could return the first n elements of an infinite series.
        # It is hard to see where n would come from, though.
        if len(self.coef) == 1 :
            try :
                quo, rem = ${nick}div(other, self.coef[0])
            except :
                return NotImplemented
        return self.__class__(quo, self.domain, self.window)

    def __rfloordiv__(self, other) :
        try :
            quo, rem = ${nick}div(other, self.coef)
        except :
            return NotImplemented
        return self.__class__(quo, self.domain, self.window)

    def __rmod__(self, other) :
        try :
            quo, rem = ${nick}div(other, self.coef)
        except :
            return NotImplemented
        return self.__class__(rem, self.domain, self.window)

    def __rdivmod__(self, other) :
        try :
            quo, rem = ${nick}div(other, self.coef)
        except :
            return NotImplemented
        quo = self.__class__(quo, self.domain, self.window)
        rem = self.__class__(rem, self.domain, self.window)
        return quo, rem

    # Enhance me
    # some augmented arithmetic operations could be added here

    def __eq__(self, other) :
        res = isinstance(other, self.__class__) \
                and self.has_samecoef(other) \
                and self.has_samedomain(other) \
                and self.has_samewindow(other)
        return res

    def __ne__(self, other) :
        return not self.__eq__(other)

    #
    # Extra methods.
    #

    def copy(self) :
        """Return a copy.

        Return a copy of the current $name instance.

        Returns
        -------
        new_instance : $name
            Copy of current instance.

        """
        return self.__class__(self.coef, self.domain, self.window)

    def degree(self) :
        """The degree of the series.

        Notes
        -----
        .. versionadded:: 1.5.0

        """
        return len(self) - 1

    def cutdeg(self, deg) :
        """Truncate series to the given degree.

        Reduce the degree of the $name series to `deg` by discarding the
        high order terms. If `deg` is greater than the current degree a
        copy of the current series is returned. This can be useful in least
        squares where the coefficients of the high degree terms may be very
        small.

        Parameters
        ----------
        deg : non-negative int
            The series is reduced to degree `deg` by discarding the high
            order terms. The value of `deg` must be a non-negative integer.

        Returns
        -------
        new_instance : $name
            New instance of $name with reduced degree.

        Notes
        -----
        .. versionadded:: 1.5.0

        """
        return self.truncate(deg + 1)

    def trim(self, tol=0) :
        """Remove small leading coefficients

        Remove leading coefficients until a coefficient is reached whose
        absolute value greater than `tol` or the beginning of the series is
        reached. If all the coefficients would be removed the series is set to
        ``[0]``. A new $name instance is returned with the new coefficients.
        The current instance remains unchanged.

        Parameters
        ----------
        tol : non-negative number.
            All trailing coefficients less than `tol` will be removed.

        Returns
        -------
        new_instance : $name
            Contains the new set of coefficients.

        """
        coef = pu.trimcoef(self.coef, tol)
        return self.__class__(coef, self.domain, self.window)

    def truncate(self, size) :
        """Truncate series to length `size`.

        Reduce the $name series to length `size` by discarding the high
        degree terms. The value of `size` must be a positive integer. This
        can be useful in least squares where the coefficients of the
        high degree terms may be very small.

        Parameters
        ----------
        size : positive int
            The series is reduced to length `size` by discarding the high
            degree terms. The value of `size` must be a positive integer.

        Returns
        -------
        new_instance : $name
            New instance of $name with truncated coefficients.

        """
        isize = int(size)
        if isize != size or isize < 1 :
            raise ValueError("size must be a positive integer")
        if isize >= len(self.coef) :
            coef = self.coef
        else :
            coef = self.coef[:isize]
        return self.__class__(coef, self.domain, self.window)

    def convert(self, domain=None, kind=None, window=None) :
        """Convert to different class and/or domain.

        Parameters
        ----------
        domain : array_like, optional
            The domain of the converted series. If the value is None,
            the default domain of `kind` is used.
        kind : class, optional
            The polynomial series type class to which the current instance
            should be converted. If kind is None, then the class of the
            current instance is used.
        window : array_like, optional
            The window of the converted series. If the value is None,
            the default window of `kind` is used.

        Returns
        -------
        new_series_instance : `kind`
            The returned class can be of different type than the current
            instance and/or have a different domain.

        Notes
        -----
        Conversion between domains and class types can result in
        numerically ill defined series.

        Examples
        --------

        """
        if kind is None:
            kind = $name
        if domain is None:
            domain = kind.domain
        if window is None:
            window = kind.window
        return self(kind.identity(domain, window=window))

    def mapparms(self) :
        """Return the mapping parameters.

        The returned values define a linear map ``off + scl*x`` that is
        applied to the input arguments before the series is evaluated. The
        map depends on the ``domain`` and ``window``; if the current
        ``domain`` is equal to the ``window`` the resulting map is the
        identity.  If the coeffients of the ``$name`` instance are to be
        used by themselves outside this class, then the linear function
        must be substituted for the ``x`` in the standard representation of
        the base polynomials.

        Returns
        -------
        off, scl : floats or complex
            The mapping function is defined by ``off + scl*x``.

        Notes
        -----
        If the current domain is the interval ``[l_1, r_1]`` and the window
        is ``[l_2, r_2]``, then the linear mapping function ``L`` is
        defined by the equations::

            L(l_1) = l_2
            L(r_1) = r_2

        """
        return pu.mapparms(self.domain, self.window)

    def integ(self, m=1, k=[], lbnd=None) :
        """Integrate.

        Return an instance of $name that is the definite integral of the
        current series. Refer to `${nick}int` for full documentation.

        Parameters
        ----------
        m : non-negative int
            The number of integrations to perform.
        k : array_like
            Integration constants. The first constant is applied to the
            first integration, the second to the second, and so on. The
            list of values must less than or equal to `m` in length and any
            missing values are set to zero.
        lbnd : Scalar
            The lower bound of the definite integral.

        Returns
        -------
        integral : $name
            The integral of the series using the same domain.

        See Also
        --------
        ${nick}int : similar function.
        ${nick}der : similar function for derivative.

        """
        off, scl = self.mapparms()
        if lbnd is None :
            lbnd = 0
        else :
            lbnd = off + scl*lbnd
        coef = ${nick}int(self.coef, m, k, lbnd, 1./scl)
        return self.__class__(coef, self.domain, self.window)

    def deriv(self, m=1):
        """Differentiate.

        Return an instance of $name that is the derivative of the current
        series.  Refer to `${nick}der` for full documentation.

        Parameters
        ----------
        m : non-negative int
            The number of integrations to perform.

        Returns
        -------
        derivative : $name
            The derivative of the series using the same domain.

        See Also
        --------
        ${nick}der : similar function.
        ${nick}int : similar function for integration.

        """
        off, scl = self.mapparms()
        coef = ${nick}der(self.coef, m, scl)
        return self.__class__(coef, self.domain, self.window)

    def roots(self) :
        """Return list of roots.

        Return ndarray of roots for this series. See `${nick}roots` for
        full documentation. Note that the accuracy of the roots is likely to
        decrease the further outside the domain they lie.

        See Also
        --------
        ${nick}roots : similar function
        ${nick}fromroots : function to go generate series from roots.

        """
        roots = ${nick}roots(self.coef)
        return pu.mapdomain(roots, self.window, self.domain)

    def linspace(self, n=100, domain=None):
        """Return x,y values at equally spaced points in domain.

        Returns x, y values at `n` equally spaced points across domain.
        Here y is the value of the polynomial at the points x.  This is
        intended as a plotting aid.

        Parameters
        ----------
        n : int, optional
            Number of point pairs to return. The default value is 100.

        Returns
        -------
        x, y : ndarrays
            ``x`` is equal to linspace(self.domain[0], self.domain[1], n)
            ``y`` is the polynomial evaluated at ``x``.

        .. versionadded:: 1.5.0

        """
        if domain is None:
            domain = self.domain
        x = np.linspace(domain[0], domain[1], n)
        y = self(x)
        return x, y



    @staticmethod
    def fit(x, y, deg, domain=None, rcond=None, full=False, w=None,
        window=$domain):
        """Least squares fit to data.

        Return a `$name` instance that is the least squares fit to the data
        `y` sampled at `x`. Unlike `${nick}fit`, the domain of the returned
        instance can be specified and this will often result in a superior
        fit with less chance of ill conditioning. See `${nick}fit` for full
        documentation of the implementation.

        Parameters
        ----------
        x : array_like, shape (M,)
            x-coordinates of the M sample points ``(x[i], y[i])``.
        y : array_like, shape (M,) or (M, K)
            y-coordinates of the sample points. Several data sets of sample
            points sharing the same x-coordinates can be fitted at once by
            passing in a 2D-array that contains one dataset per column.
        deg : int
            Degree of the fitting polynomial.
        domain : {None, [beg, end], []}, optional
            Domain to use for the returned $name instance. If ``None``,
            then a minimal domain that covers the points `x` is chosen.  If
            ``[]`` the default domain ``$domain`` is used. The default
            value is $domain in numpy 1.4.x and ``None`` in later versions.
            The ``'[]`` value was added in numpy 1.5.0.
        rcond : float, optional
            Relative condition number of the fit. Singular values smaller
            than this relative to the largest singular value will be
            ignored. The default value is len(x)*eps, where eps is the
            relative precision of the float type, about 2e-16 in most
            cases.
        full : bool, optional
            Switch determining nature of return value. When it is False
            (the default) just the coefficients are returned, when True
            diagnostic information from the singular value decomposition is
            also returned.
        w : array_like, shape (M,), optional
            Weights. If not None the contribution of each point
            ``(x[i],y[i])`` to the fit is weighted by `w[i]`. Ideally the
            weights are chosen so that the errors of the products
            ``w[i]*y[i]`` all have the same variance.  The default value is
            None.
            .. versionadded:: 1.5.0
        window : {[beg, end]}, optional
            Window to use for the returned $name instance. The default
            value is ``$domain``
            .. versionadded:: 1.6.0

        Returns
        -------
        least_squares_fit : instance of $name
            The $name instance is the least squares fit to the data and
            has the domain specified in the call.

        [residuals, rank, singular_values, rcond] : only if `full` = True
            Residuals of the least-squares fit, the effective rank of the
            scaled Vandermonde matrix and its singular values, and the
            specified value of `rcond`. For more details, see
            `linalg.lstsq`.

        See Also
        --------
        ${nick}fit : similar function

        """
        if domain is None:
            domain = pu.getdomain(x)
        elif domain == []:
            domain = $domain

        if window == []:
            window = $domain

        xnew = pu.mapdomain(x, domain, window)
        res = ${nick}fit(xnew, y, deg, w=w, rcond=rcond, full=full)
        if full :
            [coef, status] = res
            return $name(coef, domain=domain, window=window), status
        else :
            coef = res
            return $name(coef, domain=domain, window=window)

    @staticmethod
    def fromroots(roots, domain=$domain, window=$domain) :
        """Return $name instance with specified roots.

        Returns an instance of $name representing the product
        ``(x - r[0])*(x - r[1])*...*(x - r[n-1])``, where ``r`` is the
        list of roots.

        Parameters
        ----------
        roots : array_like
            List of roots.

        Returns
        -------
        object : $name
            Series with the specified roots.

        See Also
        --------
        ${nick}fromroots : equivalent function

        """
        if domain is None :
            domain = pu.getdomain(roots)
        rnew = pu.mapdomain(roots, domain, window)
        coef = ${nick}fromroots(rnew)
        return $name(coef, domain=domain, window=window)

    @staticmethod
    def identity(domain=$domain, window=$domain) :
        """Identity function.

        If ``p`` is the returned $name object, then ``p(x) == x`` for all
        values of x.

        Parameters
        ----------
        domain : array_like
            The resulting array must be if the form ``[beg, end]``, where
            ``beg`` and ``end`` are the endpoints of the domain.
        window : array_like
            The resulting array must be if the form ``[beg, end]``, where
            ``beg`` and ``end`` are the endpoints of the window.

        Returns
        -------
        identity : $name object

        """
        off, scl = pu.mapparms(window, domain)
        coef = ${nick}line(off, scl)
        return $name(coef, domain, window)
'''.replace('REL_IMPORT', rel_import))
python-numpy 1:1.6.1-6ubuntu1 / usr / share / pyshared / numpy / polynomial / polytemplate.py
