Maclaurin series of Black Box Function Given a Black Box function ff that takes a single floating-point parameter, write a function M to generate the Maclaurin series for that function, according to the following specifications:


take a parameter n, and generate all terms from x^0 to x^n where 0 &lt;= n &lt;= 10
coefficients for each term should be given to at least 5 decimal places for non-integral values; underflows can simply be ignored or replaced with 0 and left out if you wish, or you can use scientific notation and more decimal places
assume ff is a continuous function and everywhere/infinitely differentiable for x &gt; -1, and returns values for input x in the range: -10^9 &lt;= ff(x) &lt;= 10^9
programming languages with functions as first-class data types can take ff as a parameter; all other languages assume a function called ff has already been defined.
output should be formatted like # + #x + #x^2 + #x^5 - #x^8 - #x^9 ..., i.e., terms having 0 coefficient should be left out, and negative coefficients should be written - #x^i instead of + -#x^i; if the output terms are all zero, you can output either 0 or null output
output can be generated either to stdout or as a return value of type string from the function


Examples (passing ff as a first parameter):

&gt; M(x =&gt; sin(x), 8)
x - 0.16667x^3 + 0.00833x^5 - 0.00020x^7
&gt; M(x =&gt; e^x, 4)
1 + x + 0.5x^2 + 0.16667x^3 + 0.04167x^4
&gt; M(x =&gt; sqrt(1+x), 4)
1 + 0.5x - 0.125x^2 + 0.0625x^3 - 0.03906x^4
&gt; M(x =&gt; sin(x), 0)
0 (or blank is also accepted)
Python, 188 chars

D=lambda f:lambda x:100*(f(x+.005)-f(x-.005))
def M(f,n):
 s="";d=1
 for i in xrange(n+1):
  c=f(0)/d;f=D(f);d*=i+1
  if c:s+='%+.5f'%c+['','x','x^%d'%i][(i&gt;0)+(i&gt;1)]
 print s[s[:1]=='+':]


This code takes lots of liberties with accuracy and numerical stability, but it generates something close to the right answer for the given examples.

Running it:

import math
M(lambda x:math.sin(x),8)
M(lambda x:math.exp(x),4)
M(lambda x:math.sqrt(x+1),4)
M(lambda x:math.sin(x),0)


gives:

1.00000x-0.16666x^3+0.00833x^5-0.00020x^7
1.00000+1.00000x+0.50000x^2+0.16667x^3+0.04167x^4
1.00000+0.50000x-0.12500x^2+0.06251x^3-0.03907x^4
&lt;a blank line&gt;