(defun arr*k (arr k) (let* ((arr-row-num (array-dimension arr 0)) (result-arr (make-array arr-row-num))) (dotimes (i arr-row-num result-arr) (setf (aref result-arr i) (* (aref arr i) k))))) (defun add-array (arr1 arr2) (let ((arr1-row-num (array-dimension arr1 0)) (arr2-row-num (array-dimension arr2 0))) (if (/= arr1-row-num arr2-row-num) (progn (format t "arr1-row-num[~A]:arr2-row-num[~A]~%" arr1-row-num arr2-row-num) (return-from add-array))) (let ((result-arr (make-array arr1-row-num))) (dotimes (i arr1-row-num result-arr) (setf (aref result-arr i) (+ (aref arr1 i) (aref arr2 i))))))) (defun set-row (arr1 arr2 i) (let ((arr1-col-num (array-dimension arr1 1)) (arr2-col-num (array-dimension arr2 0))) (if (<= arr1-col-num i) (progn (format t "arr1-col-num[~A]:i[~A]~%" arr1-col-num i) (return-from set-row))) (if (/= arr1-col-num arr2-col-num) (progn (format t "arr1-col-num[~A]:arr2-col-num[~A]~%" arr1-col-num arr2-col-num) (return-from set-row))) (dotimes (j arr1-col-num arr1) (setf (aref arr1 i j) (aref arr2 j))))) (defun aref-row (arr col) (let* ((arr-col (array-dimension arr 1)) (result-row (make-array (list arr-col)))) (dotimes (i arr-col result-row) (setf (aref result-row i) (aref arr col i))))) (defun swap-row (arr row1 row2) (let ((temp-row (aref-row arr row1))) (set-row arr (aref-row arr row2) row1) (set-row arr temp-row row2))) (defun gauss-jordan (arr) (let ((arr-row-num (array-dimension arr 0)) (arr-col-num (array-dimension arr 1))) ;; 入力チェック (if (/= (+ arr-row-num 1) arr-col-num) (progn (format t "arr-row-num[~A]:arr-col-num[~A]~%" arr-row-num arr-col-num) (return-from gauss-jordan))) (let ((result-arr (make-array (list arr-row-num arr-col-num)))) ;; 返り値行列領域初期化 (dotimes (i arr-row-num) (set-row result-arr (aref-row arr i) i)) ;; 返り値行列の対角成分が0でないように行を入れ換える (dotimes (i arr-row-num) (if (= 0 (aref result-arr i i)) (do* ((j (+ i 1) (+ j 1)) (swap-target (mod j arr-row-num) (mod j arr-row-num ))) ((/= (aref result-arr swap-target i) 0) (swap-row result-arr i swap-target)) (if (= (mod j arr-row-num) i) (return-from gauss-jordan))))) ;; 対角成分を除く(i,j)成分を0にする (dotimes (i arr-row-num) (dotimes (j arr-row-num) (if (/= i j) (let* ((flag (if (> (* (aref result-arr j i) (aref result-arr i i)) 0) -1 1))) (if (and (/= (aref result-arr i i) 0) (/= (aref result-arr j i) 0)) (set-row result-arr (add-array (arr*k (aref-row result-arr j) (abs (aref result-arr i i))) (arr*k (aref-row result-arr i) (* (abs (aref result-arr j i)) flag))) j)))))) ;; 対角成分を1にする (dotimes (i arr-row-num result-arr) (if (/= (aref result-arr i i) 0) (set-row result-arr (arr*k (aref-row result-arr i) (/ 1 (aref result-arr i i))) i))) ;; 単位行列か否かをチェックして,解が出たかどうかをチェック (dotimes (i arr-row-num) (dotimes (j arr-row-num) (if (not (or (and (= i j) (= (aref result-arr i j) 1)) (and (/= i j) (= (aref result-arr i j) 0)))) (return-from gauss-jordan)))) (return-from gauss-jordan result-arr)))) (defun test1 () (let ((a1 #2a((1 2 3) (4 5 6))) (a2 #2a((-1 2 3) (4 5 6))) (a3 #2a((1 2 3) (-4 5 6))) (a4 #2a((-1 2 3) (-4 5 6)))) (mapcar #'gauss-jordan (list a1 a2 a3 a4)))) (defun test2 () (let ((a1 #2a((1 2 3 4) (5 6 7 8) (9 10 11 12))) (a2 #2a((1 0 0 1) (0 1 0 2) (0 0 1 3))) (a3 #2a((1 0 1 1) (1 1 0 1) (0 1 1 1))) (a4 #2a((1 4 7 10) (2 4 6 8) (3 4 7 9)))) (mapcar #'gauss-jordan (list a1 a2 a3 a4)))) (defun test3 () (let ((a1 #2a((1 0 0 1) (0 1 0 2) (0 0 1 3))) (a2 #2a((1 0 0 1) (0 0 1 3) (0 1 0 2))) (a3 #2a((0 1 0 2) (1 0 0 1) (0 0 1 3))) (a4 #2a((0 1 0 2) (0 0 1 3) (1 0 0 1))) (a5 #2a((0 0 1 3) (1 0 0 1) (0 1 0 2))) (a6 #2a((0 0 1 3) (0 1 0 2) (1 0 0 1)))) (mapcar #'gauss-jordan (list a1 a2 a3 a4 a5 a6)))) (defun test4 () (let ((a1 #2a((1 0 0 0 1) (0 1 0 0 2) (0 0 1 0 3) (0 0 0 1 4))) (a2 #2a((1 0 0 0 1) (0 1 0 0 2) (0 0 0 1 4) (0 0 1 0 3))) (a3 #2a((1 0 0 0 1) (0 0 1 0 3) (0 1 0 0 2) (0 0 0 1 4))) (a4 #2a((1 0 0 0 1) (0 0 1 0 3) (0 0 0 1 4) (0 1 0 0 2))) (a5 #2a((1 0 0 0 1) (0 0 0 1 4) (0 0 1 0 3) (0 1 0 0 2))) (a6 #2a((1 0 0 0 1) (0 0 0 1 4) (0 1 0 0 2) (0 0 1 0 3))) (a7 #2a((0 1 0 0 2) (1 0 0 0 1) (0 0 1 0 3) (0 0 0 1 4))) (a8 #2a((0 1 0 0 2) (1 0 0 0 1) (0 0 0 1 4) (0 0 1 0 3))) (a9 #2a((0 1 0 0 2) (0 0 1 0 3) (0 0 0 1 4) (1 0 0 0 1))) (a10 #2a((0 1 0 0 2) (0 0 1 0 3) (1 0 0 0 1) (0 0 0 1 4))) (a11 #2a((0 1 0 0 2) (0 0 0 1 4) (0 0 1 0 3) (1 0 0 0 1))) (a12 #2a((0 1 0 0 2) (0 0 0 1 4) (1 0 0 0 1) (0 0 1 0 3))) (a13 #2a((0 0 1 0 3) (1 0 0 0 1) (0 1 0 0 2) (0 0 0 1 4))) (a14 #2a((0 0 1 0 3) (1 0 0 0 1) (0 0 0 1 4) (0 1 0 0 2))) (a15 #2a((0 0 1 0 3) (0 1 0 0 2) (1 0 0 0 1) (0 0 0 1 4))) (a16 #2a((0 0 1 0 3) (0 1 0 0 2) (0 0 0 1 4) (1 0 0 0 1))) (a17 #2a((0 0 1 0 3) (0 0 0 1 4) (0 1 0 0 2) (1 0 0 0 1))) (a18 #2a((0 0 1 0 3) (0 0 0 1 4) (1 0 0 0 1) (0 1 0 0 2))) (a19 #2a((0 0 0 1 4) (1 0 0 0 1) (0 1 0 0 2) (0 0 1 0 3))) (a20 #2a((0 0 0 1 4) (1 0 0 0 1) (0 0 1 0 3) (0 1 0 0 2))) (a21 #2a((0 0 0 1 4) (0 1 0 0 2) (1 0 0 0 1) (0 0 1 0 3))) (a22 #2a((0 0 0 1 4) (0 1 0 0 2) (0 0 1 0 3) (1 0 0 0 1))) (a23 #2a((0 0 0 1 4) (0 0 1 0 3) (1 0 0 0 1) (0 1 0 0 2))) (a24 #2a((0 0 0 1 4) (0 0 1 0 3) (0 1 0 0 2) (1 0 0 0 1)))) (mapcar #'gauss-jordan (list a1 a2 a3 a4 a5 a6 a7 a8 a9 a10 a11 a12 a13 a14 a15 a16 a17 a18 a19 a20 a21 a22 a23 a24)))) (defun test5 () (let ((a1 #2a((1 0 0 0 1) (0 1 0 0 2) (0 0 1 0 3) (1 1 1 1 4))) (a2 #2a((1 0 0 0 1) (0 1 0 0 2) (0 0 0 1 4) (1 1 1 1 3))) (a3 #2a((1 0 0 0 1) (0 0 1 0 3) (0 1 0 0 2) (1 1 1 1 4))) (a4 #2a((1 0 0 0 1) (0 0 1 0 3) (0 0 0 1 4) (1 1 1 1 2))) (a5 #2a((1 0 0 0 1) (0 0 0 1 4) (0 0 1 0 3) (1 1 1 1 2))) (a6 #2a((1 0 0 0 1) (0 0 0 1 4) (0 1 0 0 2) (1 1 1 1 3))) (a7 #2a((0 1 0 0 2) (1 0 0 0 1) (0 0 1 0 3) (1 1 1 1 4))) (a8 #2a((0 1 0 0 2) (1 0 0 0 1) (0 0 0 1 4) (1 1 1 1 3))) (a9 #2a((0 1 0 0 2) (0 0 1 0 3) (0 0 0 1 4) (1 1 1 1 1))) (a10 #2a((0 1 0 0 2) (0 0 1 0 3) (1 0 0 0 1) (1 1 1 1 4))) (a11 #2a((0 1 0 0 2) (0 0 0 1 4) (0 0 1 0 3) (1 1 1 1 1))) (a12 #2a((0 1 0 0 2) (0 0 0 1 4) (1 0 0 0 1) (1 1 1 1 3))) (a13 #2a((0 0 1 0 3) (1 0 0 0 1) (0 1 0 0 2) (1 1 1 1 4))) (a14 #2a((0 0 1 0 3) (1 0 0 0 1) (0 0 0 1 4) (1 1 1 1 2))) (a15 #2a((0 0 1 0 3) (0 1 0 0 2) (1 0 0 0 1) (1 1 1 1 4))) (a16 #2a((0 0 1 0 3) (0 1 0 0 2) (0 0 0 1 4) (1 1 1 1 1))) (a17 #2a((0 0 1 0 3) (0 0 0 1 4) (0 1 0 0 2) (1 1 1 1 1))) (a18 #2a((0 0 1 0 3) (0 0 0 1 4) (1 0 0 0 1) (1 1 1 1 2))) (a19 #2a((0 0 0 1 4) (1 0 0 0 1) (0 1 0 0 2) (1 1 1 1 3))) (a20 #2a((0 0 0 1 4) (1 0 0 0 1) (0 0 1 0 3) (1 1 1 1 2))) (a21 #2a((0 0 0 1 4) (0 1 0 0 2) (1 0 0 0 1) (1 1 1 1 3))) (a22 #2a((0 0 0 1 4) (0 1 0 0 2) (0 0 1 0 3) (1 1 1 1 1))) (a23 #2a((0 0 0 1 4) (0 0 1 0 3) (1 0 0 0 1) (1 1 1 1 2))) (a24 #2a((0 0 0 1 4) (0 0 1 0 3) (0 1 0 0 2) (1 1 1 1 1)))) (mapcar #'gauss-jordan (list a1 a2 a3 a4 a5 a6 a7 a8 a9 a10 a11 a12 a13 a14 a15 a16 a17 a18 a19 a20 a21 a22 a23 a24)))) (defun test6 () (let ((a1 #2a((1 0 0 0 1) (0 1 0 0 2) (1 1 1 1 3) (0 0 0 1 4))) (a2 #2a((1 0 0 0 1) (0 1 0 0 2) (1 1 1 1 4) (0 0 1 0 3))) (a3 #2a((1 0 0 0 1) (0 0 1 0 3) (1 1 1 1 2) (0 0 0 1 4))) (a4 #2a((1 0 0 0 1) (0 0 1 0 3) (1 1 1 1 4) (0 1 0 0 2))) (a5 #2a((1 0 0 0 1) (0 0 0 1 4) (1 1 1 1 3) (0 1 0 0 2))) (a6 #2a((1 0 0 0 1) (0 0 0 1 4) (1 1 1 1 2) (0 0 1 0 3))) (a7 #2a((0 1 0 0 2) (1 0 0 0 1) (1 1 1 1 3) (0 0 0 1 4))) (a8 #2a((0 1 0 0 2) (1 0 0 0 1) (1 1 1 1 4) (0 0 1 0 3))) (a9 #2a((0 1 0 0 2) (0 0 1 0 3) (1 1 1 1 4) (1 0 0 0 1))) (a10 #2a((0 1 0 0 2) (0 0 1 0 3) (1 1 1 1 1) (0 0 0 1 4))) (a11 #2a((0 1 0 0 2) (0 0 0 1 4) (1 1 1 1 3) (1 0 0 0 1))) (a12 #2a((0 1 0 0 2) (0 0 0 1 4) (1 1 1 1 1) (0 0 1 0 3))) (a13 #2a((0 0 1 0 3) (1 0 0 0 1) (1 1 1 1 2) (0 0 0 1 4))) (a14 #2a((0 0 1 0 3) (1 0 0 0 1) (1 1 1 1 4) (0 1 0 0 2))) (a15 #2a((0 0 1 0 3) (0 1 0 0 2) (1 1 1 1 1) (0 0 0 1 4))) (a16 #2a((0 0 1 0 3) (0 1 0 0 2) (1 1 1 1 4) (1 0 0 0 1))) (a17 #2a((0 0 1 0 3) (0 0 0 1 4) (1 1 1 1 2) (1 0 0 0 1))) (a18 #2a((0 0 1 0 3) (0 0 0 1 4) (1 1 1 1 1) (0 1 0 0 2))) (a19 #2a((0 0 0 1 4) (1 0 0 0 1) (1 1 1 1 2) (0 0 1 0 3))) (a20 #2a((0 0 0 1 4) (1 0 0 0 1) (1 1 1 1 3) (0 1 0 0 2))) (a21 #2a((0 0 0 1 4) (0 1 0 0 2) (1 1 1 1 1) (0 0 1 0 3))) (a22 #2a((0 0 0 1 4) (0 1 0 0 2) (1 1 1 1 3) (1 0 0 0 1))) (a23 #2a((0 0 0 1 4) (0 0 1 0 3) (1 1 1 1 1) (0 1 0 0 2))) (a24 #2a((0 0 0 1 4) (0 0 1 0 3) (1 1 1 1 2) (1 0 0 0 1)))) (mapcar #'gauss-jordan (list a1 a2 a3 a4 a5 a6 a7 a8 a9 a10 a11 a12 a13 a14 a15 a16 a17 a18 a19 a20 a21 a22 a23 a24))))