function f1(n) = (
	if n <= 0 then return null;
	set(`count,count+1);
	f1(n-1);
);

function f2(n) = (
	if n <= 0 then 
		null
	else (
	      set(`count,count+1);
	      f2(n-1)
	)
);

function randtest1() = (
  for n=1 to 500 do (
    a = rand () * 5 + 0.01;
    b = rand () * 5 + 0.01;
    xx = ((log(a) + log(b) <= 2*log(a+b)));
    if not IsBoolean(xx) then return false
  );
  true
);

# check randomness a bit, standard deviation should be low
function randtest2() = (
  a = [0,0,0,0,0];
  tests = 3000;
  for n=1 to tests do (
    k = randint (5);
    a@(k+1) = a@(k+1) + 1
  );
  sd = float (StandardDeviation(a) / (tests/5));
  if (sd > 0.11) then (
    print (a);
    print (sd);
    false
  ) else
    true
);

function derivtest() = (
	epsilon = 0.0001;
	for n in [`asin,`asinh,`acos,`acosh,`acsc,`acsch,`asec,`asech,`atan,`atanh,`acot,`acoth,`sin,`sinh,`cos,`cosh,`csc,`csch,`sec,`sech,`tan,`tanh,`cot,`coth,`sqrt,`ln,`log2,`log10,`exp,`cis] do (
		d = SymbolicDerivative(n);
		for x in [0.1,0.2,0.15,0.3,0.4] do (
			y = (d call (x)) - NDerivative (n,x);
			if (not IsValue (y)) or abs(y) > epsilon then (
				error ("Bad derivative on " + n + " at " + x + "!");
				return false
			)
		)
	);
	f = SymbolicNthDerivative(atan,4);
	d = SymbolicDerivative(f);
	for x in [0.1,0.2,0.15,0.3,0.4] do (
		y = (d call (x)) - NDerivative (f,x);
		if (not IsValue (y)) or abs(y) > epsilon then (
			error ("Bad derivative on 5th atan deriv at " + x + "!");
			return false
		)
	);
	true
);

function roottestcube() = (
  for n=1 to 100 do (
    if n <= 50 then
    	p = randint(10,4)-5 + 1i* (randint(10,4)-5)
    else
        p = randint(10,4)-5;
    if p@(4) == 0 then
	p@(4) = 1;
    #print (p);
    r = CubicFormula (p);
    #print (r);
    f = PolyToFunction (p);
    if (|f(r@(1))| >= 0.01 or
	|f(r@(2))| >= 0.01 or
	|f(r@(3))| >= 0.01) then (
	#print (n);
	#print (p);
	#print (r);
	#print (f(r@(1)));
	#print (f(r@(2)));
	#print (f(r@(3)));
      return false;
    )
  );
  true
);

function roottestquart() = (
  for n=1 to 100 do (
    if n <= 50 then
    	p = randint(10,5)-5 + 1i* (randint(10,5)-5)
    else
        p = randint(10,5)-5;
    if p@(5) == 0 then
	p@(5) = 1;
    #print(p);
    r = QuarticFormula (p);
    #print (r);
    f = PolyToFunction (p);
    if (|f(r@(1))| >= 0.01 or
	|f(r@(2))| >= 0.01 or
	|f(r@(3))| >= 0.01 or
	|f(r@(4))| >= 0.01) then (
	#print (n);
	#print (p);
	#print (r);
	#print (f(r@(1)));
	#print (f(r@(2)));
	#print (f(r@(3)));
	#print (f(r@(4)));
      return false;
    )
  );
  true
);

function nthroottest() = (
	for n=1 to 500 do (
		a=(rand()*10 - 5) ;
		ee = randint(10)+1;
		if Chop(( a^(1/ee) )^ee - a) != 0 then
			return false;
	);
	true
);

function polydivtest() = (
	for n=1 to 500 do (
		sz = randint (20)+5;
	   	p=randint(20,sz)-10;
		sz = randint (15)+5;
		do (
			q=randint(20,sz)-10
		) while IsZero(q);
		d = DividePoly(p,q,&r);
		if not AddPoly(MultiplyPoly(d,q),r) == TrimPoly (p) then (
			error ("Error on Dividepoly p: " + p + " q: " + q);
			return false
		)
	);
	true
);

function LongTest() = (
	errors = 0;

	#recursion test 1
	set(`count,0);
	f1(100);
	if count != 100 then (error("error on recursion test 1");errors = errors + 1);

	#recursion test 2
	set(`count,0);
	f1(500);
	if count != 500 then (error("error on recursion test 2");errors = errors + 1);

	#recursion test 3
	set(`count,0);
	f2(100);
	if count != 100 then (error("error on recursion test 3");errors = errors + 1);

	#random test 1
	if not randtest1() then (error("error on random test 1");errors = errors + 1);

	#random test 2
	if not randtest2() then (error("error on random test 2");errors = errors + 1);

	#derivative test
	if not derivtest() then (error("error on deriv test");errors = errors + 1);

	#roots
	if not roottestcube() then (error("error on root test cubic");errors = errors + 1);
	if not roottestquart() then (error("error on root test quartic");errors = errors + 1);
	if not nthroottest() then (error("error on nth root test");errors = errors + 1);

	# polynomial division
	if not polydivtest() then (error("error on poly div test");errors = errors + 1);

	# We print instead of return since we normally use load to run this
	# and load doesn't output the return value
	if errors > 0 then
		print(false)
	else
		print(true);
);

LongTest ()