IPB
ЛогинПароль:

> ВНИМАНИЕ!

Прежде чем задать вопрос, смотрите FAQ.
Рекомендуем загрузить DRKB.

Наладить общение поможет, если вы подпишитесь по почте на новые темы в этом форуме.

 
 Ответить  Открыть новую тему 
> LU разложение.
сообщение
Сообщение #1


Профи
****

Группа: Пользователи
Сообщений: 652
Пол: Мужской
Реальное имя: Алексей

Репутация: -  20  +


Собственно есть задача, полную ее суть писать не стану,но проблема заключается в следующем нужно с помощью LU разложения решить систему.Рабочий код на С++ у меня был и я его попытался перенести на Delphi.Но к несчастью моя затея обернулась неудачей.

type
TVec = array of Extended; // векторный тип
TMatr = array of TVec; // матричный тип

function LUrazl(matr: Tmatr): Tvec;
var
i, j,k, row, col: Word;
res: TVec;
_matr_: TMatr;
L,U: TMatr;
Sum: double;
begin
_matr_ := CopyMatr(matr); // Copy не работает как надо
row := Length(_matr_);
col := Length(_matr_[0]);
SetLength(L,row,col);
SetLength(U,row,col);
if row <> col-1 then
ShowMessage('Ошибка: некорректные размерности матрицы! (МГ)'); //error
for i:=0 to row-1 do
for j:= 0 to col-1 do
begin
L [i][j] := 0;
U [i][j] := 0;

if i = j then
L [i][j] := 1;
end;



//==============================================

//находим первый столбец L[][] и первую строку U[][]

for i:= 0 to row-1 do
begin
L [i][0]:= _matr_ [i][0]/ _matr_ [0][0];
U [0][i]:= _matr_ [0][i] / L [0][0];
end;

for i:= 0 to row-1 do
begin
for j:= 0 to col-1 do
begin
U [0][ i] := _matr_[0][ i];
L [i][ 0] := _matr_[i][ 0] / U[0][ 0];
sum := 0;
for k:= 0 to i do
begin
sum :=sum+ L[i][ k] * U[k][ j];
end;
U[i][ j] := _matr_[i][ j] - sum;
if (i > j) then
L[j][ i] := 0
else
begin
sum := 0;
for k:= 0 to i do
begin
sum:=sum+L[j][ k] * U[k][ i];
end;
L[j][ i] := (_matr_[j][ i] - sum) / U[i][ i];//В этой строке выдает деление на ноль так же в дебаге видно,что значение i=65000+
end;

end;
res[i]:= U[i][ i];
end;
end;


Непонятна причина такого поведения,ведь на С++ все это работало.Если надо могу добавить файл целиком.
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
сообщение
Сообщение #2


Гуру
*****

Группа: Пользователи
Сообщений: 1 013
Пол: Мужской
Ада: Разработчик
Embarcadero Delphi: Сторонник
Free Pascal: Разработчик

Репутация: -  627  +


Добавь лучше С++-ную реализацию. Посмотрим, что сделано не так, почему не работает...
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
сообщение
Сообщение #3


Профи
****

Группа: Пользователи
Сообщений: 652
Пол: Мужской
Реальное имя: Алексей

Репутация: -  20  +


Вот С++,так как оно мне требовалось, чтобы посчитать стопочку типовых расчетов,особого энтузиазма в оформлении нет.

#include <iostream>
using namespace std;

int main ()
{
const int n=4;
double sum = 0;
double A [ n ][ n ]= {
{ 2, -5,6,1},
{ -16, 48,-55,-5},
{ 0, 24,-25,10},
{8,21,30,4 }};
double L [ n ][ n ];
double U [ n ][ n ];

//задаем матрицу A[][] ...

for (int i = 0; i < n; i++)
{
for (int j = 0; j < n; j++)
{
// cout << "\na[" << i+1 << "][" << j+1 << "] = ";
// cin >> A [i][j];
L [i][j] = 0;
U [i][j] = 0;

if (i == j)
L [i][j] = 1;
}
}

//==============================================

//находим первый столбец L[][] и первую строку U[][]

for (int i = 0; i < n; i++)
{
L [i][0] = A [i][0]/ A [0][0];
U [0][i] = A [0][i] / L [0][0];
}

for (int i = 0; i < n; i++)
{
for (int j = 0; j < n; j++)
{
U [0][ i] = A[0][ i];
L [i][ 0] = A[i][ 0] / U[0][ 0];
double sum = 0;
for (int k = 0; k < i; k++)
{
sum += L[i][ k] * U[k][ j];
}
U[i][ j] = A[i][ j] - sum;
if (i > j)
{
L[j][ i] = 0;
}
else
{
sum = 0;
for (int k = 0; k < i; k++)
{
sum += L[j][ k] * U[k][ i];
}
L[j][ i] = (A[j][ i] - sum) / U[i][ i];
}
}
}


//====================================================
cout << "\n\n";
for (int i = 0; i < n; i++)
{
for (int j = 0; j < n; j++)
{
if (i == j)
L [i][j] = 1;
cout << " " << L [i][j] << " ";
}
cout << "\n\n";
}


cout << "\n\n";

for (int i = 0; i < n; i++)
{
for (int j = 0; j < n; j++)
cout << " " << U [i][j] << " ";
cout << "\n\n";
}

return 0;
}


