Skip to content

Instantly share code, notes, and snippets.

@shirok
Created September 20, 2023 22:48
Show Gist options
  • Save shirok/5e50e8f5146ebd4bb94934ca8ae82c54 to your computer and use it in GitHub Desktop.
Save shirok/5e50e8f5146ebd4bb94934ca8ae82c54 to your computer and use it in GitHub Desktop.
;; f(z) = z^3 + 2z^2 + 1
(define (f z) (+ (expt z 3) (* 2 (square z)) 1))
;; f'(z) = 3z^2 + 4z
(define (df z) (+ (* 3 (square z)) (* 4 z)))
;; ニュートン法。初期値zから出発。変化が充分少なくなったらその時のzがf(z)=0の解。
;; 微分が0に近くなった場合はエラーとする。
;; 循環してしまうケースは検出しない (止まらなくなる)
(define (newton z)
(let ([fz (f z)]
[dfz (df z)])
(when (< (magnitude dfz) 1.0e-15)
(error "Unstable"))
(if (< (magnitude (/. fz dfz)) 1.0e-15)
z
(newton (- z (/. fz dfz))))))
#|
;; 実行例。
gosh> (newton 1)
-2.2055694304005904
gosh> (newton 1+i)
0.10278471520029513+0.6654569511528134i
gosh> (newton -1+i)
0.10278471520029517+0.6654569511528134i
;; 解であることの確認 (誤差があるので完全に0にはならない)
gosh> (f -2.2055694304005904)
0.0
gosh> (f 0.10278471520029513+0.6654569511528134i)
1.1102230246251565e-16-1.1102230246251565e-16i
gosh> (f 0.10278471520029517+0.6654569511528134i)
1.1102230246251565e-16
|#
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment