Assuming convergence (so the formula works at least for polynomial $y$), the formula can be seen by linear algebra. One has part of the infinite dimensional matrix as follows: $$\left[\begin{array}{cccccccc}*&y(0)&y'(0)&y''(0)&y'''(0)&y^{(4)}(0)&y^{(5)}(0)&y^{(6)}(0)\\
A&1&\frac 1 2&\frac 1 6&\fr...
They are just the coefficients in $A,B,\cdots$. Note that in $B$, they have a typo as well: It should read $B=y(1)-y(0)=y'(0)+\frac 1 2y''(0)+\cdots.$ The formulas for $A,B,C,\cdots$ follow from the Fundamental Theorem of Calculus.
@MaximusSu, It depends on the way you write the formula, but definitely the formula is related to Bernoulli numbers: oeis.org/… (I was too lazy to check it earlier, but when you time, you can try to derive the general terms using Bernoulli numbers.) I will double check everything within a day, depending on availability.
"I am confused on that part only, cuz it must match the exact formula stated in Wiki :))": Note that the formula in Wiki is NOT an identity (equality), it is just an asymptotic approximation, so it is irrelevant to drop or change a few terms.
First of all, the formula in your question is not the same formula as in the Wiki, so there is no reason to match the terms (some other terms disagree as well as you can check).
As for your question regarding asymptotic expansion, if f~g in the sense when b goes to infinity, then changing the term f(a) will not affect the result.
This Maclaurin formula derivation from taylor series should match the wikipedia one, the wikipedia one can help to get good closed form approximation of any sums
Discussion between Pythagoras and Max…
Discussion between Pythagoras and Max…