PDA

Просмотр полной версии : Прога на Delphi.C++ или Pascal


Rastamanka
02.02.2009, 00:31
В общем ищу софт с исходниками на данных языках.
Тема софта:
1.Решение нелинейных уравнений методом итерации
2.Решение простых линейных уравнений методом итерации

blednii
03.02.2009, 19:41
Тебе надо чтоб решение выводил или только твет?

Rastamanka
04.02.2009, 11:27
Только ответ.

Flame of Soul
04.02.2009, 13:43
Решение системы линейных алгебраических уравнений ( СЛАУ ) методом простых итераций (методом Якоби).

Входные данные - матрица коэффициентов СЛАУ, вектор правых частей, размерность системы.
Выходные данные - вектор решения.

Вспомогательная функция getDiagonal( double **coefficients, double *rightPart, int currRowAndColumn, int numberOfEquation ) преобразует матрицу СЛАУ к диагональному виду с паралелльным преобразованием вектора правых частей.

Не всякая система может быть приведена к системе с матрицей с ненулевой диагональю. Если такое приведение не найдено ( в данном случае оно ищется перестановкой строк ) – система не решается.

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

Исходный код.
#include "stdio.h"
#include "iostream.h"
#include "math.h"

// приведение матрицы коэффициентов к виду с ненулевой диагональю и соответствующее изменение вектора правых частей
// возвращает true - если приведение - успешно
bool getDiagonal( double **coefficients, double *rightPart, int currColumn, int numberOfEquation ) {
bool result = false;
int i, row = currColumn;
double tempItem;

// для матрицы 1x1 ответом является ненулевость ее единственного элемента
if ( currColumn == numberOfEquation - 1 ) {
result = coefficients[currColumn][currColumn] != 0;
}
else {
// пока не найдена перестановка приводящая матрицу к виду с ненулевой диагональю и пока мы можем еще просматривать строки ищем перестановку
while ( !result && row < numberOfEquation ) {
// если текущий элемент на диагонали нулевой ищем в столбце дальше ненулевой
if ( coefficients[row][currColumn] != 0 ) {
// если элемент ненулевой и не лежит на диагонали меняем местами сотвествующие строки в матрице и элементы в векторе прваых частей
if ( row != currColumn ) {
for ( i = 0; i < numberOfEquation; i++ ) {
tempItem = coefficients[currColumn][i];
coefficients[currColumn][i] = coefficients[row][i];
coefficients[row][i] = tempItem;
}
tempItem = rightPart[currColumn];
rightPart[currColumn] = rightPart[row];
rightPart[row] = tempItem;
}
// рекурсивный вызов фактически для подматрицы с размерностью на 1 меньше
result = getDiagonal( coefficients, rightPart, currColumn + 1, numberOfEquation );
if ( result ) {
break;
}
}
row ++;
}
}

return result;
}

// было ли найдено решение, если да - итог в параметре solution
int iteration( double **coefficients, double *rightPart, int numberOfEquation, double *solution, double precision ) {
bool result;
int i, j, k, step = 1;
double temp;
double* tempSolution;

tempSolution = new double[numberOfEquation];

// приведение матрицы коэффициентов к виду с ненулевой диагональю и соответствующее изменение вектора правых частей
result = getDiagonal( coefficients, rightPart, 0, numberOfEquation );

// если приведение успешно - работаем дальше
if ( result ) {
double fault = precision + 1;

// преобразуем матрицу и правую часть для дальнейшего решения
for ( i = 0; i < numberOfEquation; i ++ ) {
for ( j = 0; j < numberOfEquation; j ++ ) {
if ( i != j ) {
coefficients[i][j] = - coefficients[i][j] / coefficients[i][i];
}
}
rightPart[i] = rightPart[i] / coefficients[i][i];
coefficients[i][i] = 0;
}

// первое приближение решения - преобразованный вектор правых частей
for ( i = 0; i < numberOfEquation; i ++ ) {
solution[i] = rightPart[i];
}

// пока не найдено решение с допустимй погрешнстью или пока не исчерпан лимит шагов... если все расходится например
while ( fault > precision && step <= 1000 ) {

// поиск новой итерации с ее "самоиспользованием" при расчетах
for ( j = 0; j < numberOfEquation; j ++ ) {
tempSolution[j] = 0;
}
for ( i = 0; i < numberOfEquation; i ++ ) {
for ( j = 0; j < numberOfEquation; j ++ ) {
tempSolution[i] = tempSolution[i] + coefficients[i][j] * solution[j];
}
tempSolution[i] = tempSolution[i] + rightPart[i];
}

// расчет погрешности
fault = 0.0;
for ( j = 0; j < numberOfEquation; j ++ ) {
fault = fault + (solution[j] - tempSolution[j])*(solution[j] - tempSolution[j]);
}
fault = sqrt( fault );

// сохранение полученной новой итерации
for ( j = 0; j < numberOfEquation; j ++ ) {
solution[j] = tempSolution[j];
}
step++;
}
}
else {
result = 1001;
}


return step;
}


