diff options
Diffstat (limited to 'root.scm')
| -rw-r--r-- | root.scm | 66 |
1 files changed, 66 insertions, 0 deletions
@@ -149,3 +149,69 @@ df+sqrt-H df-sqrt-H))))) (loop next-z (f next-z)))))))) + +(define (secant:find-root-1 f x0 x1 prec must-bracket?) + (letrec ((stop? + (cond ((procedure? prec) prec) + ((and (integer? prec) (negative? prec)) + (lambda (x0 x1 fmax count) + (>= count (- prec)))) + (else + (lambda (x0 f0 x1 f1 count) + (and (< (abs f0) prec) + (< (abs f1) prec)))))) + (bracket-iter + (lambda (xlo flo glo xhi fhi ghi count) + (define (step xnew fnew) + (cond ((or (= xnew xlo) + (= xnew xhi)) + (let ((xmid (+ xlo (* 1/2 (- xhi xlo))))) + (if (= xnew xmid) + xmid + (step xmid (f xmid))))) + ((positive? fnew) + (bracket-iter xlo flo (if glo (* 0.5 glo) 1) + xnew fnew #f + (+ count 1))) + (else + (bracket-iter xnew fnew #f + xhi fhi (if ghi (* 0.5 ghi) 1) + (+ count 1))))) + (if (stop? xlo flo xhi fhi count) + (if (> (abs flo) (abs fhi)) xhi xlo) + (let* ((fflo (if glo (* glo flo) flo)) + (ffhi (if ghi (* ghi fhi) fhi)) + (del (- (/ fflo (- ffhi fflo)))) + (xnew (+ xlo (* del (- xhi xlo)))) + (fnew (f xnew))) + (step xnew fnew)))))) + (let ((f0 (f x0)) + (f1 (f x1))) + (cond ((<= f0 0 f1) + (bracket-iter x0 f0 #f x1 f1 #f 0)) + ((<= f1 0 f0) + (bracket-iter x1 f1 #f x0 f0 #f 0)) + (must-bracket? #f) + (else + (let secant-iter ((x0 x0) + (f0 f0) + (x1 x1) + (f1 f1) + (count 0)) + (cond ((stop? x0 f0 x1 f1 count) + (if (> (abs f0) (abs f1)) x1 x0)) + ((<= f0 0 f1) + (bracket-iter x0 f0 #f x1 f1 #f count)) + ((>= f0 0 f1) + (bracket-iter x1 f1 #f x0 f0 #f count)) + ((= f0 f1) #f) + (else + (let* ((xnew (+ x0 (* (- (/ f0 (- f1 f0))) (- x1 x0)))) + (fnew (f xnew)) + (fmax (max (abs f1) (abs fnew)))) + (secant-iter x1 f1 xnew fnew (+ count 1))))))))))) + +(define (secant:find-root f x0 x1 prec) + (secant:find-root-1 f x0 x1 prec #f)) +(define (secant:find-bracketed-root f x0 x1 prec) + (secant:find-root-1 f x0 x1 prec #t)) |