Так как типовые были сдаты, и пара вариантов были посчитаны руками и результаты были одинаковыми,поэтому я счел,что он работает верно.

Да и забыл упомянуть для текущей задачи мне требуеться,чтобы оно считалось для максимальных размерностей 101х100.А то у меня уже был случай,когда Гаусс верно считал для любой матрицы меньше 20,а для больших переставал работать,пришлось полностью пересмотреть концепцию)))

Сообщение отредактировано: Krjuger -
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
сообщение
Сообщение #4


Гуру
*****

Группа: Пользователи
Сообщений: 1 013
Пол: Мужской
Ада: Разработчик
Embarcadero Delphi: Сторонник
Free Pascal: Разработчик

Репутация: -  627  +


У тебя как минимум в двух местах неправильно сделан перенос с С++... Оба - циклы по K. И там и там надо ходить до (i - 1), а не до i... Тестировал программу на той же матрице, что и Сишный вариант?

P.S.

Матрицы вот такие должны получаться (это результат работы твоего чуть-чуть подпиленного дельфийского кода):
L-matrix:
1.0000 0.0000 0.0000 0.0000
-8.0000 1.0000 0.0000 0.0000
0.0000 3.0000 1.0000 0.0000
4.0000 5.1250 -10.4688 1.0000
U-matrix:
2.0000 -5.0000 6.0000 1.0000
0.0000 8.0000 -7.0000 3.0000
0.0000 0.0000 -4.0000 1.0000
0.0000 0.0000 0.0000 -4.9063

?

Да, кстати, у тебя там балаган с возвращаемыми значениями, ты под Res память-то не выделил. Да и как собрался возвращать и L- и U-матрицы? Проще передать через параметры...

Сообщение отредактировано: IUnknown -
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
сообщение
Сообщение #5


Профи
****

Группа: Пользователи
Сообщений: 652
Пол: Мужской
Реальное имя: Алексей

Репутация: -  20  +


Нет не тестировал,я смысла не вижу,если оно даже не отрабатывает до конца,то говорить правильности решения пока рано.
Цитата
У тебя как минимум в двух местах неправильно сделан перенос с С++... Оба - циклы по K. И там и там надо ходить до (i - 1), а не до i...

К несчастью так нельзя,когда мы самый первый раз заходим в цикл по i у нас i=0, и i-1 в цикле по K мне даст большую беду,в Сях то такой цикл просто проскочился бы не отробатывая,а тут так не катит(((

Цитата

Да, кстати, у тебя там балаган с возвращаемыми значениями, ты под Res память-то не выделил. Да и как собрался возвращать и L- и U-матрицы? Проще передать через параметры...

Мне не нужно возвращать эти матрицы,мне нужно вернуть вектор, который состоит из элементов диагонали матрици U это и будут решения системы,такой уж метод.

Цитата
Матрицы вот такие должны получаться (это результат работы твоего чуть-чуть подпиленного дельфийского кода):

В Сях то он таким и получаеться.

Сообщение отредактировано: Krjuger -
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
сообщение
Сообщение #6


Гуру
*****

Группа: Пользователи
Сообщений: 1 013
Пол: Мужской
Ада: Разработчик
Embarcadero Delphi: Сторонник
Free Pascal: Разработчик

Репутация: -  627  +


Цитата
К несчастью так нельзя,когда мы самый первый раз заходим в цикл по i у нас i=0, и i-1 в цикле по K мне даст большую беду
Бла-бла-бла... Кроме Word-а в Дельфи есть еще и Integer, который способен принимать отрицательные значения. Чуешь, откуда ошибка?

В общем, вот, прогнал на FPC твой дельфийский код, ничего не вылетает, печатает то, что приведено выше:

type
TVec = array of Extended;
TMatr = array of TVec;

procedure LUrazl(_matr_: Tmatr; var L, U : TMatr);
var
i, j, row, col: Word;
k : integer;
res: TVec;
Sum: double;
begin
// _matr_ := CopyMatr(matr);
row := Length(_matr_);
col := Length(_matr_[0]);
SetLength(L,row,col);
SetLength(U,row,col);
if row <> col-1 then
writeln('error');

for i:=0 to row-1 do
for j:= 0 to col-1 do
begin
L [i][j] := 0;
U [i][j] := 0;
if i = j then L [i][j] := 1;
end;

//==============================================

for i:= 0 to row-1 do
begin
L [i][0]:= _matr_ [i][0]/ _matr_ [0][0];
U [0][i]:= _matr_ [0][i] / L [0][0];
end;

for i:= 0 to row-1 do
begin
for j:= 0 to col-1 do
begin
U [0][ i] := _matr_[0][ i];
L [i][ 0] := _matr_[i][ 0] / U[0][ 0];
sum := 0;
for k:= 0 to i-1 do
begin
sum :=sum+ L[i][ k] * U[k][ j];
end;
U[i][ j] := _matr_[i][ j] - sum;
if (i > j) then
L[j][ i] := 0
else
begin
sum := 0;
for k:= 0 to i-1 do
begin
sum:=sum+L[j][ k] * U[k][ i];
end;
L[j][ i] := (_matr_[j][ i] - sum) / U[i][ i];
end;
end;
end;
end;

var
A : tmatr;
L, U : tMatr;
i, j : integer;

begin
setlength(a, 4, 4);

a[0, 0] := 2; a[0, 1] := -5; a[0, 2] := 6; a[0, 3] := 1;
a[1, 0] := -16; a[1, 1] := 48; a[1, 2] := -55; a[1, 3] := -5;
a[2, 0] := 0; a[2, 1] := 24; a[2, 2] := -25; a[2, 3] := 10;
a[3, 0] := 8; a[3, 1] := 21; a[3, 2] := 30; a[3, 3] := 4;

LURazl(A, L, U);

writeln('L-matrix:');
for i := 0 to Length(L) - 1 do
begin
for j := 0 to Length(L[0]) - 1 do
write(L[i, j]:10 :4);
writeln;
end;
writeln('U-matrix:');
for i := 0 to Length(U) - 1 do
begin
for j := 0 to Length(U[0]) - 1 do
write(U[i, j]:10 :4);
writeln;
end;

end.
Дельфи запускать лень smile.gif

Сообщение отредактировано: IUnknown -
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
сообщение
Сообщение #7


Гость






Цитата(Krjuger @ 1.11.2011 17:52) *

Собственно есть задача, полную ее суть писать не стану,но проблема заключается в следующем нужно с помощью LU разложения решить систему.Рабочий код на С++ у меня был и я его попытался перенести на Delphi.Но к несчастью моя затея обернулась неудачей.

type
TVec = array of Extended; // векторный тип
TMatr = array of TVec; // матричный тип

function LUrazl(matr: Tmatr): Tvec;
var
i, j,k, row, col: Word;
res: TVec;
_matr_: TMatr;
L,U: TMatr;
Sum: double;
begin
_matr_ := CopyMatr(matr); // Copy не работает как надо
row := Length(_matr_);
col := Length(_matr_[0]);
SetLength(L,row,col);
SetLength(U,row,col);
if row <> col-1 then
ShowMessage('Ошибка: некорректные размерности матрицы! (МГ)'); //error
for i:=0 to row-1 do
for j:= 0 to col-1 do
begin
L [i][j] := 0;
U [i][j] := 0;

if i = j then
L [i][j] := 1;
end;
//==============================================

//находим первый столбец L[][] и первую строку U[][]

for i:= 0 to row-1 do
begin
L [i][0]:= _matr_ [i][0]/ _matr_ [0][0];
U [0][i]:= _matr_ [0][i] / L [0][0];
end;

for i:= 0 to row-1 do
begin
for j:= 0 to col-1 do
begin
U [0][ i] := _matr_[0][ i];
L [i][ 0] := _matr_[i][ 0] / U[0][ 0];
sum := 0;
for k:= 0 to i do
begin
sum :=sum+ L[i][ k] * U[k][ j];
end;
U[i][ j] := _matr_[i][ j] - sum;
if (i > j) then
L[j][ i] := 0
else
begin
sum := 0;
for k:= 0 to i do
begin
sum:=sum+L[j][ k] * U[k][ i];
end;
L[j][ i] := (_matr_[j][ i] - sum) / U[i][ i];//В этой строке выдает деление на ноль так же в дебаге видно,что значение i=65000+
end;

end;
res[i]:= U[i][ i];
end;
end;


Непонятна причина такого поведения,ведь на С++ все это работало.Если надо могу добавить файл целиком.




А ВЫ не могли бы сбросить вашу задачу пожалуйста, она мне очень нужна...(lovetanya.15@mail.ru) Спасибо
 К началу страницы 
+ Ответить 

 Ответить  Открыть новую тему 
1 чел. читают эту тему (гостей: 1, скрытых пользователей: 0)
Пользователей: 0

 





- Текстовая версия 29.03.2024 14:23
500Gb HDD, 6Gb RAM, 2 Cores, 7 EUR в месяц — такие хостинги правда бывают
Связь с администрацией: bu_gen в домене octagram.name