Программы для решения систем линейных

Алгебраических уравнений

 

Файл GAUSS.CPP.

 

// Подключение математики:

#include <math.h>

// Математика подключена.

gauss(a,b,x,n)

float a[][nmax]; // Исходная матрица коэффицентов:a

float b[]; // Правые части

float x[]; // Результат (будет получен при решении)

int n; // Количество элементов в векторах и размерность матрицы

{

int i,j,k,imax;

int nm1=n-1; // Используется для удобства вместо n, так как

// элементы массивов отсчитываются с 0, а не с 1

int err=0; // Начальное значение кода ошибки =0 (т.е. ошибки нет),

// перед окончанием функции переменная err передается

// в качестве значения функции

float amax,t,mi,sum;

 

// Проверка корректности переменной n и константы nmax

// в случае противоречий выход из функции:

// endofproc -метка конца функции

 

if (n<2) { err=2; goto endofproc; } // В матрице <2 элементов

if (n>nmax) { err=3; goto endofproc; }

// n>реальной размерности =nmax

 

// Теперь все корректно.

 

// Далее следует...

// Цикл из n-1 шагов, на каждом шаге находится главный элемент и

// преобразовывается матрица.

// Главный элемент определяется только из очередного столбца,

// а не из всей матрицы.

// n-й шаг цикла не нужен так как остается один элемент a[n-1][n-1],

// который сам по себе главный (рассматривается после цикла).

 

for (k=0; k<=nm1-1; k++)

{

// Поиск главного элемента в текущем столбце (k):

 

for (i=k,amax=0.0; i<=nm1; i++)

if (fabs(a[i][k]) > fabs(amax)) amax=a[(imax=i)][k];

// Главный элемент найден.

// Если он =0, то ошибка, так как матрица вырождена и выход из функции

// (метка endofproc стоит в конце функции):

if (amax==0.0) { err=1; goto endofproc; }

// Замена строки с найденным главным элементом на строку с номером k

for (j=k;j<=nm1;j++)

{

t=a[imax][j];

a[imax][j]=a[k][j];

a[k][j]=t;

}

// Строки заменены.

// Замена элементов правой части для тех же строк:

t=b[imax];

b[imax]=b[k];

b[k]=t;

// Заменены.

// Определение mi для каждой строки (ниже главной)

// и ее преобразование для найденного mi:

for (i=k+1;i<=nm1;i++)

{

mi=-a[i][k]/amax;

a[i][k]=0.0;

for (j=k+1;j<=nm1;j++)

a[i][j]+=a[k][j]*mi;

b[i]+=b[k]*mi;

}

// На текущем шаге k матрица и правые части преобразованы.

} // Конец основного цикла (k увеличивается на единицу и все с начала)

// Проверка единственного оставшегося (и поэтому главного) элемента

// если он=0, то ошибка, так как матрица вырождена и выход из функции

if (a[nm1][nm1]==0.0) { err=1; goto endofproc; }

// Теперь матрица a и вектор b преобразованы (подготовлены к поиску x)

// Решение полученной системы (обратный ход):

for (i=nm1;i>=0;i--)

{

sum=0.0;

for (j=i+1;j<=nm1;j++)

sum+=a[i][j]*x[j];

X[i]=(b[i]-sum)/a[i][i];

}

// Решение получено.

endofproc: // Конец процедуры: метка для аварийного выхода

return(err); // Функция будет иметь значение ошибки err

}

 

 

Файл MITER.CPP.

 

// Подключение математики:

#include <math.h>

// Математика подключена.

// Функция получения мантиссы и порядка float числа:

mantpor(x,mant,por)

float *x,*mant;

double *por;

{

int st;

st=(int)floor(log10(*x));

*por=pow10(st);

*mant=*x/(*por);

}

// Конец функции mantpor.

// Функция, округляющая float число до n разрядов:

float round(x,n)

float *x;

int n;

{

float xabs,mant,dpown,xr;

double por;

int signx;

if (*x!=0.0)

{

xabs=fabs(*x);

signx=(*x>=0)? 1:-1;

mantpor(&xabs,&mant,&por);

dpown=pow10(n);

xr=signx*(floor(mant*dpown+0.5)/dpown)*por;

}

else xr=0.0;

return (xr);

}

// Конец функции round.

// Функция, решающая приведенную систему методом итераций:

miter(a,b,x,n,ntch,it)

float a[][nmax]; // Приведенная матрица коэффициентов.

float b[]; // Вектор свободных членов.

float x[]; // Вектор неизвестных.

int n, // Количество элементов в векторах и размерность матрицы.

ntch, // Количество знаков после запятой в мантиссе результата,

// которые должны остаться после округления.

*it; // Количество итераций, потребовавшихся для получения x.

