Вчера делал свою курсовую работу по методу конечных элементов и, в частности, нужно было решить систему линейных алгебраических уравнений (СЛАУ).
Так как число неизвестных было не велико (21), то я выбрал обычный метод Гаусса.
Описание алгоритма:
Алгоритм решения СЛАУ методом Гаусса подразделяется на два этапа.
1) На первом этапе осуществляется так называемый прямой ход, когда путём элементарных преобразований над строками систему приводят к ступенчатой или треугольной форме, либо устанавливают, что система несовместна. А именно, среди элементов первого столбца матрицы выбирают ненулевой, перемещают его на крайнее верхнее положение перестановкой строк и вычитают получившуюся после перестановки первую строку из остальных строк, домножив её на величину, равную отношению первого элемента каждой из этих строк к первому элементу первой строки, обнуляя тем самым столбец под ним. После того, как указанные преобразования были совершены, первую строку и первый столбец мысленно вычёркивают и продолжают пока не останется матрица нулевого размера. Если на какой-то из итераций среди элементов первого столбца не нашёлся ненулевой, то переходят к следующему столбцу и проделывают аналогичную операцию.
2) На втором этапе осуществляется так называемый обратный ход, суть которого заключается в том, чтобы выразить все получившиеся базисные переменные через небазисные и построить фундаментальную систему решений, либо, если все переменные являются базисными, то выразить в численном виде единственное решение системы линейных уравнений. Эта процедура начинается с последнего уравнения, из которого выражают соответствующую базисную переменную (а она там всего одна) и подставляют в предыдущие уравнения, и так далее, поднимаясь по «ступенькам» наверх. Каждой строчке соответствует ровно одна базисная переменная, поэтому на каждом шаге, кроме последнего (самого верхнего), ситуация в точности повторяет случай последней строки.
Алгоритмическая сложность:
Метод Гаусса требует порядка O(n3) действий.
Исходный код на Delphi:
Только кодкопировать в буфер обменапечать
type TMatrixE = array of array of Extended;
type TVectorE = array of Extended;
function SolveSLAE(A: TMatrixE; B: TVectorE): TVectorE;
var i,j,k,n:integer;
h:real;
begin
n:=Length(A)-1;
SetLength(Result, Length(A));
//Прямой ход — исключение переменных
for i:=0 to n-1 do
for j:=i+1 to n do
begin
A[j,i]:=-A[j,i]/A[i,i];
for k:=i+1 to n do
A[j,k]:=A[j,k]+A[j,i]*A[i,k];
B[j]:=B[j]+A[j,i]*B[i]
end;
Result[n]:=B[n]/A[n,n];
//Обратный ход — нахождение корней
for i:=n-1 downto 0 do
begin
h:=B[i];
for j:=i+1 to n do h:=h-Result[j]*A[i,j];
Result[i]:=h/A[i,i]
end;
end;
Где A – матрица коэффициентов СЛАУ, B – Вектор свободных членов.
На этом на сегодня всё.