diff --git a/contents/thomas_algorithm/code/clisp/thomas.lisp b/contents/thomas_algorithm/code/clisp/thomas.lisp new file mode 100644 index 000000000..600cc0b1d --- /dev/null +++ b/contents/thomas_algorithm/code/clisp/thomas.lisp @@ -0,0 +1,34 @@ +;;;; Thomas algorithm implementation in Common Lisp + +(defmacro divf (place divisor) + "Divides the value at place by divisor" + `(setf ,place (/ ,place ,divisor))) + +(defun helper (v1 v2 v3 row) + (- (svref v1 row) (* (svref v2 row) (svref v3 (1- row))))) + +(defun thomas (diagonal-a diagonal-b diagonal-c last-column) + "Returns the solutions to a tri-diagonal matrix non-destructively" + ;; We have to copy the inputs to ensure non-destructiveness + (let ((a (copy-seq diagonal-a)) + (b (copy-seq diagonal-b)) + (c (copy-seq diagonal-c)) + (d (copy-seq last-column))) + (divf (svref c 0) (svref b 0)) + (divf (svref d 0) (svref b 0)) + (loop + for i from 1 upto (1- (length a)) do + (divf (svref c i) (helper b a c i)) + (setf (svref d i) (/ (helper d a d i) (helper b a c i)))) + (loop + for i from (- (length a) 2) downto 0 do + (decf (svref d i) (* (svref c i) (svref d (1+ i))))) + d)) + +(defparameter diagonal-a #(0 2 3)) +(defparameter diagonal-b #(1 3 6)) +(defparameter diagonal-c #(4 5 0)) +(defparameter last-column #(7 5 3)) + +;; should print 0.8666667 1.5333333 -0.26666668 +(format t "~{~f ~}~%" (coerce (thomas diagonal-a diagonal-b diagonal-c last-column) 'list)) diff --git a/contents/thomas_algorithm/thomas_algorithm.md b/contents/thomas_algorithm/thomas_algorithm.md index 042293fa1..5af6c3334 100644 --- a/contents/thomas_algorithm/thomas_algorithm.md +++ b/contents/thomas_algorithm/thomas_algorithm.md @@ -129,6 +129,8 @@ You will find this algorithm implemented [in this project](https://scratch.mit.e [import, lang:"crystal"](code/crystal/thomas.cr) {% sample lang="kotlin" %} [import, lang:"kotlin"](code/kotlin/thomas.kt) +{% sample lang="lisp" %} +[import, lang:"lisp"](code/clisp/thomas.lisp) {% sample lang="ruby" %} [import, lang="ruby"](code/ruby/thomas.rb) {% sample lang="js" %}