void main() {
int i, j;
int size;
double **coefficients, *rightPart, *solution, precision;

cout << "Simple iteration method.\nEnter system dimension: ";
cin >> size;

coefficients = new double*[size];
for ( i = 0; i < size; i++ ) {
coefficients[i] = new double[size];
}
rightPart = new double[size];
solution = new double[size];

for ( i = 0; i < size; i ++ ){
cout << "Enter " << i + 1 << " row: ";
for ( j = 0; j < size; j ++ ){
cin >> coefficients[i][j];
}
}

cout << "Enter right part: ";
for ( j = 0; j < size; j ++ ){
cin >> rightPart[j];
}

cout << "Enter precision: ";
cin >> precision;

int steps = iteration( coefficients, rightPart, size, solution, precision );
if ( steps > 1000 ) {
cout << "Solution for this matrix of coefficients not exist or not found";
}
else {
cout << "Solution is:\n";
for ( j = 0; j < size; j ++ ){
cout << solution[j] << "\n";
}
cout << "Number of step: " << steps;
}

cout << "\nPress \"Enter\" to continue..." << endl;
getchar();
}

Примеры использования (распечатка листинга):

Simple iteration method.
Enter system dimension: 3
Enter 1 row: 10 1 1
Enter 2 row: 2 10 1
Enter 3 row: 2 2 10
Enter right part: 12 13 14
Enter precision: 0.1
Solution is:
0.9946
0.9934
0.9916
Number of step: 4
Press "Enter" to continue...

Flame of Soul
04.02.2009, 13:54
Метод простых итераций довольно медленно сходится. Для его ускорения существует метод Зейделя

Входные данные - матрица коэффициентов СЛАУ, вектор правых частей, размерность системы, заданная точность.
Выходные данные - вектор решения.

Вспомогательная функция getDiagonal( double **coefficients, double *rightPart, int currColumn, int numberOfEquation ) преобразует матрицу СЛАУ к виду с ненулевыми диагональными элементами с параллельным преобразованием вектора правых частей.

Не всякая система может быть приведена к системе с матрицей с ненулевой диагональю. Если такое приведение не найдено ( в данном случае оно ищется перестановкой строк ) – система не решается.
Быстрее метода простых итераций.

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

Исходный код:
#include "stdio.h"
#include "iostream.h"
#include "math.h"

// приведение матрицы коэффициентов к виду с ненулевой диагональю и соответствующее изменение вектора правых частей
// возвращает true - если приведение - успешно
bool getDiagonal( double **coefficients, double *rightPart, int currColumn, int numberOfEquation ) {
bool result = false;
int i, row = currColumn;
double tempItem;

// для матрицы 1x1 ответом является ненулевость ее единственного элемента
if ( currColumn == numberOfEquation - 1 ) {
result = coefficients[currColumn][currColumn] != 0;
}
else {
// пока не найдена перестановка приводящая матрицу к виду с ненулевой диагональю и пока мы можем еще просматривать строки ищем перестановку
while ( !result && row < numberOfEquation ) {
// если текущий элемент на диагонали нулевой ищем в столбце дальше ненулевой
if ( coefficients[row][currColumn] != 0 ) {
// если элемент ненулевой и не лежит на диагонали меняем местами сотвествующие строки в матрице и элементы в векторе прваых частей
if ( row != currColumn ) {
for ( i = 0; i < numberOfEquation; i++ ) {
tempItem = coefficients[currColumn][i];
coefficients[currColumn][i] = coefficients[row][i];
coefficients[row][i] = tempItem;
}
tempItem = rightPart[currColumn];
rightPart[currColumn] = rightPart[row];
rightPart[row] = tempItem;
}
// рекурсивный вызов фактически для подматрицы с размерностью на 1 меньше
result = getDiagonal( coefficients, rightPart, currColumn + 1, numberOfEquation );
if (result) {
break;
}
}
row ++;
}
}

return result;
}

// было ли найдено решение, если да - итог в параметре solution
int iteration( double **coefficients, double *rightPart, int numberOfEquation, double *solution, double precision ) {
bool result;
int i, j, k, step = 1;
double temp;
double* tempSolution;

tempSolution = new double[numberOfEquation];

// приведение матрицы коэффициентов к виду с ненулевой диагональю и соответствующее изменение вектора правых частей
result = getDiagonal( coefficients, rightPart, 0, numberOfEquation );

// если приведение успешно - работаем дальше
if ( result ) {
double fault = precision + 1;

// преобразуем матрицу и правую часть для дальнейшего решения
for ( i = 0; i < numberOfEquation; i ++ ) {
for ( j = 0; j < numberOfEquation; j ++ ) {
if ( i != j ) {
coefficients[i][j] = - coefficients[i][j] / coefficients[i][i];
}
}
rightPart[i] = rightPart[i] / coefficients[i][i];
coefficients[i][i] = 0;
}

// первое приближение решения - преобразованный вектор правых частей
for ( i = 0; i < numberOfEquation; i ++ ) {
solution[i] = rightPart[i];
}

// пока не найдено решение с допустимй погрешнстью или пока не исчерпан лимит шагов... если все расходится например
while ( fault > precision && step <= 1000 ) {

// поиск новой итерации
for ( j = 0; j < numberOfEquation; j ++ ) {
tempSolution[j] = 0;
}
for ( i = 0; i < numberOfEquation; i ++ ) {
for ( j = 0; j < i; j ++ ) {
tempSolution[i] = tempSolution[i] + coefficients[i][j] * tempSolution[j];
}
for ( j = i; j < numberOfEquation; j ++ ) {
tempSolution[i] = tempSolution[i] + coefficients[i][j] * solution[j];
}
tempSolution[i] = tempSolution[i] + rightPart[i];
}

// расчет погрешности
fault = 0;
for ( j = 0; j < numberOfEquation; j ++ ) {
fault = fault + (solution[j] - tempSolution[j])*(solution[j] - tempSolution[j]);
}
fault = sqrt( fault );

// сохранение полученной новой итерации
for ( j = 0; j < numberOfEquation; j ++ ) {
solution[j] = tempSolution[j];
}
step++;
}
}
else {
result = 1001;
}


return step;
}


void main() {
int i, j;
int size;
double **coefficients, *rightPart, *solution, precision;

cout << "Simple iteration method.\nEnter system dimension: ";
cin >> size;

coefficients = new double*;
for ( i = 0; i < size; i++ ) {
coefficients[i] = new double[size];
}
rightPart = new double[size];
solution = new double[size];

for ( i = 0; i < size; i ++ ){
cout << "Enter " << i + 1 << " row: ";
for ( j = 0; j < size; j ++ ){
cin >> coefficients[i][j];
}
}

cout << "Enter right part: ";
for ( j = 0; j < size; j ++ ){
cin >> rightPart[j];
}

cout << "Enter precision: ";
cin >> precision;

int steps = iteration( coefficients, rightPart, size, solution, precision );
if ( steps > 1000 ) {
cout << "Solution for this matrix of coefficients not exist or not found";
}
else {
cout << "Solution is:\n";
for ( j = 0; j < size; j ++ ){
cout << solution[j] << "\n";
}
cout << "Number of step: " << steps;
}

cout << "\nPress \"Enter\" to continue..." << endl;
getchar();
}
[SIZE=1]
Примеры использования (распечатка листинга):

Simple iteration method.
Enter system dimension: 3
Enter 1 row: 10 1 1
Enter 2 row: 2 10 1
Enter 3 row: 2 2 10
Enter right part: 12 13 14
Enter precision: 0.1
Solution is:
1.00068
0.997944
1.00028
Number of step: 3
Press "Enter" to continue...



PS: Больше не чем помочь немогу