fun polyadd(p1,p2:real list)=
if length(p1)=0 then p2
else if length(p2)=0 then p1
else (hd(p1)+hd(p2))::polyadd(tl(p1),tl(p2));

fun scalarmult(l :real list, k:real)=
if length(l)=0 then nil
else (k*hd(l))::scalarmult(tl(l),k);

fun polydiff(p1,p2 :real list)=
if length(p1)=0 then scalarmult(p2,~1.0)
else if length(p2)=0 then p1
else (hd(p1)-hd(p2))::polydiff(tl(p1),tl(p2));

fun polymult(l1,l2 :real list)=
if length(l1)=0 then nil 
else 
polyadd(scalarmult(l2,hd(l1)),0.0::polymult(tl(l1),l2)); 

 fun polypower(l:real list,n :int)=
	let fun polypower_iter(l,n,c)=
		if n=0 then [1.0]
		else if n=1 then polymult(l,c)
		else if n mod 2 = 1 
		     then polypower_iter(polymult(l,l),n div 2, 
			polymult(c,l))
		else polypower_iter(polymult(l,l),n div 2, c) 
	in polypower_iter(l,n,[1.0]) end;  

fun computepoly(l:real list,x :real)=
	let fun eval_iter(s :real list,a,p,xpn :real)=
 	if length(s)=0 then p 
	else eval_iter(tl(s),a,p+hd(s)*xpn,a*xpn)
	in eval_iter(l,x,0.0,1.0) end;


fun differen(l :real list)=
	let fun d1(a :real list,n:int)=if 
	length(a)=0 then [] else real(n)*hd(a)::d1(tl(a),n+1)
	in d1(tl(l),1) end;

fun solve(l :real list,epsilon)=
	let
fun mod1(x :real)=if x <0.0 then ~x else x;
fun 
solve_iter(m:real list,assum :real,n :int)=if 
mod1(computepoly(m,assum))<epsilon 
then 
assum else  if n>100 then assum else  
solve_iter(m,assum-computepoly(m,assum)/computepoly(differen(m),assum),n+1)
in
solve_iter(l,0.0,1) end;

fun definteg(p :real list, a,b :real)=
let fun defin1(q :real list, c,d :real, n :int)=
if length(q)=0 then 0.0 else  
hd(q)*(Math.pow(d,real(n+1))-Math.pow( 
c,real(n+1)))/real(n+1)+defin1(tl(q),c,d,n+1)
in defin1(p,a,b,0) end;