Phase 2 of the coding period has started. This week has gone in wrapping up the left-over work of FormalPowerSeries. I had a meeting with Sartaj on Tuesday 25th of June, about the work left to be done on FormalPowerSeries module. We agreed that some minor changes need to be done on the PRs that have already been published. And also the class PR needs to be published. Here are the deliverables that have been completed this week.

  • PR #17011, which returns the connvoluted formal power series has been refactored. According to Kalevi’s comment here, and Sartaj’s follow-up discussion here, we felt that all the various forms of convolutions are not really that applicable to formal power series, backed by the evidence that there is no prior literature describing these operations for power series. Thus I restricted the function to only the linear infinite cycle case, and removed all other methods of convolution. The function name was changed to product, and now it returns something like this –
>>> from sympy import fps, sin, exp, convolution
>>> from sympy.abc import x
>>> f1 = fps(sin(x))
>>> f2 = fps(exp(x))

>>> f1.product(f2, x, 4)
x + x**2 + x**3/3 + O(x**4)

  • PR #17064, which was published to create compose and inverse functions, were changed a bit. Initially, all the sequences were being created and used in the functions, which were not at all optimal. After Sartaj’s valuable suggestions, I created all the possibe sequences in the __init__ function of the FormalPowerSeries class itself.
def __init__(self, *args):
    ak = args[4][0]
    k = ak.variables[0]
    self.ak_seq = sequence(ak.formula, (k, 1, oo))
    self.fact_seq = sequence(factorial(k), (k, 1, oo))
    self.bell_coeff_seq = self.ak_seq * self.fact_seq
    self.sign_seq = sequence((-1, 1), (k, 1, oo))

As we can see here, fact_seq, which is a factorial sequence, sign_seq, which is a sign alternating sequence, and a custom bell_coeff_seq have been defined in __init__. Since sequences are lazy, defining them at object creation time won’t affect performance.

Also, the code in bothe the compose and inverse required some optimization changes. More can be understood by following the comments and commits in the PR.

  • Till now, all the operations have been returning truncated terms. The problem is, it breaks the functionality of the FormalPowerSeries class, i.e a function defined in the FormalPowerSeries class should return a related object. But the problem with that was, it was very difficult to represent the coefficient sequence of the resultant composed or inverted formal power series. So, after a huge amount of discussion, I decided to create a FiniteFormalPowerSeries class, which will actually return a list of the coefficients of the coefficient sequence, based on the order n which the user provides.

Based on that, I creates a PR #17134. This PR consists of the implementation of FiniteFormalPowerSeries class, in which the coefficient sequence is a list of real numbers, computed based on the algorithm of the function calling it. Presently, compose and inverse, which are functions of FormalPowerSeries class, return a FiniteFormalPowerSeries object.

An example would be –

>>> from sympy import *
>>> x = symbols('x')
>>> f1, f2 = fps(exp(x)), fps(sin(x))

>>> f1.compose(f2)
FiniteFormalPowerSeries(exp(sin(x)), x, 0, 1, ([1, 1/2, 0, -1/8, -1/15], SeqFormula(x**_k, (_k, 0, oo)), 1))

>>> f1.compose(f2).truncate(6)
1 + x + x**2/2 - x**4/8 - x**5/15 + O(x**6)

>>> f1.inverse()
FiniteFormalPowerSeries(exp(-x), x, 0, 1, ([-1, 1/2, -1/6, 1/24, -1/120], SeqFormula(x**_k, (_k, 0, oo)), 1))

>>> f1.inverse().truncate(6)
1 - x + x**2/2 - x**3/6 + x**4/24 - x**5/120 + O(x**6)

The implementation of the FiniteFormalPowerSeries is somewhat like this :-

class FiniteFormalPowerSeries(FormalPowerSeries):
    def __new__(cls, *args):
        args = map(sympify, args)
        return Expr.__new__(cls, *args)

    def __init__(self, *args):
        pass

    def _eval_term(self, pt):
        try:
            pt_xk = self.xk.coeff(pt)
            pt_ak = S.Zero
            if pt is not S.Zero:
                pt_ak = self.ak[pt-1]  # Simplify the coefficients
        except IndexError:
            term = S.Zero
        else:
            term = (pt_ak * pt_xk)

        if self.ind:
            ind = S.Zero
            for t in Add.make_args(self.ind):
                pow_x = self._get_pow_x(t)
                if pt == 0 and pow_x < 1:
                    ind += t
                elif pow_x >= pt and pow_x < pt + 1:
                    ind += t
            term += ind

        return term.collect(self.x)

    def polynomial(self, n=6):
        """Truncated series as polynomial.
        Returns series expansion of ``f`` upto order ``O(x**n)``
        as a polynomial(without ``O`` term).
        """
        if (n > len(self.ak)+1):
            raise ValueError("The order of truncate should be less than or equal to the"
                "order specified in the compose function.")
        terms = []
        for i, t in enumerate(self):
            xp = self._get_pow_x(t)
            if xp >= n:
                break
            elif xp.is_integer is True and i == n + 1:
                break
            elif t is not S.Zero:
                terms.append(t)

        return Add(*terms)

As you can see, only the _eval_term and the polynomial functions have been overrriden, as compared to the FormalPowerSeries, from which it is inherited. Documentation and tests remain to be added, so have been added in the TODO section.

That’s it for Week 5 !! See you in the next week !!