Решение системы уравнений в фортране

Настоящее пособие представляет собой вторую часть методических материалов, предназначенных для пользователей Комплекса программ PARALG для решения задач линейной алгебры с использованием распределенной памяти, реализуемого в НИВЦ МГУ (http://www.srcc.msu.su/num_anal/par_prog/).

Указанный Комплекс программ ориентирован на реализацию программ на ФОРТРАН’е — 77 в стиле SPMD (Single Program Multiple Data) с использованием для обмена информацией между параллельными процессами передачи сообщений с помощью примитивов MPI (Massage Passing Interface).

В настоящей второй части пособия рассматривается решение задач из раздела «Решение систем линейных алгебраических уравнений» с помощью программ, включенных к настоящему времени в указанный Комплекс.

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

Повторим здесь лишь вкратце основные действия, требующиеся от пользователя программ Комплекса PARALG.
Прежде чем обратиться к какой — либо из целевых параллельных программ пользователь должен выбрать для ее выполнения и описать так называемую «решетку процессов», т.е. «выстроить» параллельные процессы в общем случае в двумерное упорядоченное множество, состоящее из NPROW строк, в каждой из которых содержится по NPCOL процессов. Таким образом общее число используемых процессов NP = NPROW * NPCOL.

Поскольку обмен информацией между процессами происходит с помощью передачи сообщений, и для задач линейной алгебры используется специально разработанный для этой цели пакет BLACS, то необходимо в начале вычислений написать операторы инициализации этого пакета и определения контекста, соответствующего выбранной решетке процессов (см. пп.3,4).

По окончании вычислений необходимо также сделать выход из определенного ранее контекста (т.е. освободить использованные процессы решетки).

Кроме выбора решетки процессов, необходимо распределить исходные матрицы (векторы), обрабатываемые целевыми программами по параллельным процессам выбранной решетки.

Процесс распределения состоит из разбиения исходной матрицы (вектора), называемой глобальной, на прямоугольные (квадратные) блоки, содержащие по NB * MB элементов матрицы, с последующим «циклическим» распределением их по процессам решетки (см. пп.6). Это необходимо ввиду использования базовыми подпрограммами высокоэффективных блочных алгоритмов обработки матриц.

Таким образом при вычислениях в локальной памяти каждого из процессов содержится только некоторая (называемая локальной) часть исходной глобальной матрицы (вектора). Поэтому в головной фортранной программе, вызывающей целевую программу Комплекса, для размещения матрицы достаточно заказать память, равную по величине самой большой локальной части глобальной матрицы (т.к. размеры локальных частей на разных процессах могут несколько отличаться по величине) (см. пп.5.1, 6).

В связи с вышесказанным, при вызове целевых программ Комплекса помимо адреса локальной части матрицы (вектора) им необходимо передать в качестве фактического параметра адрес «дескрипторного массива» или «дескриптора». Информация, содержащаяся в дескрипторе, позволяет всегда однозначно определить в локальной памяти, какого процесса и на каком именно месте располагается каждый элемент глобальной матрицы A ( I, J ) (см. п.5.2).

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

Распределение блоков исходных глобальных матриц (векторов) по параллельным процессам значительно упрощается посредством использования соответствующих служебных подпрограмм (см. п.6.3).

Вопросы обработки ошибок и выдачи диагностических сообщений достаточно подробно рассматриваются в п.7 третьей части пособия.

2. Краткое описание содержания раздела Комплекса для решения систем линейных уравнений.

Одна из задач, решаемых Комплексом, — это линейные системы (системы линейных алгебраических уравнений)

Для решения системы (1) применяются специальные блочные реализации (блочные алгоритмы) для прямого и обратного исключения Гаусса.

Решение системы разбито на два этапа:

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

1. Виды факторизации, зависящие от свойств матрицы A, представлены ниже.

1.1. Для матриц общего вида применяется LU — факторизация с выбором ведущего элемента по столбцам:

1.2. Для симметричных или эрмитовых положительно определенных матриц применяется разложение Холецкого:

1.3. Для ленточных матриц общего вида применяется LU — факторизация с выбором ведущего элемента по столбцам:

1.4. Для ленточных матриц общего вида с диагональным преобладанием, включая трехдиагональные матрицы общего вида, применяется LU — факторизация без выбора ведущего элемента:

1.5. Для симметричных и эрмитовых положительно определенных ленточных матриц применяется разложение Холецкого:

1.6. Для симметричных и эрмитовых положительно определенных трехдиагональных матриц применяется LDL T — факторизация:

Указанные факторизации выполняются разными базовыми подпрограммами Комплекса:

• для плотных (вещественных или комплексных) матриц общего вида и симметричных (или эрмитовых) положительно определенных — PDGETRF (PZGETRF) и PDPOTRF (PZPOTRF);
• для ленточных матриц общего вида (вещественных или комплексных) — PDGBTRF (PZGBTRF) (с выбором ведущего элемента по столбцу) и PDDBTRF (PZDBTRF) (без выбора ведущего элемента);
• для симметричных/эрмитовых ленточных положительно определенных матриц — PDPBTRF (PZPBTRF);
• для трехдиагональных матриц общего вида (вещественных или комплексных) — PDDTTRF (PZDTTRF) (без выбора ведущего элемента);
• для симметричных/эрмитовых трехдиагональных положительно определенных матриц — PDPTTRF (PZPTTRF).

2. Второй этап — решение системы с использованием полученной на предыдущем этапе факторизации матрицы A выполняется следующими базовыми подпрограммами Комплекса:

• для плотных (вещественных или комплексных) матриц общего вида и симметричных (или эрмитовых) положительно определенных — PDGETRS (PZGETRS) и PDPOTRS (PZPOTRS);
• для ленточных матриц общего вида (вещественных или комплексных) — PDGBTRS (PZGBTRS) и PDDBTRS (PZDBTRS);
• для симметричных/эрмитовых ленточных положительно определенных матриц — PDPBTRS (PZPBTRS);
• для трехдиагональных матриц общего вида (вещественных или комплексных) — PDDTTRS (PZDTTRS);
• для симметричных/эрмитовых трехдиагональных положительно определенных матриц — PDPTTRS (PZPTTRS).

Приведем, наконец, список всех целевых программ, включенных в настоящее время в раздел Комплекса «Решение систем линейных алгебраических уравнений».

— подраздел «Решение систем с невырожденными матрицами общего вида»:

— подраздел «Решение систем с симметричными или эрмитовыми невырожденными матрицами»:

— подраздел «Решение систем с невырожденными матрицами специального вида»:

Видео:Решение систем уравнений методом подстановкиСкачать

Системы линейных алгебраических уравнений

Программа LA10 осуществляет разложение матрицы A в произведение двух треугольных матриц L и U по компактной схеме Гаусса с выбором главного элемента по столбцу:
TA = LU,
где L — нижняя треугольная матрица с единичной диагональю, U — верхняя треугольная матрица, T — матрица перестановок.
Матрица TA есть матрица A с переставленными строками. Поскольку в матрице T порядка N в каждой строке и в каждом столбце содержится ровно одна единица и N-1 нулей, то в программе вместо матрицы сохраняется только вектор, содержащий информацию о перестановках строк матрицы A .
Процесс разложения матрицы включает в себя N-1 шаг. На каждом шаге с номером p вычисляются промежуточные величины a * _ip :
_p-1
a * _ip = a_ip — ∑ l_iku_kj (i = p, p+1, . N) ,
k=1
где a_ip, l_ik, u_kp — элементы матриц A , L и U соответственно. Далее находится такой номер строки m , при котором величина |a * _ip| максимальна (для i от p до N ). Если m≠p , то строки с номерами m и p переставляются. Если a * _m,p = 0, то вычисления прекращаются.
Элементы p -го столбца матрицы L и p -ой строки матрицы U вычисляются по формулам:
l_ip = a * _ip / a_pp (i = p+1, p+2, . N) ,
u_pp = a * _pp ,
_N-1
u_pj = a_pj — ∑ l_pku_kj (j = p+1, p+2, . N) .
k=1
На последнем, дополнительном шаге, вычисляется элемент u_N,N :
_N-1
u_NN = a_NN — ∑ l_Nku_kN .
k=1
На этом процесс разложения заканчивается.

CALL LA10(A, AF, T, Det, Error)

A — входной параметр;
AF , T , Det , Error — выходные параметры.

Real A(1:N, 1:N) — массив, содержащий элементы исходной матрицы;
Real AF(1:N, 1:N) — массив, содержащий элементы треугольных матриц, единичная диагональ нижней треугольной матрицы не хранится. При вызове массивы A и AF можно совмещать в памяти, в этом случае результат работы программы записывается в массив A ;
Integer T(1:N) — массив, содержащий информацию о перестановках строк матрицы A ;
Real Det — переменная, содержащая значение определителя матрицы A ;
Integer Error — индикатор ошибки. Error =0, если ошибок нет и Error =65 или 66, если разложение матрицы не закончено и работа программы завершилась аварийно.

Программа LA11 решает систему линейных уравнений Ax=b с матрицей коэффициентов A , преобразованной при помощи программы LA10 . Программа LA10 должна закончиться без ошибок ( Error =0).
Решение системы Ax=b разбивается на решение двух треугольных систем Ly=b и Ux=y и выполняется по формулам:
_i-1
y_i = b_i — ∑ l_ijy_j (i = 1, 2, . N)
j=1
для первой системы (прямой ход исключения неизвестных) и
_N
x_i = (y_i — ∑ u_ijx_j)/u_jj (i = N, N-1, . 1)
j=i+1
для второй системы (обратный ход исключения неизвестных).
Перед решением первой системы массив свободных членов b переставляется согласно массиву T , который содержит информацию о перестановках матрицы A при работе программы LA10 .

CALL LA11(AF, T, B, X)

AF, T, B — входные параметры;
X — выходной параметр.

Real AF(1:N, 1:N) — массив, содержащий элементы треугольных матриц L и U ;
Integer T(1:N) — массив, содержащий информацию о перестановках строк матрицы A ;
Real B(1:N) — массив, содержащий вектор свободных членов;
Real X(1:N) — массив, содержащий вектор решения системы. При вызове массивы B и X можно совмещать в памяти, в этом случае результат работы программы записывается в массив B .

Программа LA12 предназначена для итерационного решения системы линейных уравнений Ax=b с матрицей коэффициентов A , преобразованной при помощи программы LA10 . Программа LA10 должна закончиться без ошибок ( Error =0).
Уточненное решение x (p+1) получается из приближенного x p следующим образом. Сначала вычисляется вектор невязки
r (p) = b (p) — Ax (p) .
Затем, решая систему
Ad (p) = r (p)
с использованием LU разложения матрицы A (прогр. LA10 ), получаем вектор поправки d (p) , который прибавляется к приближенному решению
x (p+1) = x (p) + d (p) .
Этот процесс повторяется до получения заданной точности или до достижения заданного количества итераций, начиная с нулевого приближения x (0) =0 .
На каждом шаге итерационного процесса, начиная с p>0 , проверяется выполнение неравенства
║d (p) ║ ≤ ½║d (p-1) ║ , (1)
где ║d (p) ║ = ∑|d_i (p) | — норма вектора поправок. Если неравенство (1) справедливо, то итерационный процесс продолжается.
Итерационный процесс заканчивается при появлении одной из следующих ситуаций:
1. Начиная с некоторого p для заданного ε>0 компоненты вектора поправок удовлетворяют неравенству
|d_i (p) | ≤ ε•|x_i (p) | (i = 1, 2, . N) , (2)
в этом случае x=x (p) выдается в качестве решения.
2. Начиная с некоторого p неравенство (1) перестает выполняться. В этом случае, если выполняется неравенство
║d (p) ║ ≤ ε•║x (p) ║ , (3)
то заданной точности удовлетворяет только норма полученного решения и в качестве решения выдается x=x (p) .
3. Начиная с некоторого p неравенство (1) перестает выполняться и неравенство (3) не выполняется также. Это свидетельствует о том, что итерационный процесс не сходится.
4. Неравенство (1) выполняется для всех p , но достигнуто максимальное количество итераций. Тогда в качестве решения выдается x=x (p+1) .

CALL LA12(A, AF, T, B, Eps, Max, X, Error)

A, AF, B, T, Eps, — входные параметры;
Max — входной и выходной параметр;
X, Error — выходные параметры.

Real A(1:N, 1:N) — массив, содержащий элементы исходной матрицы;
Real AF(1:N, 1:N) — массив, содержащий элементы треугольных матриц L и U ;
Integer T(1:N) — массив, содержащий информацию о перестановках строк матрицы A ;
Real B(1:N) — массив, содержащий вектор свободных членов;
Real Eps — заданная точность решения;
Integer Max — на входе: максимальное количество итераций, которое разрешено выполнить; на выходе: количество итераций, выполненное программой;
Real X(1:N) — массив, содержащий вектор решения системы;
Integer Error — индикатор ошибки;
Error=0 , если каждая компонента вектора x удовлетворяет заданной точности;
Error=1 , если только норма x удовлетворяет заданной точности;
Error=2 , если приближенное решение получено, но норма решения меньше заданной точности;
Error=3 , если достигнуто максимальное количество итераций, но норма решения не удовлетворяет заданной точности;
Error=65 , если итерационный процесс расходится.

Программа осуществляет разложение симметричной матрицы A с ненулевыми главными минорами в произведение диагональной и двух треугольных матриц с единичными диагоналями согласно уравнению
T T AT = L T DL ,
где L — верхняя треугольная матрица, L T — транспонированная к ней нижняя треугольная матрица, D — диагональная матрица, T — матрица перестановок.
С целью ограничения роста элементов в процессе преобразования в программе применяется выбор главного элемента на диагонали с одновременной перестановкой строк и столбцов. Такая перестановка сохраняет симметрию всех промежуточных матриц. Матрица перестановок T в программе не вычисляется, а вместо нее информация о перестановках сохраняется в линейном массиве.
Процесс преобразования осуществляется по формулам:
_i-1
d_i = ∑ a 2 _iid_k (i = 1, 2, . N),
k=1
_i-1
a_ij = ∑ a_ika_jk/d_j (i = 2, 3, . N; j = 1, 2, . N-1).
k=1

CALL LA15(A, D, Trp, T, Error)

A — входной и выходной параметр;
Trp — входной параметр;
D , T и Error — выходные параметры.

Real A(1:N, 1:N) — массив, содержащий элементы исходной матрицы. На выходе: ниже главной диагонали располагается матрица L , главная диагональ и элементы выше главной диагонали остаются без изменений;
Real D(1:N) — массив, который на выходе содержит элементы диагональной матрицы;
Logical Trp — логическая переменная. Если при вызове программы Trp=true, то программа осуществляет поиск главного элемента с соответствующими перестановками строк и столбцов, если Trp=false, то поиск главного элемента не производится и программа работает без перестановок;
Integer T(1:N) — массив, который на выходе содержит информацию о произведенных перестановках. Если при вызове программы Trp=false, то массив заполняется последовательно числами от 1 до N.
Integer Error — индикатор ошибки. Error =0, если ошибок нет и равно 65, если программа завершилась аварийно.

Программа LA16 решает систему линейных уравнений Ax=b с симметричной матрицей коэффициентов A , разложенной в произведение двух треугольных и диагональной матриц при помощи программы LA15 . Программа LA15 должна закончиться без ошибок ( Error =0).
Решение системы Ax=b разбивается на решение двух треугольных и одной линейной систем L T z=b , Dy=z , Lx=y и осуществляется по формулам:
_i-1
z_i = b_i — ∑ l_ijz_j (i = 1, 2, . N) ,
j=1

y_i = z_i/d_ii (i = 1, 2, . N) ,
_N
x_i = y_i — ∑ l_ijy_j (i = N, N-1, . 1) .
j=i+1
Перед решением первой системы переставляется массив b , а после решения последней системы переставляется массив x согласно массиву T , который содержит информацию о перестановках строк и столбцов матрицы A при работе программы LA15 .

CALL LA16(A, D, T, B, X)

A, D, T, B — входные параметры;
X — выходной параметр.

Real A(1:N, 1:N) — массив, содержащий матрицу L (ниже главной диагонали);
Real D(1:N) — массив, содержащий главную диагональ матрицы D ;
Integer T(1:N) — массив, содержащий информацию о перестановках строк и столбцов;
Real B(1:N) — массив, содержащий вектор свободных членов;
Real X(1:N) — массив, содержащий вектор решения системы. При вызове массивы B и X можно совмещать в памяти, в этом случае результат работы программы записывается в массив B .

Программа LA17 предназначена для итерационного решения системы линейных уравнений Ax=b с симметричной матрицей коэффициентов A , преобразованной при помощи программы LA15 . Программа LA15 должна закончиться без ошибок ( Error =0).
Описание алгоритма программы полностью совпадает с описанием алгоритма программы LA12 .

CALL LA17(A, D, T, B, Eps, Max, X, Error)

A, D, B, T, Eps — входные параметры;
Max — входной и выходной параметр;
X, Error — выходные параметры.

Real A(1:N, 1:N) — массив, содержащий матрицу L (ниже главной диагонали);
Real D(1:N) — массив, содержащий главную диагональ матрицы D ;
Integer T(1:N) — массив, содержащий информацию о перестановках строк и столбцов;
Real B(1:N) — массив, содержащий вектор свободных членов;
Real Eps — заданная точность решения;
Integer Max — на входе: максимальное количество итераций, которое разрешено выполнить; на выходе: количество итераций, выполненное программой;
Real X(1:N) — массив, содержащий вектор решения системы. При вызове массивы B и X можно совмещать в памяти, в этом случае результат работы программы записывается в массив B .
Integer Error — индикатор ошибки;
Error=0 , если каждая компонента вектора x удовлетворяет заданной точности;
Error=1 , если только норма x удовлетворяет заданной точности;
Error=2 , если приближенное решение получено, но норма решения меньше заданной точности;
Error=3 , если достигнуто максимальное количество итераций, но норма решения не удовлетворяет заданной точности;
Error=65 , если итерационный процесс расходится.

Программа решает систему линейных уравнений с трехдиагональной матрицей коэффициентов методом прогонки. Исходная система может быть представлена в виде

где a_i — элементы главной диагонали, b_i — верхней диагонали, c_i — нижней диагонали.
Решение системы проводится в два этапа. В процессе прямого хода прогонки диагональные элементы в каждом уравнении выражаются через внедиагональные и последовательно исключается нижняя диагональ исходной матрицы коэффициентов. Получаем следующую систему уравнений

Исключая из них x₁ , получаем:

Подставляя найденное решение в систему (**), получаем окончательное решение исходной системы.

CALL LA30(A, B, C, D, X, Error)

A, B, C, D — входные параметры;
X, Error — выходные параметры.

Real A(1:N) — массив, содержащий элементы главной диагонали;
Real B(1:N) — массив, содержащий элементы верхней диагонали;
Real C(1:N) — массив, содержащий элементы нижней диагонали;
Real D(1:N) — массив, содержащий вектор свободных членов;
Real X(1:N) — массив, содержащий вектор решения;
Integer Error — индикатор ошибки
Error = 0 , если решение системы получено;
Error = 1 , если при вызове программы N Error = 65 , если матрица системы вырожденная.

Программа решает систему линейных уравнений с пятидиагональной матрицей коэффициентов методом прогонки. Исходная система может быть представлена в виде

где a_i — элементы главной диагонали, b_i и c_i — элементы верхних диагоналей, e_i и d_i — элементы нижних диагоналей.
Для нахождения коэффициентов прямого хода диагональные элементы в каждом уравнении выражаются через внедиагональные и последовательно исключаются две нижние диагонали исходной матрицы коэффициентов. Получаем следующую систему уравнений

Вычислением коэффициентов p_i , q_i , r_i , s_i и t_i по этим рекуррентным формулам заканчивается прямой ход прогонки. Для нахождения коэффициентов обратного хода в системе (*) начиная с x_n-2 последовательно подставляем нижнее уравнение в два верхних и тем самым исключаем две верхние диагонали исходной матрицы коэффициентов. Система уравнений будет преобразована к виду

На этом заканчивается обратный ход прогонки.
Значение переменных x_n-1 и x_n можно найти из совместно решаемых первых двух уравнений системы (**) и уравнений

Исключая переменные x₁ и x₂ , получаем

Решая эту систему по методу Крамера, имеем

CALL LA31(A, B, C, D, E, G, X, Error)

A, B, C, D, E, G — входные параметры;
X, Error — выходные параметры.

Real A(1:N) — массив, содержащий элементы главной диагонали;
Real B(1:N), C(1:N) — массивы, содержащие элементы верхних диагоналей;
Real D(1:N), E(1:N) — массивы, содержащие элементы нижних диагоналей;
Real G(1:N) — массив, содержащий вектор свободных членов;
Real X(1:N) — массив, содержащий вектор решения;
Integer Error — индикатор ошибки
Error = 0 , если решение системы получено;
Error = 1 если при вызове программы N Error = 65 , если матрица системы вырожденная.

Программа LA44 предназначена для итерационного решения системы линейных уравнений Ax=b с симметричной полжительно-определенной матрицей коэффициентов A методом сопряженных градиентов.
Итерационный процесс строится по следующей схеме. Исходя из начального приближения x (0) полагают s (1) = r (0) = Ax (0) — b и последовательно вычисляют рекуррентные соотношения
r (k) = r (k-1) — α (k) As (k) , (1)
x (k) = x (k-1) — α (k) s (k) , (2)
s (k+1) = r (k) + β (k) s (k) , (3)
k=1, 2, . ,
где
α (k) = (r (k-1) , r (k-1) )/(As (k) , s (k) )
β (k) = (r (k) , r (k) )/(r (k-1) , r (k-1) ) .
Теоретически векторы r (k) и x (k) должны удовлетворять условиям
r (k) = Ax (k) — b, (4)
(r (i) , r (k) ) = 0 для i ≠ k , (5)
r (n) = 0, x (n) = x ∞ , (6)
где x ∞ — точное решение, а n — размерность матрицы системы. Однако из-за наличия ошибок округления при вычислениях реальный результат отличается от теории и соотношения (4), (5) и (6) не выполняются. Поэтому полученный после n итераций вектор x (n) может значительно отличаться от точного решения. Чтобы получить достаточно близкое приближение к точному решению в общем случае число итераций k должно превосходить размерность системы n . При этом рано или поздно наступит момент, когда невязка r (k) начнет резко уменьшаться и в конечном итоге точность решения будет определяться числом обусловленности матрицы A .
Вычисления прекращаются, как только выясняется, что продолжение итераций не улучшает решения. Критерием окончания процесса служит оценка малости величины ║r (k) ║ 2 = (r (k) , r (k) ) . Поскольку вектор невязок r (k) в алгоритме вычисляется рекуррентно через r (k-1) и s (k) , то из-за наличия ошибок округления полученное значение r (k) несколько отличается от величины истинной невязки r *(k) , полученной из вектора Ax (k) — b . По этой разнице можно судить о точности вычислений.
В рассматриваемом алгоритме вычисляется величина
σ = ║r *(k) — r (k) ║ 2 , (7)
где вектор r (k) определяется из выражения (1), а вектор r *(k) из выражения (4). Величина σ возрастает с ростом k . Как только ║r (k) ║ 2 и σ становятся величинами одного порядка, это означает, что невязки r (k) и r *(k) достигли уровня ошибок округления.
Однако непосредственно использовать неравенство ║r (k) ║ 2 нельзя, т.к. может оказаться, что величина ║r (k) ║ 2 все время чуть-чуть больше σ . Поэтому в алгоритме используется неравенство
║r (k) ║ 2 , (8)
в котором ω определяется выражением
ω = exp((k/n) 2 ) , (9)
где k — номер итерации, n — порядок системы уравнений.
Постоянно возрастающий коэффициент ω в выражении (8) обеспечивает окончание процесса даже при очень плохо обусловленной матрице A и, кроме того, позволяет завершить вычисления, когда ║r (k) ║ 2 и σ быстро приближаются друг к другу.
Т.к. определение величины σ требует значительного времени из-за необходимости вычислять истинную невязку, то значение σ вычисляется через √ — n шагов и остается постоянным до вычисления нового значения σ . Неравенство (8) проверяется на каждом шаге алгоритма.

CALL LA44(Amv, B, X, Max, Error)

Amv, B — входные параметры;
X, Max — входные и выходные параметры;
Error — выходной параметр.

Amv — процедура, которая вычисляет произведение Av для произвольного вектора v ;
Real B(1:N) — массив, содержащий вектор свободных членов;
Real X(1:N) — массив, который на входе содержит вектор начального приближения, если оно известно. В этом случае параметр Max должен быть больше нуля. Если начальное приближение неизвестно, то параметр Max должен быть отрицательным. На выходе массив содержит вектор решения системы;
Max — параметр, который на входе задает максимально возможное количество итераций. Параметр должен быть больше нуля, если в массиве X содержится начальное приближение и должен быть меньше нуля, если оно не известно. Для определения максимально возможного количества итераций программа использует абсолютное значение параметра. На выходе параметр содержит число итераций, которое потребовалось программе для получения приближения к решению;
Integer Error — индикатор ошибки
Error = 0 , если для получения решения потребовалось не более Max итераций;
Error = 1 , если для получения решения Max итераций оказалось недостаточно;
Error = 65 , если матрица А либо плохо обусловлена и ее наименьшее собственное значение имеет величину порядка ошибки округления, либо не является положительно-определенной. В обоих случаях решение получить нельзя.

Видео:Решение системы уравнений методом ГауссаСкачать

Мак-Кракен Д., Дорн У. Численные методы и программирование на Фортране

ПЕРЕВОД С АНГЛИЙСКОГО Б. Н. КАЗАКА

ПОД РЕДАКЦИЕЙ И С ДОПОЛНЕНИЕМ Б. М. НАЙМАРКА

Издание второе, стереотипное

ИЗДАТЕЛЬСТВО «МИР» Москва 1977

Книга является руководством по структуре и использованию алгоритмического языка ФОРТРАН при решении вычислительных задач на современных электронных машинах.

Специфика и простота трансляторов для ФОРТРАНа, эффективность оттранслированных программ и методика выявления и оценки ошибок выгодно отличают этот язык от других алгоритмических языков. Все это обусловило широкое внедрение ФОРТРАНа в технику программирования за рубежом.

Одновременно в книге подробно излагаются тщательно отобранные численные методы, применение которых иллюстрируется на многочисленных практических примерах.

Объединение численных методов и основ программирования на ФОРТРАНе делает эту книгу полезной для широкого круга читателей, как для студентов и аспирантов втузов, так и для инженеров и специалистов по теории программирования.

Редакция литературы по математическим наукам

©Перевод на русский язык, «Мир», 1977 041(01) – 77

Оглавление книги Численные методы и программирование на Фортране
Предисловие редактора перевода
К второму изданию
Предисловие
Глава 1. Основы программирования на Фортране
1.1 Применение цифровых вычислительных машин
1.2. Последовательные этапы в «решении задачи» с помощью ЭЦВМ
1.3. Программа на Фортране
1.4. Константы. Упражнения
1.5. Переменные и наименование переменных. Упражнения
1.6. Операции и выражения. Упражнения
1.7. Математические функции
1.8. Арифметические операторы. Упражнения
1.9. Ввод и вывод. Упражнения
1.10. Передача управления. Операторы GO ТО и IF. Упражнения
1.11. Операторы PAUSE, STOP и END. Упражнения
1.12. Написание программы, ее перфорация на перфокартах и постановка ее на ЭЦВМ
1.13. Практический пример 1: Площадь треугольника

Глава 2. Ошибки
2.1. Введение
2.2. Относительные и абсолютные ошибки
2.3. Ошибки, содержащиеся в исходной информации
2.4. Ошибки ограничения
2.5. Ошибки округления
2.6. Распространение ошибок
2.7. Графы вычислительных процессов
2.8. Примеры
2.9. Памятка программисту
Упражнения

Глава 3. Практическое вычисление функций
3.1. Введение
3.2. Степенные ряды
3.3. Полиномы Чебышева
3.4. Экономизация степенных рядов
3.5. Вычисление ряда
3.6. Рациональные приближения и непрерывные дроби
3.7. Элементарные функции
3.8. Практический пример 3: Ошибки при прямом вычислении синуса по ряду Тейлора
Упражнения

Глава 4. Некоторые простые программы
4.1. Введение
4.2. Практический пример 4: Расчет колонны
4.3. Частотная характеристика сервомеханизма. Отладка программы
4.4. Практический пример 6: Интеграл вероятностей

Глава 5. Численное решение уравнений
5.1. Введение
5.2. Метод последовательных приближений
5.3. Усовершенствованный метод последовательных приближений
5.4. Метод Ньютона — Рафсона
5.5. Случай почти равных корней
5.6. Сравнение методов и их ошибок округления
5.7. Корни многочленов
5.8. Влияние неточности коэффициентов многочлена
5.9. Системы уравнений
5.10. Комплексные корни
5.11. Нахождение исходного приближения
5.12. Практический пример 7: Процесс роста монокристалла из пара
Упражнения

Глава 6. Численное интегрирование
6.1. Введение
6.2. Правило трапеций
6.3. Ошибка ограничения для метода трапеций
6.4. Ошибки округления при использовании метода трапеций
6.5. Экстраполяционный переход к пределу
6.6. Правило Симпсона
6.7. Метод Гаусса
6.8. Численные примеры и сравнение методов
6.9. Практический пример 8: Светимость электрической лампочки
Упражнения

Глава 7. Переменные с индексами и оператор DO
7.1. Определения
7.2. Примеры использования переменных с индексами
7.3. Для чего нужны переменные с индексами?
7.4. Оператор DIMENSION
7.5. Допустимые формы индексов
7.6. Оператор DO
7.7. Дальнейшие определения
7.8. Правила использования оператора DO
7.9. Дальнейшие примеры использования оператора DO
7.10. Практический пример 9: Линейная интерполяция
Упражнения

Глава 8. Системы линейных алгебраических уравнений
8.1. Введение
8.2. Метод исключения (метод Гаусса)
8.3. Ошибки округления
8.4. Уточнение решения
8.5. Влияние погрешностей коэффициентов. Достижимая точность решения
8.6. Итерационные методы решения систем линейных уравнений
8.7. Сравнение методов
8.8. Практический пример 10: Проведение кривой методом наименьших квадратов. Упражнения

Глава 9. Функции, подпрограммы и вспомогательные операторы
9.1. Введение
9.2. Функции, предусмотренные в программе-трансляторе
9.3. Арифметический оператор-функция
9.4. Подпрограммы FUNCTION и SUBROUTINE
9.5. Таблица основных характеристик функций и подпрограмм
9.6. Операторы EQUIVALENCE и COMMON
9.7. Практический пример 11: Решение квадратных уравнений с помощью подпрограмм. Упражнения

Глава 10. Обыкновенные дифференциальные уравнения
10.1. Введение
10.2. Решение с помощью рядов Тейлора
10.3. Методы Рунге — Кутта
10.4. Анализ ошибок, возникающих при использовании методов Рунге — Кутта
10.5. Методы прогноза и коррекции
10.6. Анализ ошибок при использовании методов прогноза и коррекции
10.7. Достижимая точность
10.8. Сравнение методов
10.9. Практический пример 12: Полет сверхзвукового самолета. Упражнения

Глава 11. Уравнения в частных производных
11.1. Введение и некоторые определения
11.2. Разностные уравнения
11.3. Эллиптические уравнения
11.4. Решение эллиптического разностного уравнения
11.5. Гиперболические уравнения
11.6. Решение гиперболического разностного уравнения
11.7. Параболические уравнения
11.8. Решение параболического разностного уравнения
11.9. Практический пример 13: Распределение температуры в трубе квадратного сечения. Упражнения

Приложение 1. Сводка методов ввода и вывода информации в ФОРТРАНе
П.1.1. Основные сведения
П.1.2. Список переменных в операторе ввода — вывода
П.1.3. Оператор FORMAT
П.1.4. Дополнительные приемы построения оператора FORMAT
П.1.5. Операции с магнитной лентой

Приложение 2. Некоторые употребительные математические формулы
Ответы к упражнениям

Дополнение. Сводка основных правил программирования на языке ФОРТРАН. Б. М. Наймарк
Д.1. Основные символы языка ФОРТРАН
Д.2. Числа
Д.З. Переменные без индексов
Д.4. Индексы
Д.5. Переменные с индексами
Д.6. Выражения
Д.7. Функции
Д.8. Операторы
Д.9. Описательные операторы
Д.10. Исполнимые операторы
Д.11. Примечания

Эта книга должна дать читателю, специализирующемуся в области техники или науки, основные понятия о численном решении различных задач на современных электронных цифровых вычислительных машинах.

Основное внимание в книге уделено тщательно подобранным и весьма полезным практическим методам исследования различных численных задач, встречающихся в современной науке и технике. Изложенные в книге методы полезны не только сами по себе, но они к тому же подготавливают читателя для изучения более сложных методов и помогают понять их во всей глубине. Точно так же, как физику необходимо иметь интуитивное понятие об инерции до того, как он сможет глубоко изучить механику, так и изучение численных методов должно предшествовать тщательному изучению численного анализа. Цель книги, таким образом, двойная: дать читателю необходимые знания для практических вычислений и заложить основы для дальнейшего изучения предмета.

Объединение численных методов и программирования на ФОРТРАНЕ — это гораздо больше, нежели просто издание двух обычных книг в общей обложке. Сведения по программированию на ФОРТРАНЕ появляются в книге в том порядке, в котором это требуется для различных классов задач и для соответствующих им методов численного решения. Основные примеры программирования также в большинстве случаев иллюстрируют те или иные численные методы. Очень важно также, что во всей книге используются современные понятия. Например, весьма важный вопрос об ошибках округления трактуется в терминах арифметики с плавающей запятой (арифметики нормализованных чисел с постоянным количеством значащих цифр в мантиссе), введенные там же графические представления о вычислительных процессах облегчают полное понимание.

Численные методы, включенные в книгу, были тщательно проанализированы с точки зрения их пригодности для использования на цифровых вычислительных машинах.

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

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

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

Тринадцать практических примеров включены в книгу по трем причинам. Во-первых, необходимо было проиллюстрировать тот материал, который относится к численным методам и к программированию. Во-вторых, нужно было показать, как различные аспекты численных методов и программирования взаимодействуют и сочетаются между собой. (В связи с этим в нескольких практических примерах показаны последствия легкомысленного подхода к ошибкам вычисления.) В-третьих, необходимо было дать понятие о диапазоне областей науки и техники, в которых могут использоваться вычислительные машины. Содержание практических примеров, однако, таково, что понимание материала не зависит от степени знакомства читателя с данной областью науки или техники. Во многих случаях в книге воспроизведена печать на выходном бланке вычислительной машины.

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

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

На основе этой книги можно построить лекционный курс на один семестр (по четыре часа занятий в неделю), если у студентов имеется возможность работать на вычислительной машине и если ожидается, что студенты выполнят достаточно большое число упражнений. С другой стороны, книгу можно использовать и в качестве пособия при гораздо более кратком курсе; для этого можно ограничиться меньшим числом упражнений или пропустить некоторые главы. Более того, главы, посвященные численным методам, не очень жестко связаны между собой. Преподаватель, который пожелает, например, уделить меньше внимания главе 3 о вычислении функций и сконцентрировать внимание аудитории на системах линейных уравнений, довольно легко может это сделать. Возможны и многие другие способы построения учебного курса на основе этой книги.

Мы с большой признательностью отмечаем сотрудничество тех, кто помогал нам в написании и подготовке этой книги: мисс Агнесы Кульке и Джеймса Грисмера из IBM, Уильяма Найбека из инженерного колледжа Милуоки, Фреда Гринбергера из РЭНД Корпорейшен, Джека Хол-лингсуорта и Поля Макслойна из политехнического института Реннселера, Юджина Голуба и Гарольда Р. ван Зорена из Стэнфордского университета и Чарльза Дэвидсона из Висконсинского университета. Помощь этих специалистов была для нас очень ценной: своим участием они помогли сформировать эту книгу во многих очень важных деталях. Наконец, мы выражаем свою признательность фирме IBM, ив частности доктору Герману X. Голдстайну, без поддержки и участия которого эта книга не была бы написана.

Даниель Д. Мак-Кракен
Уильям С. Дорн

Скачать книгу Мак-Кракен Д., Дорн У. Численные методы и программирование на Фортране. Издательство «Мир», Москва, 1977

📹 Видео

Матричный метод решения систем уравненийСкачать

Метод Крамера за 3 минуты. Решение системы линейных уравнений - bezbotvyСкачать

ФСР. Система однородных уравнений. Общее решениеСкачать

Cистемы уравнений. Разбор задания 6 и 21 из ОГЭ. | МатематикаСкачать

Решение систем уравнений второго порядка. 8 класс.Скачать

Решение системы уравнений методом Гаусса 4x4Скачать

Метод Ньютона (метод касательных) Пример РешенияСкачать

9 класс, 11 урок, Методы решения систем уравненийСкачать

Математика | Система уравнений на желтую звездочку (feat Золотой Медалист по бегу)Скачать

Решение системы уравнений методом обратной матрицы - bezbotvyСкачать

Решение систем уравнений методом сложенияСкачать

Решение системы уравнений методом Крамера 2x2Скачать

Решение системы линейных уравнений с двумя переменными способом подстановки. 6 класс.Скачать

МЕТОД ПОДСТАНОВКИ 😉 СИСТЕМЫ УРАВНЕНИЙ ЧАСТЬ I#математика #егэ #огэ #shorts #профильныйегэСкачать

2 минуты на формулы Крамера ➜ Решение систем уравнений методом КрамераСкачать

Решение систем уравнений. Методом подстановки. Выразить YСкачать

Система уравнений. Метод алгебраического сложенияСкачать

Решение системы линейных уравнений графическим методом. 7 класс.Скачать