Форум программистов, компьютерный форум, киберфорум
Lisp
Войти
Регистрация
Восстановить пароль
Карта форума Темы раздела Блоги Сообщество Поиск Заказать работу  
 
Рейтинг 4.79/29: Рейтинг темы: голосов - 29, средняя оценка - 4.79
0 / 0 / 0
Регистрация: 24.03.2016
Сообщений: 4
1

Метод Крамера или метод Гаусса. Реализация

24.03.2016, 14:29. Показов 5751. Ответов 3
Метки нет (Все метки)

Author24 — интернет-сервис помощи студентам
Доброго времени суток. Стоит задача написать метод Гаусса или метод Крамера для решения СЛАУ на lisp, как это сделать даже на уровне алгоритма, я не знаю. На вход подается матрица, нужно будет работать с большими матрицами, например 10 на 10. На Lisp никогда не приходилось кодить) Помогите, пожалуйста, очень срочно.
0
Лучшие ответы (1)
Programming
Эксперт
94731 / 64177 / 26122
Регистрация: 12.04.2006
Сообщений: 116,782
24.03.2016, 14:29
Ответы с готовыми решениями:

Метод Гаусса или метод Крамера. Реализация
Доброго времени суток. Стоит задача написать метод Гаусса или метод Крамера для решения СЛАУ на F#,...

СЛАУ. Метод обратной матрицы, метод Гаусса, метод Крамера, метод Зейделя
Помогите ребят. Не могу построить алгоритмы для этих методов Язык C++

Метод Крамера и Гаусса для n неизвестных
Здравствуйте.Мне надо написать визуальную программу по решению слау методом крамера и гаусса.Метод...

Метод Крамера или обратной матрцы!
Всем привет! мне надо написать программу для решения уранений метод Крамера или обратной матрцы, но...

3
Модератор
Эксперт функциональных языков программированияЭксперт Python
37302 / 20736 / 4272
Регистрация: 12.02.2012
Сообщений: 34,127
Записей в блоге: 14
24.03.2016, 16:08 2
Основой метода Крамера является вычисление определителей. Вот один из способов:

Lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
(defun mksig (n &optional (res '(1)))
  (cond ((= n 1) res)
        (t (mksig (- n 1) (append res (list (- (car (last res)))))))))
 
==> mksig
 
(defun minor (matr n m)
  (let ((tmp (del-el matr n)))
       (mapcar #'(lambda (x) (del-el x m)) tmp)))
 
==> minor
 
(defun del-el (lst n)
  (cond ((= n 1) (cdr lst))
        (t (cons (car lst) (del-el (cdr lst) (- n 1)))))) 
 
==> del-el
 
(defun det (matr)
 (let ((n (length matr)))
  (cond ((= 2 n) (- (* (caar matr) (cadr (cadr matr))) (* (cadar matr) (caadr matr))))
        (t (let ((s (car matr)) (z (mksig n)) (w (range 1 n)))
            (apply '+ (mapcar #'(lambda (x y p) (* x y (det (minor matr 1 p)))) s z w)))))))
            
==> det
 
(det '((1 3  3  1)
       (1 4  6  4)
       (1 5 10 10)
       (1 6 15 20)))
 
==> 1
 
(det '((10.96597 35.10765 96.72356)
       (2.35765 -84.11256 0.87932)
       (18.24689 22.13579 1.11123)))
 
==> 152731.360015724
Добавлено через 6 минут
А вот метод Гаусса:

Lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
(defun SLE! (matr rs &optional (n (length matr)) (k 0) (z 1))
  (cond ((= k (- n 1))
         (rev-walk matr rs))
        (t (let* ((umatr (if (zerop k) nil (subseq matr 0 k)))
                  (urs   (if (zerop k) nil (subseq rs 0 k)))
                  (row (nth k matr))
                  (flgz (chk-zero row))
                  (f (if flgz 1 (+ 0.0 (nth k (nth k matr)))))
                  (row-k (if flgz row (mapcar #'(lambda (x) (/ x f)) row))) 
                  (dmatr (list row-k))
                  (drs (list (/ (nth k rs) f))))
                 (if flgz 
                     0
                     (progn                      
                        (dolist (i (range! (+ k 1) (- n 1)) dmatr)
                          (let ((v (get-elt matr i k)) (row-i (nth i matr)))
                            (push (mapcar #'(lambda (x y) (+ x (* y (- v)))) row-i row-k)
                                  dmatr)
                            (push (+ (nth i rs) (/ (* (- v) (nth k rs)) f)) drs)))   
                        (SLE! (append umatr (reverse dmatr)) (append urs (reverse drs))  n (+ k 1) (* z f))))))))                        
 
==> SLE!
 
(defun rev-walk (matr y)
  (let* ((n (length y)) 
         (r nil)
         (x (replicate 1.0 n)))
        (dotimes (j n t)
          (let* ((i (- n j 1))
                 (row (nth i matr))
                 (s   (apply '+ (mapcar #'(lambda (a x) (* a x))  (cdr (subseq row i)) (cdr (subseq x i))))))
                      (setf (nth i x) (/ (- (nth i y) s) (nth i row))))) x ))  
 
==> rev-walk
 
(defun chk-zero (row) 
  (reduce #'(lambda (ac x) (and ac (<= (abs x) 1.0E-11))) row :initial-value t)) 
 
==> chk-zero
 
 
(defun range! (n m)
  (cond ((= n m) (list n))
        ((> n m) (range! m n))
        (t (let ((res (list m)))
                 (dotimes (i (- m n) res) (push (- (car res) 1) res))))))
 
==> range!
 
(defun get-elt (matr r c)
  (nth c (nth r matr)))  
  
==> get-elt
 
(sle! '((1 1 1 1) (1 2 2 2) (1 2 3 3) (1 2 3 4)) '(4 6 7 8))
 
==> (2.0 1.0 0.0 1.0)
3
4527 / 3521 / 358
Регистрация: 12.03.2013
Сообщений: 6,038
24.03.2016, 23:44 3
Лучший ответ Сообщение было отмечено Catstail как решение

Решение

Лисп лиспом, но линейную алгебру надо знать!

Десять уравнений ― это маленькие системы. В численных методах счёт уравнений может идти на десятки и сотни тысяч. Но это отдельная наука, я плохо её знаю. Тем не менее, по Крамеру всё равно не хочется решать ничего выше второго порядка. Нерациональный он.

Вообще, у систем может быть разное поведение. Они могут не иметь решения. Если у них есть решения, то либо ровно одно, либо бесконечно много. Во втором случае все решения можно получить, выразив так называемые главные неизвестные через свободные, и давая свободным произвольные значения. Какое API вы предполагаете для решения этих вопросов, трудно сказать. Мне кажется, имеет смысл написать две функции: приведение матрицы к ступенчатому виду и к «приведённому» ступенчатому виду (не припомню, как он по-русски, а по-английски ― reduced row echelon form). В частности, по приведённому ступенчатому виду расширенной матрицы системы невооружённым глазом видно решение системы. Неприведённый ступенчатый вид матрицы можно использовать и для других целей: вычисление ранга, например, или определителя, если при получении этого вида не используется умножение строк на число.

Начнём программу на лиспе так:
Lisp
1
2
3
4
5
6
7
8
9
10
11
(defpackage #:gauss
  (:use #:cl)
  (:export #:row-echelon
           #:reduced-row-echelon
           #:row-dimension
           #:column-dimension
           #:switch-rows
           #:multiply-row
           #:add-row))
 
(in-package #:gauss)
Пакет ― это навроде пространства имён. Когда вы открываете лисп, пакетом по умолчанию является CL-USER. Однако когда вы пишете что-то, вы обычно определяете свои пакеты и работаете в них. В данном случае мы определяем пакет GAUSS, импортируем в него пакет CL (это никнейм для CL-USER) и экспортируем несколько символов.

Вообще, чтобы получить доступ к символу mysym из пакета mypackage, можно использовать префиксную запись mypackage::mysym. Пакет может экспортировать некоторые символы. Такие символы называются внешними, и в префиксной записи можно ставить одно двоеточие: mypackage:my-external-sym. Внешние символы можно рассматривать как API пакета. Когда вы видите в коде префикс с двойным двоеточием, это всегда подозрительно: почему пользуются неэспортированным символом? Не деталь ли это реализации, которая может измениться в следующем релизе?

Если один пакет импортирует другой, то все внешние символы импортируемого пакета можно использовать без префиксов. Это простая, но удобная система.

Мы импортируем пакет CL, чтобы можно было без префиксов использовать обычные лисповые функции и т. п. Экспортируем мы имена функций, представляющие ценность для внешнего мира. Это row-echelon и reduced-row-echelon, вычисляющие соответственно ступенчатую и приведённую ступенчатую формы, а также некоторые утилиты: число строк и столбцов матрицы и элементарные операции со столбцами.

Потом мы пишем форму in-package, которая делает пакет GAUSS текущим.

Обычно все формы defpackage собирают в отдельном файле, по традиции называемом package.lisp, а файлы с кодом начинают сразу строчкой in-package. Однако если вас интересует организация библиотек в виде файлов, лучше в отдельной теме обсудить. Пока и так сойдёт.

Теперь переходим к коду.

В лиспе есть многомерные массивы, поэтому, не мудрствуя, будем работать с двумерными массивами.

Размеры многомерного массива можно достать функцией array-dimension. Раз мы работаем с матрицами, удобно завести функции row-dimension и column-dimension, возвращающие число строк и столбцов. Функции маленькие, поэтому заинлайним их (предварительная оптимизация, так сказать).

Lisp
1
2
3
4
5
6
7
8
9
(declaim (inline row-dimension column-dimension))
 
(defun row-dimension (a)
  "Return the number of rows of a matrix A."
  (array-dimension a 0))
 
(defun column-dimension (a)
  "Return the number of columns of a matrix A."
  (array-dimension a 1))
Теперь ― метод Гаусса.

Мы имеем право применять три вида элементарных преобразований строк: перестановка двух строк, умножение строки на число, и прибавление к строке другой строки, умноженной на число. Нетрудно понять, что соответствующие операции для систем уравнений не влияют на множество решений (если только во втором случае не умножать на 0).

Напишем функции, которые выполняют эти преобразования.

Поскольку создавать на каждый чих новую матрицу накладно, мы будем преобразовывать матрицы in place, то есть изменять их сами. Чтобы не забывать, что функции изменяют свои аргументы (что для функций, вообще говоря, нетипично), мы в докстринг включим слово destructively.

При этом функциям неплохо и возвращать что-нибудь. Договоримся, что они будут возвращать саму матрицу.
Lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
(defun switch-rows (a i j)
  "Destructively swap rows I and J of a matrix A, return A."
  (dotimes (k (column-dimension a) a)
    (psetf (aref a i k) (aref a j k)
           (aref a j k) (aref a i k))))
 
(defun multiply-row (a i alpha)
  "Destructively multiply the Ith row of A by alpha, return A."
  (dotimes (k (column-dimension a) a)
    (setf (aref a i k) (* (aref a i k) alpha))))
 
(defun add-row (a i j alpha)
  "Destructively add to Ith row of A its Jth row multiplied by alpha, return A."
  (dotimes (k (column-dimension a) a)
    (incf (aref a i k) (* (aref a j k) alpha))))
Реализация очевидна: счётчик цикла проходит от 0 до числа числа столбцов минус 1 (нумерация массивов в лиспе всегда с нуля) и выполняются соответствующие присваивания.

dotimes ― один из немудрёный операторов цикла в лиспе, который запускает счётчик всегда с нуля и перебирает указанное число значений счётчика. То есть dotimes i n эквивалентно for (i = 0; i < n; ++i). В лиспе есть обычай отовсюду возвращать значения, это очень удобно. По умолчанию dotimes возвращает неинтересное значение nil, но факультативный третий аргумент в шапке является выражением, которое нужно вернуть. Мы указываем вернуть a, то есть саму матрицу.

Вообще, функция возвращает значение последней вычисленной формы. В наших функциях по одной форме ― dotimes; их значение и будет возвращено.

Доступ к элементам массива осуществляется функцией aref: aref a i j эквивалентно a[i, j].

Когда мы меняем местами два элемента, мы можем использовать временную переменную для хранения одной из них. Вместо этого мы воспользуемся оператором параллельного присваивания psetf. Обычный оператор присваивания ― setf.

incf эквивалентно +=. К сожалению, аналога *= нет, поэтому в multiply-row честно пишем выражение полностью.

Каков алгоритм приведения матрицы к ступенчатому виду?

Берём первый нулевой столбец. Если в нём одни нули, переходим к следующему столбцу. Если в нём есть ненулевой элемент, то его можно считать первым нулевым по номеру, чего всегда можно добиться перестановкой строк. Тогда мы зануляем все элементы столбца № 0, стоящие ниже ненулевого, путём прибавления к последующим строкам строки 0, умноженной на подходящие числа. Занулив столбец, двигаемся дальше, зануляя следующие столбцы.

Я не в ударе, почитайте любой учебник линейной алгебры.

Напишем вспомогательные функции.

Сначала ― занулить все элементы под a[i,j]. Здесь нужно делать цикл по номерам строк от i+1 до числа строк. К сожалению, dotimes умеет только с нуля. Поэтому мы воспользуемся мудрёным макросом ― loop с ключевыми словами. Он представляет большие возможности ― например, больше, чем сишный for. На самом деле это специальный язык для циклов ― язык внутри языка, который надо учить отдельно. Он не без недостатков и не всем нравится, но часто удобен. Здесь мы используем его наподобие паскалевского for; однако, below означает, что счётчик принимает только значения, строго меньшие верхнего предела.
Lisp
1
2
3
4
5
(defun eliminate-column-below (a i j)
  "Assuming that a[i,j] is nonzero, destructively eliminate nonzero entries below a[i,j]; return a."
  (loop for k from (+ i 1) below (array-dimension a 0)
        do (add-row a k i (- (/ (aref a k j) (aref a i j))))
        finally (return a)))
Клауза finally выполняется в конце итерации; в нашем случае она выполняет форму (return a), возвращающую из цикла значение a. Так как цикл ― единственная форма в функции, его значение становится значением функции.

Заодно напишем зануление элементов над a[i,j], пригодится для обратного хода метода Гаусса. Его можно было бы сделать с помощью dotimes, но я просто подправил код предыдущей функции. В for нижний предел 0 можно не писать.
Lisp
1
2
3
4
5
(defun eliminate-column-above (a i j)
  "Assuming that a[i,j] is nonzero, destructively eliminate nonzero entries above a[i,j]; return a."
  (loop for k below i
        do (add-row a k i (- (/ (aref a k j) (aref a i j))))
        finally (return a)))
Теперь нам надо уметь искать ненулевые элементы в столбце j, начиная со строки i. Если все нулевые, вернём ложь ― nil. Значение nil удобно возвращать в случае «неудачи» ― конечно, если оно не может появиться в «удачном» случае.
Lisp
1
2
3
4
5
6
(defun find-pivot-row (a i j)
  "Return the first row number starting from i having a nonzero entry in the jth column, or nil if it does not exist."
  (loop for k from i below (row-dimension a)
        unless (zerop (aref a k j))
          do (return k)
        finally (return nil)))
unless ― одна из условных клауз цикла loop. В данном случае мы перебираем номера строки, если находим ненулевой (zerop проверяет, нулевой ли у неё аргумент), возвращаем номер строки из цикла; если же мы не выйдем из цикла досрочно, то нулевых элементов не было, и возвращаем nil. Вообще, насчёт loop почитайте «Практический Common Lisp» или ещё что-то; поначалу можете использовать более ортодоксальные средства.

Сестра-близнец пригодится для обратного хода. (Вообще, мне начинает не нравиться дублирование кода.)
Lisp
1
2
3
4
5
6
(defun find-pivot-column (a i j)
  "Return the first column number starting from i having a nonzero entry in the ith row, or nil if it does not exist."
  (loop for k from j below (column-dimension a)
        unless (zerop (aref a i k))
          do (return k)
        finally (return nil)))
Теперь можно написать приведение к ступенчатому виду. Опять не на лиспе, а на «лупе». А что делать, задача простая, продвинутые средства для неё не нужны, достаточно массивов и циклов. Мы не будем использовать умножение строк на числа, чтобы функцию можно было использовать для вычисления определителей.
Lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
(defun row-echelon (a)
  "Compute the row echelon form in-place without multiplying rows by numbers; return a."
  (loop with row-dimension = (row-dimension a)
        with column-dimension = (column-dimension a)
        with current-row = 0
        with current-col = 0
        while (and (< current-row row-dimension)
                   (< current-col column-dimension))
        for pivot-row = (find-pivot-row a current-row current-col)
        do (when pivot-row
             (unless (= pivot-row current-row)
               (switch-rows a pivot-row current-row))
             (eliminate-column-below a current-row current-col)
             (incf current-row))
        do (incf current-col)
        finally (return a)))
Реализовано примерно то, что я пытался написать выше. Мы начинаем со строки и столбца 0. Если в текущем столбце нет ненулевых элементов, увеличиваем номер текущего столбца. Если есть ненулевой элемент, но не в текущей строке, меняем строки; убедившись, что в текущих строке и столбце не 0, зануляем столбец под этим элементом и увеличиваем на 1 текущую строку и столбец.

Мы используем следующие возможности цикла loop. Во-первых, with ... = вводит локальные переменные для цикла. Во-вторых, мы совмещаем одновременно делаем цикл типа while, и цикл по номерам строк с ненулевыми элементами. Дело в том, что клауза for может иметь вид типа for i from 1 to 10, а может иметь вид типа for i = выражение, в последнем случае на каждой итерации переменной i присваивается значение выражения. Тело цикла у нас разбито в два do. Сначала ― что делать, если нашли-таки строку с ненулевым элементом (я воспользовался лисповым оператором while, эквивалентным if с одной веткой, а можно было воспользоваться одноимённой клаузой loop); а во втором do ― увеличение текущего столбца, которое нужно делать при любом раскладе. Да, incf по умолчанию увеличивает на 1.

Следующая функция выполняет обратный ход метода Гаусса в предположении, что аргумент уже в ступенчатом виде. Алгоритм простой: в каждой строке ищется ненулевой элемент и зануляются все элементы над ним. Здесь мы совмещаем два счётчика и цикл while (если в строке уже нет ненулевых элементов, нужно заканчивать).
Lisp
1
2
3
4
5
6
7
8
9
(defun reduce-row-echelon (a)
  "Assuming that a has the row echelon form, compute the reduced row echelon form in-place without multiplying rows by numbers; return a."
  (loop for i below (row-dimension a)
        for j = (find-pivot-column a i 0) then (find-pivot-column a i j)
        while j
        unless (= 1 (aref a i j))
        do (multiply-row a i (/ 1 (aref a i j)))
        do (eliminate-column-above a i j)
        finally (return a)))
Наконец, следующая функция возвращает приведённую ступенчатую форму. Она, в частности, и решает уравнения.
Lisp
1
2
3
(defun reduced-row-echelon (a)
  "Compute the reduced row echelon form in-place without multiplying rows by numbers; return a."
  (reduce-row-echelon (row-echelon a)))
Пример 567 из Проскурякова:
https://www.cyberforum.ru/cgi-bin/latex.cgi?<br />
3x_1 - 2x_2 - 5x_3 + x_4 = 3<br />
2x_1 - 3x_2 + x_3 + 5x_4 = -3<br />
x_1 + 2x_2 - 4 x_4 = -3<br />
x_1 - x_2 - 4x_3 + 9x_4 = 22<br />
Lisp
1
2
GAUSS> (reduced-row-echelon #2a((3 -2 -5 1 3) (2 -3 1 5 -3) (1 2 0 -4 -3) (1 -1 -4 9 22)))
#2A((1 0 0 0 -1) (0 1 0 0 3) (0 0 1 0 -2) (0 0 0 1 2))
Это значит, исходная система эквивалентна следующей:
https://www.cyberforum.ru/cgi-bin/latex.cgi?<br />
x_1 = -1<br />
x_2 = 3<br />
x_3 = -2<br />
x_4 = 2<br />
Нетрудно понять, какое решение.

Неопределённая система, № 690:
https://www.cyberforum.ru/cgi-bin/latex.cgi?<br />
2x_1 - 3x_2 + 5x_3 + 7x_4 = 1<br />
4x_1 - 6x^2 + 2x_3 + 3x_4 = 2<br />
2x_1 - 3x_2 - 11x_3 - 15x_4 = 1<br />
Lisp
1
2
GAUSS> (reduced-row-echelon #2A((2 -3 5 7 1) (4 -6 2 3 2) (2 -3 -11 -15 1)))
#2A((1 -3/2 0 1/16 1/2) (0 0 1 11/8 0) (0 0 0 0 0))
Значит, система эквивалентна такой:
https://www.cyberforum.ru/cgi-bin/latex.cgi?x_1 - 3/2 x_2 + 1/16 x_4 = 1/2<br />
x_3 + 11/8 x_4 = 0
Кто немного знает линейную алгебру, понимает, что второе и четвёртое неизвестные ― свободные, они могут принимать любые значения, и общее решение имеет ви
https://www.cyberforum.ru/cgi-bin/latex.cgi?<br />
x_1 = 1/2 + 3/2c_1 - 1/16 c_2<br />
x_2 = c_1<br />
x_3 = - 11/8 c_2<br />
x_4 = c_2<br />
Из приведённой формы всё видно невооружённым глазом.

У Проскурякова в ответе в качестве главных взяты x1 и x2. С его ответом можно сравнить, если эти столбцы записать последними в матрице системы.

Несовместная система:
https://www.cyberforum.ru/cgi-bin/latex.cgi?<br />
x_1 + x_2 = 1<br />
x_1 + x_2 = 2<br />
Lisp
1
2
GAUSS>(reduced-row-echelon #2a((1 1 1) (1 1 2)))
#2A((1 1 0) (0 0 1))
Исходная система эквивалентна такой:
https://www.cyberforum.ru/cgi-bin/latex.cgi?<br />
x_1 + x_2 = 0<br />
0 = 1<br />
несовместность которой очевидна.

Учите алгебру!
2
4527 / 3521 / 358
Регистрация: 12.03.2013
Сообщений: 6,038
25.03.2016, 15:15 4
Весь код из предыдущего сообщения:
Lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
(defpackage #:gauss
  (:use #:cl)
  (:export #:row-echelon
           #:reduced-row-echelon
           #:row-dimension
           #:column-dimension
           #:switch-rows
           #:multiply-row
           #:add-row))
 
(in-package #:gauss)
 
(declaim (inline row-dimension column-dimension))
 
(defun row-dimension (a)
  "Return the number of rows of a matrix A."
  (array-dimension a 0))
 
(defun column-dimension (a)
  "Return the number of columns of a matrix A."
  (array-dimension a 1))
 
(defun switch-rows (a i j)
  "Destructively swap rows I and J of a matrix A, return A."
  (dotimes (k (column-dimension a) a)
    (psetf (aref a i k) (aref a j k)
           (aref a j k) (aref a i k))))
 
(defun multiply-row (a i alpha)
  "Destructively multiply the Ith row of A by alpha, return A."
  (dotimes (k (column-dimension a) a)
    (setf (aref a i k) (* (aref a i k) alpha))))
 
(defun add-row (a i j alpha)
  "Destructively add to Ith row of A its Jth row multiplied by alpha, return A."
  (dotimes (k (column-dimension a) a)
    (incf (aref a i k) (* (aref a j k) alpha))))
 
(defun eliminate-column-below (a i j)
  "Assuming that a[i,j] is nonzero, eliminate nonzero entries below a[i,j]; return a."
  (loop for k from (+ i 1) below (array-dimension a 0)
        do (add-row a k i (- (/ (aref a k j) (aref a i j))))
        finally (return a))) 
 
(defun eliminate-column-above (a i j)
  "Assuming that a[i,j] is nonzero, destructively eliminate nonzero entries above a[i,j]; return a."
  (loop for k below i
        do (add-row a k i (- (/ (aref a k j) (aref a i j))))
        finally (return a))) 
 
(defun find-pivot-row (a i j)
  "Return the first row number starting from i having a nonzero entry in the jth column, or nil if it does not exist."
  (loop for k from i below (row-dimension a)
        unless (zerop (aref a k j))
          do (return k)
        finally (return nil))) 
 
(defun find-pivot-column (a i j)
  "Return the first column number starting from i having a nonzero entry in the ith row, or nil if it does not exist."
  (loop for k from j below (column-dimension a)
        unless (zerop (aref a i k))
          do (return k)
        finally (return nil))) 
 
(defun row-echelon (a)
  "Compute the row echelon form in-place without multiplying rows by numbers; return a."
  (loop with row-dimension = (row-dimension a)
        with column-dimension = (column-dimension a)
        with current-row = 0
        with current-col = 0
        while (and (< current-row row-dimension)
                   (< current-col column-dimension))
        for pivot-row = (find-pivot-row a current-row current-col)
        do (when pivot-row
             (unless (= pivot-row current-row)
               (switch-rows a pivot-row current-row))
             (eliminate-column-below a current-row current-col)
             (incf current-row))
        do (incf current-col)
        finally (return a))) 
 
(defun reduce-row-echelon (a)
  "Assuming that a has the row echelon form, compute the reduced row echelon form in-place without multiplying rows by numbers; return a."
  (loop for i below (row-dimension a)
        for j = (find-pivot-column a i 0) then (find-pivot-column a i j)
        while j
        unless (= 1 (aref a i j))
        do (multiply-row a i (/ 1 (aref a i j)))
        do (eliminate-column-above a i j)
        finally (return a)))
 
(defun reduced-row-echelon (a)
  "Compute the reduced row echelon form in-place; return a."
  (reduce-row-echelon (row-echelon a)))
Добавлено через 15 часов 30 минут
В последнем докстринге надо убрать, конечно, without multiplying. Здесь пофикшено.
2
25.03.2016, 15:15
IT_Exp
Эксперт
87844 / 49110 / 22898
Регистрация: 17.06.2006
Сообщений: 92,604
25.03.2016, 15:15
Помогаю со студенческими работами здесь

Переделать программу, использующую метод Гаусса в метод Барейса
Всем ДВС! подскажите пожалуйста, как переделать эту программу использующую метод Гаусса в метод...

Базис или нет (метод Гаусса)
Здравствуйте! помогите пожалуста: Я сделала прогу решение системы уравнений методом Гаусса. Она...

Метод Гаусса-Зейделя,метод Якоби, LU разложения
есть у кого то примеры решения???в Wolfram Mathematica

Метод итераций Якоби и метод Гаусса-Зейделя
Подскажите что-нибудь,вот у мну лабораторная: РЕШЕНИЕ СИСТЕМ ЛИНЕЙНЫХ И НЕЛИНЕЙНЫХ УРАВНЕНИЙ 1....


Искать еще темы с ответами

Или воспользуйтесь поиском по форуму:
4
Ответ Создать тему
КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin
Copyright ©2000 - 2024, CyberForum.ru