[SciPy-dev] Generic polynomials class (was Re: Volunteer for Scipy Project)
Charles R Harris
Tue Oct 13 19:15:37 CDT 2009
On Tue, Oct 13, 2009 at 5:55 PM, Anne Archibald
> 2009/10/13 Charles R Harris <firstname.lastname@example.org>:
> >> I'm not sure I see why returning all those NotImplementedErrors is a
> >> problem - yes, the code needs to be written to override the methods
> >> that return NotImplementedError, but that code needs to be written no
> >> matter how we set up an interface to it.
> > It's not trivial to do a complete workup. My template class is 698 lines
> > long, but it exists once in one place, and that does the trick for *all*
> > classes.
> Aha. I think I see where the confusion is. In my proposed interface,
> the coefficient-oriented functions would *be* the ones in the Basis
> object. So, for example,
> chebadd(c1,c2,interval) -> ChebyshevBasis(interval).add(c1, c2)
> or, more likely, you'd have somewhere earlier
> cheb = ChebyshevBasis(interval)
> and then it's
> chebadd(c1,c2,interval) -> cheb.add(c1,c2)
> So the coefficient-oriented interface that we need to have is what
> lives in the basis objects, and it is what the Polynomial object uses
> to implement (e.g.) p1*p2. There would be no separate chebmul,
> polymul, lagrangemul, ....
But the separate functions are what we want. They provide the most flexible
low level implementation and make the fewest assumptions as to what the
programmer might want to use them for. You don't really have a base class,
you essentially have classes derived from different base classes, which
isn't so different from what I am proposing.
> >> And having the base class
> >> allows us to take any code that is the same between bases - like the
> >> companion matrix code - and put it in the base class. So the only
> >> "wasted" code here is the one-line wrappers that return
> >> NotImplementedError, and these serve as documentation for users - any
> >> such method is something they're going to have to implement for their
> >> own polynomial representations. If you look at the three different
> >> representations I have currently implemented, there isn't really much
> >> repeated code.
> > But it is also incomplete. And you have to have the separate functions
> > anyway, so it is easier if you just need to write them and then you are
> > done. If you want to change or fix the class interface, you don't have to
> > touch the implementations of the functions. If you want to change the
> > implementations, then you don't need to touch the interface. The two
> > become mostly independent, which is generally a good thing.
> The Basis object is effectively just a namespace holding the
> collection of functions, and possibly any data they need (interval,
> for example). So the separation you're describing is still present in
> my proposed interface; it's the separation between Basis objects
> (which contain the functions but no particular plumbing) and the
> Polynomial objects (which are all plumbing).
Why not let the module hold the functions? It is the more natural namespace
> >> > I think your method of getting the companion matrix is clever, and
> >> > having
> >> > definitions for zero, one, and x is a useful touch that should
> >> > be
> >> > part of any interface. I'm less convinced by the use of an empty list
> >> > for
> >> > zero.
> >> This choice is a natural one in some ways, because it means you never
> >> normally have a leading coefficient of zero for polynomials of any
> >> degree. It also means degree = len(coefficients)-1, apart from
> >> polynomials that are represented in a degree-inflated manner. And the
> >> places you need to special-case it are almost all places you would
> >> need to special-case zero polynomials anyway. I realize it may not be
> >> so convenient for a vector-of-coefficients interface, but I think even
> >> there it's the natural approach.
> > But does it avoid special casing? Can you add it transparently to another
> > array?
> > In : x = arange(3)
> > In : x + 
> > ValueError Traceback (most recent call
> > /home/charris/<ipython console> in <module>()
> > ValueError: shape mismatch: objects cannot be broadcast to a single shape
> > In : x + 
> > Out: array([0, 1, 2])
> This is a special case anyway - it's one of the very few instances in
> which you can add coefficients the way you're describing. Normally you
> have to extend one or the other so that they're the same length (and
> in some bases this is a non-trivial operations, rather than simply
> padding with zeros). If what you mean here is to add a scalar, there's
> no reason to wrap it in a list;
Consistency. Usually these things are thought of as algebras with a unit
over the scalars.
> In : np.arange(3)+0
> Out: array([0, 1, 2])
> In : np.arange(3)+2
> Out: array([2, 3, 4])
> If what you're adding is actually a vector of coefficients, you should
> be surprised if it works:
> In : np.arange(3)+np.arange(2)
> ValueError Traceback (most recent call last)
> /home/peridot/<ipython console> in <module>()
> ValueError: shape mismatch: objects cannot be broadcast to a single shape
> In : np.arange(3)+np.arange(1)
> Out: array([0, 1, 2])
> > I think it is a complication.
> It's a question of trading one set of special cases for another. If
> the zero polynomial is the only polynomial with zero leading
> coefficient, that's going to force a fair number of special cases as
> well. (Polynomials with zero leading coefficient that are generated as
> a result of calculations are trouble of a different sort, since it's
> entirely likely that the leading coefficient will be corrupted by
> roundoff error and come out to something near zero, where the notion
> of "near" is unclear.)
Yes, that is always going to be a problem.
> In any case, it's only really important for people who want to use the
> coefficient-array interface; the polynomial class can hide this detail
> under the rug. I'm not terribly attached to this representation.
-------------- next part --------------
An HTML attachment was scrubbed...
More information about the Scipy-dev