3
$\begingroup$

I have a function that can evaluate ordinary generating function, and need to construct an approximation to the exponential generating function by using a few calls to the ordinary generating function. What functionality in Mathematica would make it easy?

A typical example is below. In this case I can EGF from OGF by using Laplace transform, but in general case this may be too expensive.

ogf = (159.774 - 1330.3 z + 5024.07 z^2 - 11359.6 z^3 + 17084.1 z^4 - 
   17943.2 z^5 + 13428.2 z^6 - 7159.87 z^7 + 2665.27 z^8 - 
   659.611 z^9 + 97.6646 z^10 - 6.55327 z^11)/(26.629 - 246.114 z + 
   1040.46 z^2 - 2660.31 z^3 + 4581.59 z^4 - 5598.57 z^5 + 
   4977.06 z^6 - 3243. z^7 + 1537.03 z^8 - 516.715 z^9 + 
   116.943 z^10 - 15.9963 z^11 + 1. z^12);
egfTrue = -InverseLaplaceTransform[ogf, z, -t];
Plot[egfTrue, {t, 1, 10}]

enter image description here

Motivation: this would allow efficient powers of diagonal+rank1 matrix by relying on $(I-A)^t\approx \exp(-tA)$ which is related to EGF. Meanwhile corresponding OGF takes $O(d)$ time to compute at a single point using trick with Woodbury formula,

$\endgroup$

1 Answer 1

2
$\begingroup$

Maybe this help:

Using Heaviside expansion formula we can speed up to 20X

ClearAll["`*"]; Remove["`*"];

ogf = (159.774 - 1330.3 z + 5024.07 z^2 - 11359.6 z^3 + 17084.1 z^4 - 
 17943.2 z^5 + 13428.2 z^6 - 7159.87 z^7 + 2665.27 z^8 - 
 659.611 z^9 + 97.6646 z^10 - 6.55327 z^11)/(26.629 - 246.114 z + 
 1040.46 z^2 - 2660.31 z^3 + 4581.59 z^4 - 5598.57 z^5 + 
 4977.06 z^6 - 3243. z^7 + 1537.03 z^8 - 516.715 z^9 + 
 116.943 z^10 - 15.9963 z^11 + 1. z^12);
ogfZ = Rationalize[ogf, 0];

INV = -Total[(Numerator[ogfZ]/D[Denominator[ogfZ], z]*Exp[-z*t] /. 
  NSolve[Denominator[ogfZ] == 0, z,WorkingPrecision -> 20])];(*Heaviside formula*)

egfTrue = -InverseLaplaceTransform[ogf, z, -t];

Plot[{egfTrue, INV}, {t, 1, 10}, PlotStyle -> {Red, {Dashed, Black}}]
Plot[egfTrue - INV // Re // Evaluate, {t, 1, 10}, WorkingPrecision -> 20] // Quiet
$\endgroup$
2
  • 1
    $\begingroup$ If we don't Rationalize and use Solve we can speed up to 70X. $\endgroup$ Mar 21, 2023 at 12:20
  • 1
    $\begingroup$ Neat! That has some curious resemblance to harder to implement solution here -- write $e^{-x}$ as $1/p(x)$ for polynomial $p$, express it as sum of terms like $r(t y_i)=(1-t y_i)^{-1}$, putting OGF in place of $r$ then somehow gives EGF $\endgroup$ Mar 21, 2023 at 14:54

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.