{

float sum,xst[nmax];

int err=0; // Значение ошибки (сейчас=0, т.е. ее нет)

// в конце функции в качестве ее значения возвращается

// эта переменная.

int nm1=n-1; // Используется вместо n, так как

// массивы отсчитываются с 0, а не с 1

int ntchp1=ntch+1;

int i,j,tochnok;

// Проверка корректности переменной n и константы nmax

// в случае противоречий ошибка и выход из функции:

// endofproc - метка в конце функции

if (n<2) {err=2; goto endofproc;} // В матрице <2 элементов

if (n>nmax) {err=3; goto endofproc;} // n>реальной размерности =nmax

// Теперь все корректно.

// Проверка на корректность количества разрядов мантиссы для округления:

// Если n выходит за ограничения, то ошибка и выход из функции:

if (! ((0<=ntch)&&(ntch<=6))) {err=1; goto endofproc;}

// Проверка закончена.

// Проверка соответствия матрицы a условиям сходимости процесса:

// По строкам:

for (i=0,sum=0.0; (i<=nm1)&&(sum<1.0); i++)

for (j=0,sum=0.0; j<=nm1; j++) sum+=fabs(a[i][j]);

// По строкам проверено.

// Если по строкам не подходит, то по столбцам:

if (sum>=1)

{

// По столбцам:

for (j=0,sum=0.0; (j<=nm1)&&(sum<1.0); j++)

for (i=0,sum=0.0; i<=nm1;i++) sum+=fabs(a[i][j]);

// По столбцам проверено.

// Если и по столбцам не подходит, то ошибка и выход из функции:

if (sum>=1) {err=1; goto endofproc;}

}

// Проверка на сходимость закончена.

// Начальная установка массива x (в принципе может быть любая):

for (j=0;j<=nm1;j++) x[j]=b[j];

// Закончена.

*it=0; // Обнуляется счетчик итераций.

// Основной цикл: на каждом шаге находятся новые значения x,

// выход при достижении необходимой точности x:

do

{

(*it)++; // Увеличение счетчика итераций.

for (j=0;j<=nm1;j++) xst[j]=x[j]; // Запоминает старое x в xst

// Получает новые x:

for (i=0;i<=nm1;i++)

{

for (j=0,sum=0.0; j<=nm1; j++) sum+=a[i][j]*xst[j];

X[i]=b[i]+sum;

}

// Новые значения x получены.

// Проверяет отличие нового x от старого,

// т. е. точность до ntch+1 знаков мантиссы:

for (j=0,tochnok=1; (j<=nm1)&&tochnok; j++)

tochnok=(round(&x[j],ntchp1)==round(&xst[j],ntchp1));

// Точность проверена, если все нормально, то tochnok=1 иначе 0

} // Конец тела основного цикла.

while (!(tochnok)); // Выход, если достигнута необходимая точность.

// Округление найденного x до ntch знаков мантиссы:

for (j=0;j<=nm1;j++) x[j]=round(&x[j],ntch);

// Округлены.

endofproc: // Метка для аварийного выхода в случае ошибок:

return (err); // Код ошибки возвращается в качестве значения функции.

}

 

Список литературы

 

1. Амосов А. А., Дубинский Ю. А., Копченова Н. В. Вычислительные методы для инженеров: Учеб. пособие. - М.: Высш. шк., 1994.

2. Бахвалов Н. С., Жидков Н. П., Кобельков Г. М. Численные методы. - М.: Наука, 1987.

3. Березин И. С., Жидков Н. П. Методы вычислений: В 2-х т.- М.: Наука, 1966. Т.1.

4. Демидович Б. П., Марон И. А. Основы вычислительной математики. - М.: Наука, 1970.

5. Копченова Н. В., Марон И. А. Вычислительная математика в примерах и задачах. - М.: Наука, 1972.

6. Марчук Г. И. Методы вычислительной математики. - М.: Наука, 1989.

7. Плис А. И., Сливина Н. А. Лабораторный практикум по высшей математике: Учеб. пособие для втузов. - 2-е изд. - М.: Высш. шк., 1994.

8. Самарский А. А., Гулин А. В. Численные методы. - М.: Наука, 1989.

9. Форсайт Дж., Малькольм М., Моулер К. Машинные методы математических вычислений/Пер. с англ. - М.: Мир, 1980.

10. Хемминг Р. Численные методы/Пер. с англ. - М.: Наука, 1972.

11. Методические указания к лабораторным работам по дисциплине “Вычислительная математика” / Сост.: В. Н. Кафтасьев, А. Р. Лисс, М. С. Титов; СПбГЭТУ.- СПб., 1997.

12. Воробьева Г. Н., Данилова А. Н. Практикум по численным методам. - М.: Высш. шк., 1979.

13. Волков Е. А. Численные методы: Учеб. пособие. - М.: Наука, 1982.

14. Егоренков Д. Л., Фрадков А. Л., Харламов В. Ю. Основы математического моделирования. Построение и анализ моделей с примерами на языке MATLAB: Учеб. пособие/ Под ред. А. Л. Фрадкова; БГТУ. - СПб., 1994.

15. Компьютерная математика: Методические указания к лабораторным работам / Сост. И.А. Назаров; СПбГЭТУ. - СПб., 1993.

 

 

СОДЕРЖАНИЕ

 

Введение…….…………………………………………………………………………….3

1. Интерполирование функций.……………………………………………………….…4

1.1. Интерполяционные формулы для неравноотстоящих узлов.………………….….5

Лабораторная работа № 9…..……………………………………………………..….8

1.2. Интерполяционные формулы для равноотстоящих узлов…………………..…...17

Лабораторная работа № 10…………………………………………………….……18

2. Решение систем линейных алгебраических уравнений……………………….…...25

2.1. Решение систем линейных алгебраических уравнений методом Гаусса……..…25

Лабораторная работа № 11.……………………………………………………..…..26

2.2. Решение систем линейных алгебраических уравнений методом простой итерации…………………………………………………………………………….……26

Лабораторная работа № 12………………………………………………………......27

2.3. Программы для решения систем линейных алгебраических уравнений….….…28

Список литературы………………………………………………..……………..………30

Редактор Н. В. Рощина Лицензия ЛР № 020617 от 10.08.92 ___________________________________________________________________ Подписано в печать. Формат 60x84^1/16. Бумага тип. № 2. Печать офсетная. Усл. печ. л. 1,86. Уч.-изд. л. 2,0. Тираж 100 экз. Заказ Издательско-полиграфический центр СПбГЭТУ ___________________________________________________________________   197376, С.-Петербург, ул. Проф. Попова, 5