cal/cal/root.cal

{ RCSid $Id: root.cal,v 1.2 2010/09/04 19:26:09 greg Exp $ }
{
        root.cal - calculate zeroes of functions using Newton's method.

        2/3/88
}

FTINY : 1e-9;

root(f, x0) = root2(f, x0, x0-f(x0)/lim(f1(f,x0)), nit);

root2(f, x0, x1, i) = if(i,
                        if(err-abs(x1-x0),
                                x1,
                                root2(f,x1,x1-f(x1)/lim(f1(f,x1)),i-1)),
                        0);

abs(x) = if(x, x, -x);
lim(x) = if(x, max(x,err), -max(-x,err));
max(a,b) = if(a-b, a, b);

f1(f,x) = (f(x+FTINY)-f(x-FTINY))/(2*FTINY);    { numerical derivative }

err = 1e-5;
nit = 100;
Revision:	1.3
Committed:	Tue Oct 11 15:40:54 2011 UTC (12 years, 6 months ago) by greg
Branch:	MAIN
CVS Tags:	rad5R2, rad4R2P2, rad5R0, rad5R1, rad4R2, rad4R1, rad4R2P1
Changes since 1.2:	+3 -3 lines
Log Message:	Changed default epsilon values
#	Content
1	{ RCSid $Id: root.cal,v 1.2 2010/09/04 19:26:09 greg Exp $ }
2	{
3	root.cal - calculate zeroes of functions using Newton's method.
4
5	2/3/88
6	}
7
8	FTINY : 1e-9;
9
10	root(f, x0) = root2(f, x0, x0-f(x0)/lim(f1(f,x0)), nit);
11
12	root2(f, x0, x1, i) = if(i,
13	if(err-abs(x1-x0),
14	x1,
15	root2(f,x1,x1-f(x1)/lim(f1(f,x1)),i-1)),
16	0);
17
18	abs(x) = if(x, x, -x);
19	lim(x) = if(x, max(x,err), -max(-x,err));
20	max(a,b) = if(a-b, a, b);
21
22	f1(f,x) = (f(x+FTINY)-f(x-FTINY))/(2*FTINY); { numerical derivative }
23
24	err = 1e-5;
25	nit = 100;