|
|
1.1 root 1: #print
2: Write a function
3: cubrt(x)
4: that returns the cube root of a floating point number.
5: Put it on a file named "cubrt.c"; compile and test it,
6: then type "ready".
7: (If you don't know how to compute cube roots, try Newton's method).
8: #once #create reldif.c
9: double reldif(a,b)
10: double a,b;
11: {
12: double c,d;
13: if (a==0. && b==0.) return(0.);
14: c = a>0 ? a : -a;
15: d = b>0 ? b : -b;
16: c = c>d ? c : d;
17: return( (a-b)/c );
18: }
19: #once #create tzaqc.c
20: main()
21: {
22: double cubrt();
23: double reldif();
24: double a, b, eps;
25:
26: a = 8.;
27: b = 2.;
28: eps = 1e-5;
29: if (reldif(cubrt(a), b) > eps)
30: exit(1);
31:
32: a = 0.;
33: b = 0.;
34: if (reldif(cubrt(a), b) > eps)
35: exit(1);
36:
37: a = 1e6;
38: b = 1e2;
39: if (reldif(cubrt(a), b) > eps)
40: exit(1);
41: exit(0);
42: }
43: #user
44: cc tzaqc.c cubrt.o reldif.c
45: a.out
46: #succeed
47: /* one way */
48: double cubrt(x)
49: double x;
50: {
51: /* Newton's method: x <- x - (x**3-a)/(3*x*x) */
52: double y, yn, dabs();
53: y = 0.;
54: yn = x;
55: while (dabs(y-yn) > y*1.e-8) {
56: y = yn;
57: yn = y - (y*y*y-x)/(3*y*y);
58: }
59: return(yn);
60: }
61:
62: double dabs(x)
63: double x;
64: {
65: return(x>0 ? x : -x);
66: }
67: #log
68: #next
69: 50.1a 10
70: 43.1b 5
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.