Void MultMatr(double[,] a, double[,] b, double[,] c)

{

int m = A.GetLength(0);

int n = A.GetLength(1);

int p = B.GetLength(1);

double S =0;

for(int i=0; i < m; i++)

for (int j = 0; j < p; j++)

{

S = 0;

for (int k = 0; k < n; k++)

S += A[i, k] * B[k, j];

C[i, j] = S;

}

}

Void FillDg()

{

for (int i = 0; i < m; i++)

for (int j = 0; j < p; j++)

dataGridView3.Rows[i].Cells[j].Value

= C[i, j].ToString();

}

Задачи (ввод, вывод и другие простые задачи с массивами)

250. Организуйте в консольном приложении ввод и вывод одномерного массива строкового типа.

251. Организуйте в Windows приложении ввод и вывод одномерного массива строкового типа.

252. Организуйте в консольном приложении ввод массива «Сотрудники», содержащего фамилии сотрудников. Введите массив «Заявка», элементы которого содержат фамилии сотрудников и, следовательно, должны содержаться в массиве сотрудников. Обеспечьте контроль корректности ввода данных.

253. Организуйте в Windows приложении ввод массива «Сотрудники», содержащего фамилии сотрудников. Введите массив «Заявка», элементы которого содержат фамилии сотрудников и, следовательно, должны содержаться в массиве сотрудников. Обеспечьте контроль корректности ввода данных.

254. Организуйте в Windows приложении ввод массива «Сотрудники», содержащего фамилии сотрудников. Создайте массив «Заявка», элементы которого должны содержаться в массиве сотрудников. Для создания массива «Заявка» постройте форму «Два списка», содержащую два элемента ListBox, источником данных для первого из них служит массив «Сотрудники». Пользователь переносит данные из первого списка во второй, формируя данные для массива «Заявка». После формирования данные переносятся в массив. Для построения формы используйте шаблон, описанный в лекции 24 учебника.

255. Организуйте в консольном приложении ввод и вывод двумерного массива строкового типа.

256. Организуйте в Windows приложении ввод и вывод двумерного массива строкового типа.

257. Организуйте в консольном приложении ввод массива «Сотрудники» из двух столбцов, содержащего фамилии и имена сотрудников. Введите массив «Заявка» той же структуры, элементы которого должны содержаться в массиве сотрудников. Обеспечьте контроль корректности ввода данных. Организуйте вывод обоих массивов.

258. Организуйте в Windows приложении ввод массива «Сотрудники» из двух столбцов, содержащего фамилии и имена сотрудников. Введите массив «Заявка» той же структуры, элементы которого должны содержаться в массиве сотрудников. Обеспечьте контроль корректности ввода данных. Организуйте вывод обоих массивов.

259. (*) Организуйте в консольном приложении ввод и вывод массива «Машины», содержащего 4 столбца: «Владелец», «Марка», «Номер», «Год Выпуска». При вводе данных обеспечьте их корректность. Поле «Владелец» должно быть строкой в формате «фамилия имя», где фамилия и имя должны начинаться с большой буквы и состоять из букв алфавита кириллицы, включая дефис. Номер машины должен соответствовать формату, принятому для номеров машин. При выводе сохраняйте структуру массива.

260. (*) Организуйте в Windows приложении ввод и вывод массива «Машины», содержащего 4 столбца: «Владелец», «Марка», «Номер», «Год Выпуска». При вводе данных обеспечьте их корректность. Поле «Владелец» должно быть строкой в формате «фамилия имя», где фамилия и имя должны начинаться с большой буквы и состоять из букв алфавита кириллицы, включая дефис. Номер машины должен соответствовать формату, принятому для номеров машин. При выводе сохраняйте структуру массива.

261. (*) В консольном приложении уже построен массив «Машины» (см. задача 3.9) . Построить массив «Цветные машины», в котором к столбцам массива «Машины» добавляется 5-й столбец «Цвет». Организуйте диалог с пользователем, выясняя цвет для каждой машины из массива «Машины».

262. (*) В Windows приложении уже построен массив «Машины» (см. задача 3.10) . Построить массив «Цветные машины», в котором к столбцам массива «Машины» добавляется 5-й столбец «Цвет». Организуйте диалог с пользователем, выясняя цвет для каждой машины из массива «Машины».

263. Организуйте в консольном приложении ввод и вывод одномерного массива арифметического типа (от byte до double).

264. Организуйте в Windows приложении ввод и вывод одномерного массива арифметического типа (от byte до double).

265. Организуйте в консольном приложении ввод массива «Сотрудники», содержащего фамилии сотрудников, и массива «Зарплата». Обеспечьте контроль корректности ввода данных о зарплате, проверяя диапазон возможных значений.

266. Организуйте в Windows приложении ввод массива «Сотрудники», содержащего фамилии сотрудников, и массива «Зарплата». Обеспечьте контроль корректности ввода данных о зарплате, проверяя диапазон возможных значений.

267. Организуйте в консольном приложении ввод и вывод матрицы - двумерного массива арифметического типа.

268. Организуйте в Windows приложении ввод и вывод матрицы - двумерного массива арифметического типа.

269. Организуйте в консольном приложении ввод массива декартовых координат n точек на плоскости. Вычислите массив полярных координат этих точек и организуйте вывод этого массива. Обеспечьте контроль вводимых значений.

270. Организуйте в Windows приложении ввод массива декартовых координат n точек на плоскости. Вычислите массив полярных координат этих точек и организуйте вывод этого массива. Обеспечьте контроль вводимых значений.

271. Организуйте в консольном приложении ввод массива полярных координат n точек на плоскости. Вычислите массив декартовых координат этих точек и организуйте вывод этого массива. Обеспечьте контроль вводимых значений.

272. Организуйте в Windows приложении ввод массива полярных координат n точек на плоскости. Вычислите массив декартовых координат этих точек и организуйте вывод этого массива. Обеспечьте контроль вводимых значений.

273. Организуйте в консольном приложении ввод массива декартовых координат n точек в трехмерном пространстве. Вычислите массив полярных координат этих точек и организуйте вывод этого массива. Обеспечьте контроль вводимых значений.

274. Организуйте в Windows приложении ввод массива декартовых координат n точек в трехмерном пространстве. Вычислите массив полярных координат этих точек и организуйте вывод этого массива. Обеспечьте контроль вводимых значений.

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

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

Массивы и классические алгоритмы математики

Полиномы

Полиномом n-й степени P_n(x) называют функцию:

(1)

Если рассматривать график этой функции на плоскости, то x и P_n(x) – это декартовы координаты точек графика функции. Значения a_k (k из интервала [0,n]) называются коэффициентами полинома. Все они принадлежат одному типу и при программной работе с полиномами представляются одномерным массивом с n+1 элементами.

Если задан массив коэффициентов полинома A, то вычислить значение полинома в точке x не представляет особой сложности. Но ни один уважающий себя программист не позволит себе вычислять значение полинома, буквально пользуясь схемой (1), требующей n-1 операций возведения в степень, n операций умножения и n операций сложения. Прекрасный пример того, как можно упростить алгоритм, дает схема Горнера, вообще не требующая возведения в степень. В этой схеме полином P_n(x) представляют в виде:

(2)

Удобнее представлять схему Горнера в рекуррентной форме:

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

Если P_n(x) – полином n-й степени с коэффициентами a_i, Q_n(x) – полином n-й степени с коэффициентами b_i и P_n(x) = Q_n(x), то из этого следует равенство соответствующих коэффициентов:

(3)

Многие задачи над полиномами связаны с определением их корней. Напомню, x₀ является корнем полинома, если P_n(x₀) = 0. У полинома n-й степени не более чем n действительных корней. Если n – нечетно, то полином имеет хотя бы один действительный корень. Все корни полинома принадлежат некоторому конечному интервалу [c, d]. Вне этого интервала поведение полинома определяется его старшим членом – a_nxⁿ. Для полинома четной степени обе ветви уходят в +∞, если a_n >0 и в -∞, если a_n<0. Для полинома нечетной степени ветви полинома вне интервала [c, d] разнонаправлены. Если a_n>0, то правая ветвь уходит в +∞, а левая ветвь - в -∞. Если a_n<0, то левая ветвь уходит в +∞, а правая ветвь в -∞.

Когда по каким-либо физическим соображениям интервал [c, d] известен хотя бы приблизительно, то задача нахождения корней полинома облегчается, в противном случае она может быть довольно трудной, особенно для случая близко расположенных корней.

Исследование интервала

Рассмотрим один из простых алгоритмов, исследующих существует ли на заданном интервале [e, f] хотя бы один корень. Один корень заведомо существует, если полином на концах исследуемого интервала имеет разные знаки или один из концов интервала уже является корнем полинома. Это условие и будет характерным признаком поиска нужного интервала. Если исходный интервал [e, f] удовлетворяет характерному признаку, то задача решена и такой интервал найден. В противном случае в цикле по k вычислим h = L/2^k, где L – длина исходного интервала (L= f-e). Затем организуем внутренний цикл, в котором проверим характерный признак на всех интервалах длины h. Если интервал будет найден, то вычисления завершаются, в противном случае переходим к следующему шагу цикла по k, производя очередное дробление h. Завершение цикла по k означает, что если исследуемый интервал [e, f] и содержит корни, то это близкие пары корней, отстоящие друг от друга на расстояние, меньшее h – заключительной длины интервала по завершении цикла по k.

Приведу несколько практических рекомендаций, полезных при реализации этой схемы. Внутренний цикл следует организовать так, чтобы не повторять вычисление полинома в тех точках, в которых это вычисление проводилось на предыдущих шагах цикла. Это означает, что когда шаг h = L/2^k, то во внутреннем цикле достаточно вычислить значение полинома не более чем в 2^k-1 точках. Внешний цикл достаточно ограничить числом в интервале от 10 до 20, поскольку уже при k=10 величина исходного интервала L уменьшится более чем в 1000 раз, что вполне достаточно в большинстве практических ситуаций. Хотя следует помнить, что в ряде ситуаций практики приходится иметь дело с резко осциллирующими функциями, где близкие корни являются правилом, а не исключением.

Алгоритмы нахождения корня полинома

Рассмотрим несколько простых схем нахождения корня полинома. Заметим, что все эти схемы применимы к нахождению корней любых функций, а не только полиномов. Как всегда в программировании речь идет не столько о точном нахождении корня, сколько о нахождении корня с заданной точностью ε. Так что, если x₀ – это точное значение корня, то нам достаточно найти x^* - такое что |x₀ – x^*| < ε.

Схема дихотомии отрезка (деление пополам):

Эта схема прекрасно подходит, когда предварительно проведено исследование интервала существования корня и найден такой интервал [e, f], на концах которого полином принимает разные знаки, так что существует корень внутри интервала. Если исходный интервал мал и сравним с заданной точностью ε, то в качестве корня можно выбрать середину этого интервала. Если же исходный интервал больше, чем значение ε, то интервал можно разделить пополам и из двух половинок выбрать ту, для которой выполняется характерный признак существования корня. Понятно, что если признак выполняется для всего интервала, то он обязательно будет выполняться для одной из его половинок. Деление отрезка пополам приводит к быстрому уменьшению длины отрезка, так что 10 - 20 делений достаточно, чтобы найти интервал длины, меньшей ε, а следовательно и корень полинома с заданной точностью.

Метод простой итерации:

Формально метод применим и в том случае, когда неизвестен интервал, в котором существует корень функции. Пусть x₀ – некоторое заданное начальное приближение к корню полинома. Тогда можно построить следующий итерационный процесс:

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

Метод простой итерации обладает полезным свойством «неподвижной точки». Корни функции являются «неподвижными точками» метода. Нетрудно заметить, что если на некотором шаге x_k-1 сошлось к корню полинома x^*, то x_k и все последующие итерации будут равны x^*, так что итерационный процесс из найденного корня не уходит.

Метод Ньютона

Этот метод чуть более сложен в реализации, но обладает лучшей сходимостью в сравнении с методом простой итерации, хотя и здесь сходимость во многом зависит от удачного выбора начального приближения x₀. Для произвольной функции f итерационный процесс метода Ньютона выглядит так:

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

Понижение степени полинома

Если для полинома P(x) n-й степени найден корень x₁, то можно понизить степень полинома, построив полином P1(x) степени n-1, у которого все корни совпадают с корнями полинома P(x) за исключением того, что у него нет корня x₁.

Запишем соотношение, связывающее полиномы:

Учитывая соотношение (3) о равенстве двух полиномов одной степени, можно выписать n+1 соотношение, связывающее коэффициенты этих полиномов. Эти соотношения нетрудно разрешить относительно неизвестных коэффициентов b_k. В результате получим:

(4)

Заметьте, неизвестных всего n, а уравнений можно построить – n+1. Но последнее уравнение (a0 + b₀x1 = 0) является следствием предыдущих и используется для контроля вычислений.

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

Нахождение коэффициентов полинома по его корням

До сих пор рассматривалась задача отыскания корней полинома с заданными коэффициентами. Иногда приходится решать обратную задачу – найти коэффициенты полинома, если известны его корни – x₁, x₂, … x_n. Полиномов с одинаковыми корнями существует бесчисленное множество. Однако среди них существует единственный полином с коэффициентом a_n, равным единице. Этот полином называется приведенным, его то и будем строить. Все остальные полиномы получаются из приведенного полинома умножением всех коэффициентов на произвольное число a_n, от которого требуется лишь, чтобы оно не было равно нулю. Поэтому для однозначного решения задачи требуется задать n корней и коэффициент при старшем члене полинома. Тогда можно записать следующее равенство:

Для нахождения коэффициентов полинома P_n(x) воспользуемся, как обычно соотношением (3). Но применить его напрямую сложно. Поэтому воспользуемся процессом, обратным к процессу понижения степени. Построим вначале P₁(x) - полином первой степени, у которого x₁ является единственным корнем. Затем повысим степень и построим полином второй степени – P₂(x), у которого появляется еще один корень – x₂. Продолжая этот процесс, дойдем до искомого полинома P_n(x). При вычислении коэффициентов нового полинома будем использовать коэффициенты уже посчитанного полинома на единицу меньшей степени. Получающееся в результате соотношения близки к тому, что приведены для случая понижения степени полинома.

Коэффициенты полинома первой степени P₁(x) выписываются явно:

Коэффициенты полинома k-й степени вычисляются через коэффициенты полинома степени k-1:

Переходя к коэффициентам, получим следующие уравнения:

(5)

В соотношении (5) через a' обозначены коэффициенты полинома степени k-1. На самом деле схема безопасна и позволяет считать коэффициенты на том же месте, не требуя дополнительной памяти. Приведу алгоритм вычисления коэффициентов полинома по его корням в виде схемы, приближенной к языку C#.

Дано:

an – коэффициент при старшем члене полинома P_n(x);

n – степень полинома;

x – массив корней полинома (x[0], x[1], …x[n]);

Вычислить:

Массив a – массив коэффициентов полинома (a[0], a[1], …a[n]).

//Вычисляем коэффициенты полинома первой степени

a[1]= 1; a[0] = -x[0];

//цикл по числу полиномов

for(int k=2;k<=n; k++)

{

//Вычисляем коэффициенты полинома степени k

//Вначале старший коэффициент

a[k]= a[k-1];

//затем остальные коэффициенты, кроме последнего

for(int i=k-1;i>0; i--)

{

a[i] = a[i-1]- a[i]*x[k-1];

}

//теперь младший коэффициент

a[0]= -a[0]*x[k-1];

}

//Последний этап – умножение коэффициентов на an

for(int i=0; i<=n; i++)

a[i] = a[i]*an;

Полином Лагранжа

Пусть на плоскости заданы n+1 точка: R₀(x₀, y₀), R₁(x₁, y₁), R₂(x₂, y₂),…, R_n(x_n, y_n). Полиномом Лагранжа P_L(x) называется полином n-й степени, проходящий через все точки R_k. Если точки R_k не образуют возвратов, то такой полином существует и является единственным. Под возвратом понимается ситуация, когда существуют две точки R_i и R_j такие, что x_i = x_j.

Как построить такой полином? Лагранж предложил следующий алгоритм. Полином P_L(x) строится как сумма n+1 полиномов n-й степени:

Каждый из полиномов P_k(x), входящих в сумму, строится следующим образом. Корнями полинома P_k(x) являются все точки R_i за исключением точки R_k. Единственность P_k(x) обеспечивается за счет того, что коэффициент при старшем члене an подбирается так, чтобы полином проходил через точку R_k. В записи Лагранжа полином P_k(x) выглядит следующим образом:

(6)

В записи (6) в числителе находится приведенный полином, построенный по корням, а y_k, деленное на знаменатель в формуле (6), задает an – старший коэффициент полинома.

Условия, накладываемые на полиномы P_k(x), обеспечивают выполнение требований к полиному Лагранжа – сумма полиномов P_k(x) будет полиномом, проходящим через все заданные точки.

Поскольку алгоритм построения приведенного полинома по его корням уже разобран, то схема построения полинома Лагранжа можно выглядеть так:

//Полином Лагранжа определяется как сумма из n+1

//полиномов Pk, для которых известны корни.

for(int k=0; k<=n; k++)

{

//Задание корней для полинома Pk

for(int i =0; i<k; i++)

roots[i] = X[i];

for(int i =k+1; i<=n; i++)

roots[i-1] = X[i];

//Вычисление коэффициентов приведенного полинома по его корням

coefk = CalcCoefFromRoots(roots);

//вычисление An - старшего коэффициента полинома.

An = Y[k] / HornerP(coefk,X[k]);

//Добавление очередного полинома Pk к PL - сумме полиномов

for(int i =0; i<=n; i++)

{

coefL[i]= coefL[i]+An*coefk[i];

}

}

В этой схеме:

X и Y – массивы, задающие декартовы координаты точек, через которые проходит полином Лагранжа,

n – степень полинома,

roots – массив корней приведенного полинома P_k,

coefk - массив его коэффициентов,

An – старший коэффициент полинома, вычисляемый из условия прохождения полинома P_k через точку с координатами X[k], Y[k],

coefL - массив коэффициентов полинома Лагранжа,

HornerP – метод, вычисляющий по схеме Горнера значение полинома по его коэффициентам и значению координаты x,

CalcCoefFromRoots – метод, вычисляющий массив коэффициентов приведенного полинома по его корням.

Сложение и умножение полиномов

При рассмотрении полинома Лагранжа возникала необходимость в нахождении суммы полиномов одинаковой степени, заданных своими коэффициентами. Пусть P(x) и Q(x) – полиномы степени n и m соответственно, заданные своими коэффициентами, и пусть для определенности n >= m. Тогда суммой полиномов называется полином R(x) степени n, коэффициенты которого вычисляются следующим образом:

Пусть полиномы P(x) и Q(x) заданы, подобно полиному Лагранжа, точками, через которые они проходят:

Тогда нетрудно найти подобное представление и для полинома R(x), представляющего сумму полиномов:

В этом случае понадобится вычислить значения полинома Q(x) в n точках.

Если полиномы P(x) и Q(x) заданы своими корнями, то определить корни полинома суммы не удается, более того у суммы вообще может не быть корней. В этом случае для каждого полинома по корням можно вычислить коэффициенты, а затем определить коэффициенты полинома суммы. Можно также рассматривать корни, как частный случай задания множества точек, через которые проходит полином и применить предыдущую схему для определения множества точек, через которые проходит полином суммы.

Рассмотрим теперь операцию умножения полиномов:

Нетрудно понять, что полином S(x) является полиномом степени n+m и имеет n+m+1 коэффициент. Как вычисляется произведение, если заданы полиномы сомножители P(x) и Q(x)? Замечу, что произведение полиномов часто встречается на практике и имеет специальное имя – свертка полиномов.

В отличие от сложения полиномов проще всего найти свертку, если заданы корни обоих полиномов. В этом случае никаких вычислений не требуется, поскольку n корней P(x) и m корней Q(x) будут n+m корнями S(x). Если у полиномов P(x) и Q(x) есть совпадающие корни, то у S(x) появятся кратные корни.

Если исходные полиномы P(x) и Q(x) заданы своими точками, то нетрудно получить набор точек для полинома произведения. Схема во многом похожа на ту, что имеет место при сложении полиномов, заданных точками:

Для получения множества точек, задающих представление полинома S(x), приходится вычислять значение полинома Q(x) в n точках и значение полинома P(x) в m точках, а затем выполнять соответствующее умножение значений двух полиномов.

Если исходные полиномы P(x) и Q(x) заданы своими коэффициентами, то имеем:

Каждый член первой суммы приходится умножать на все члены второй суммы и затем приводить подобные члены при одинаковых степенях x. Нетрудно заметить, что в результате коэффициенты полинома S(x) определяются следующими соотношениями:

Суммирование идет по всем наборам k и r, дающим в сумме значение i. Понятно, что для крайних значений (i=0 и i=n+m) сумма состоит из одного члена, поскольку подобные члены для x в нулевой степени и степени n+m отсутствуют. Число членов суммирования увеличивается при приближении к середине интервала [0, n+m].

Итоги

Подводя некоторые итоги, отметим, что полином можно задать тремя разными способами – его коэффициентами, корнями, точками, через которые проходит полином. Если заданы коэффициенты полинома, то за время, пропорциональное n², (T(n) = O(n²)) можно вычислить значения полинома в n+1 точках. Для вычисления значения полинома в одной точке применяется схема Горнера, выполняющая вычисления за линейное (пропорциональное n) время. Существует и обратное преобразование. Если заданы n+1 точки, через которые проходит полином, то алгоритм Лагранжа позволяет за время O(n²) вычислить коэффициенты полинома. Задача получения коэффициентов полинома по точкам называется задачей интерполяции, а полином Лагранжа называется интерполяционным полиномом.

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

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

Когда полиномы заданы своими коэффициентами, то вычисление значения полинома в заданной точке выполняется по схеме Горнера за линейное время. Сложение полиномов также является легкой операцией и выполняется за линейное время. Свертку полиномов в этом случае выполнить сложнее. Рассмотренный нами алгоритм требует уже квадратичного времени.

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

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

Одно замечание к задаче свертки полиномов. Приведенный алгоритм решения этой задачи для полиномов, заданных своими коэффициентами, требует квадратичного времени. Ввиду практической важности этой задачи много внимания уделялось поиску наиболее эффективного по временной сложности алгоритма. Существуют алгоритмы, решающие эту задачу за время O(n*log(n)). Эти алгоритмы используют технику быстрого преобразования Фурье и обратного к нему. Они сложнее в реализации и требуют работы с комплексными числами или выполнения операций модульной арифметики. Здесь они только упоминаются и детально не рассматриваются.

Задачи

В задачах этого раздела уже не говорится о том, какого типа проект следует строить – консольный или Windows. Предполагается, что обычной практикой является построение Windows приложений.

277. Полином P(x) задан своими коэффициентами. Дан массив координат X. Вычислить, используя схему Горнера, массив значений полинома в точках x_i.

278. Полином P(x) задан своими корнями и старшим коэффициентом an. Дан массив координат X. Вычислить массив значений полинома в точках x_i.

279. (задача интерполяции) Полином P(x) задан координатами n+1 точек, через которые он проходит. Дан массив координат X. Вычислить массив значений полинома в точках x_i.

280. Полином P(x) задан своими корнями и старшим коэффициентом an. Вычислить коэффициенты полинома.

281. (задача построения интерполяционного полинома Лагранжа) Полином P(x) задан координатами n+1 точек, через которые он проходит. Вычислить коэффициенты полинома.

282. Полином P(x) задан своими коэффициентами. Дан массив чисел X. Построить полином Q(X), имеющий своими корнями числа из массива X и корни полинома P(x).

283. Полином P(x) задан своими коэффициентами. Для полинома известны два его корня – x0 и xn. Построить полином Q(x), корни которого совпадают с корнями полинома P(x) за исключением корней x0 и xn.

284. Полиномы P(x) и Q(x) заданы своими корнями и старшими коэффициентами. Вычислить коэффициенты суммы полиномов P(x) и Q(x).

285. Полиномы P(x) и Q(x) заданы своими корнями и старшими коэффициентами. Вычислить коэффициенты произведения полиномов P(x) и Q(x).

286. Полиномы P(x) и Q(x) заданы своими коэффициентами. Вычислить коэффициенты суммы полиномов P(x) и Q(x).

287. Полиномы P(x) и Q(x) заданы своими коэффициентами. Вычислить коэффициенты произведения полиномов P(x) и Q(x).

288. Полиномы P(x) и Q(x) заданы точками, через которые они проходят. Вычислить коэффициенты суммы полиномов P(x) и Q(x).

289. Полиномы P(x) и Q(x) заданы точками, через которые они проходят. Вычислить коэффициенты произведения полиномов P(x) и Q(x).

290. Полином P(x) задан своими коэффициентами. Определить интервал, если он существует, на котором полином имеет хотя бы один корень.

291. Полином P(x) задан точками, через которые он проходит. Определить интервал, если он существует, на котором полином имеет хотя бы один корень.

292. Для полинома P(x), заданного своими коэффициентами, известен интервал, на котором полином имеет хотя бы один корень. Найти корень с заданной точностью, используя схему дихотомии.

293. Построить интерфейс пользователя, позволяющий ему находить корни полинома. В основу поиска положить схему исследования интервала и дихотомии.

294. Построить интерфейс пользователя, позволяющий ему находить корни полинома. В основу поиска положить метод простой итерации.

295. Построить интерфейс пользователя, позволяющий ему находить корни полинома. В основу поиска положить метод Ньютона.

Проект

296. Построить проект, включающий построение класса Polinom и интерфейс пользователя. Методы класса должны реализовать все алгоритмы, рассмотренные в этом разделе. Интерфейс пользователя должен позволять пользователю решать основные задачи, возникающие при работе с полиномами.

Алгоритмы линейной алгебры

Матрицей называется набор чисел, состоящий из m строк и n столбцов. Для программиста матрица – это двумерный массив. Матрица называется квадратной, если m = n и прямоугольной в противном случае. Числа m и n определяют размерность матрицы. Над прямоугольными матрицами определены операции транспонирования, сложения, умножения.

Пусть A – матрица размерности m*n (из m строк и n столбцов) с элементами a_i,j. Транспонированной матрицей B = A^T называют матрицу размерности n*m, элементы которой b_i,j = a_j,i. В транспонированной матрице строки исходной матрицы становятся столбцами.

Операция сложения определена над прямоугольными матрицами одинаковой размерности. Пусть A, B, C – прямоугольные матрицы размерности m*n. Тогда сумма матриц определяется естественным образом:

Операция умножения определена над прямоугольными матрицами, у которых число столбцов первого сомножителя равно числу строк второго сомножителя. Матрица произведения имеет число строк, равное числу строк первого сомножителя, и число столбцов, равное числу столбцов второго сомножителя. Пусть A – матрица размерности m*p, B – размерности p*n, тогда матрица C= A*B имеет размерность m*n. Элементы матрицы произведения определяются как сумма попарных произведений элементов строки первого сомножителя на элементы столбца второго сомножителя.

(7)

Умножение всегда определено для прямой и транспонированной матрицы. Если A – прямоугольная матрица размерности m*n, то всегда определена квадратная матрица B размерности m*m

B = A*A^T = B^T = (A*A^T)^T = (A^T)^T*A^T= A*A^T

Результатом такого произведения является симметричная матрица. Квадратная матрица называется симметричной, если a_i,j = a_j,i для всех i и j, или, что тоже, если A = A^T. Операции транспонирования, сложения и умножения обладают следующими свойствами:

(A^T)^T = A; (A+B)^T = A^T + B^T; (A*B)^T = B^T * A^T

Квадратные матрицы

Квадратная матрица называется диагональной, если все элементы, кроме диагональных, равны нулю, то есть a_i,j = 0 при i /=j.

Квадратная матрица называется единичной, если все элементы, кроме диагональных, равны нулю, а диагональные элементы равны единице, то есть a_i,j = 0 при i /=j и a_i,j = 1 при i = j. Единичная матрица обозначается обычно буквой E, и она играет роль единицы при умножении матриц, поскольку для любой квадратной матрицы A и единичной матрицы E той же размерности имеют место соотношения:

A*E = E*A = A

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

Функция, задающая определитель, обладает рядом важных свойств:

• Определитель диагональной матрицы равен произведению диагональных элементов. Отсюда следует, что определитель матрицы E равен 1;

• Определитель матрицы не меняется при выполнении над матрицей элементарных преобразований. Под элементарной операцией (преобразованием) понимается прибавление к любой строке матрицы линейной комбинации других ее строк. В частности, если к строке матрицы с номером j прибавить строку с номером k (k /= j), умноженную на некоторое число, то определитель матрицы не изменится;

• Если все элементы одной строки матрицы умножить на некоторое число q, то определитель матрицы изменится в q раз (умножается на q);

• Если переставить местами строки j и k, то модуль определителя не изменится, но изменится знак, если разность |k-j| является нечетным числом;

• Определитель произведения матриц равен произведению определителей: D(A*B) = D(A)*D(B)

Не приводя общего формального определения, рассмотрим ниже алгоритм вычисления определителя матрицы, основанный на его свойствах.

Если определитель квадратной матрицы A не равен нулю, то существует обратная матрица, обозначаемая как A^-1. Прямая и обратная матрицы связаны соотношением:

A*A^-1 = A^-1*A = E

Операции транспонирования, умножения и обращения матриц связаны соотношениями:

(A^T)^-1 = (A^-1)^T; (A*B)^-1 = B^-1*A^-1

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

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

В круглых скобках для клеток заданы их размерности. Пусть теперь некоторые клетки нулевые, например, таковыми являются клетки D, F и G. Тогда матрица M2 имеет вид:

Для вычисления матрицы M2 необходимо будет найти произведение трех пар матриц, но значительно меньших размеров, чем исходные матрицы. В целом объем вычислений сократится более чем в три раза.

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

Системы линейных уравнений

Рассмотрим систему из n линейных уравнений с n неизвестными:

(8)

В матричном виде эта система записывается намного элегантнее:

A*x=b (9)

Здесь вектор неизвестных x рассматривается как столбец – прямоугольная матрица размерности n*1. Аналогичный вид имеет вектор правых частей b системы уравнений. В матричном виде условие существования решения системы линейных уравнений (8) и нахождение самого решения формулируется совсем просто. Для существования решения необходимо и достаточно, чтобы определитель матрицы A был отличен от нуля. Тогда у матрицы A существует обратная матрица A^-1. Для нахождения решения системы умножим обе части уравнения (9) на A^-1. Тогда получим:

A^-1*(A*x) = A^-1*b → (A^-1*A)*x = A^-1*b → (E)*x = A^-1*b → x = A^-1*b (10)

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

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

A*X = B

задается соотношением:

X = A^-1*B

Здесь B – прямоугольная матрица размерности n*m, каждый столбец которой представляет вектор правых частей одной системы уравнений. Соответствующий столбец матрицы X дает решение этой системы. Что произойдет, если в качестве матрицы B рассмотреть единичную матрицу? Очевидно, что тогда матрица X будет представлять собой обратную матрицу A^-1. Несмотря на кажущуюся очевидность соотношения A^-1 = A^-1*E, в нем есть определенный смысл, который постараюсь сейчас прояснить. Три задачи – вычисление определителя, решение системы линейных уравнений, нахождение обратной матрицы – имеют одинаковую вычислительную сложность и требуют, если не применять специальные алгоритмы, выполнения порядка n³ операций умножения и сложения. Если посмотреть на соотношение (10), то кажется, что решить систему уравнений несколько сложнее, чем вычислить обратную матрицу, поскольку нужно вначале найти обратную матрицу, а затем умножить ее на вектор правых частей b. Однако реальный алгоритм, рассматриваемый ниже, находящий решение системы, вычислительно проще, чем тот же алгоритм, находящий обратную матрицу. Для такого алгоритма найти обратную матрицу это все равно, что решить n систем линейных уравнений с одной и той же матрицей A в левой части, используя матрицу E в качестве правых частей.

Алгоритм Гаусса

Рассмотрим сейчас алгоритм Гаусса, позволяющий найти решение всех интересующих нас задач – вычислить определитель матрицы, решить m систем линейных уравнений, найти обратную матрицу. Построим вначале расширенную матрицу AB, состоящую из двух клеток:

AB = ||A B||

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

После того, как построена расширенная матрица, вся специфика конкретной задачи теряется, - над расширенной матрицей выполняются одни и те же действия с параллельным вычислением определителя матрицы A. В чем суть этих действий? Над матрицей AB последовательно выполняются элементарные преобразования – деление элементов строки на число, что изменяет величину определителя, и вычитание из одной строки матрицы другой строки, умноженной на некоторое число. Цель наших действий состоит в том, чтобы в расширенной матрице AB клетку A преобразовать в единичную матрицу E. Поскольку каждое элементарное действие можно рассматривать, как умножение слева на некоторую матрицу, то совокупность преобразований, переводящая A в E, эквивалентна умножению слева на матрицу A^-1. Но это означает также, что эти преобразования переводят клетку B в матрицу A^-1*B, что и дает решение исходных задач. Поскольку в результате преобразования A переходит в единичную матрицу, определитель которой известен и равен 1, а для каждого преобразования известно как меняется величина определителя, то параллельно вычисляется и величина определителя исходной матрицы A.

Рассмотрим на простом примере матричный вид элементарных операций. Пусть элементарная операция состоит в том, что к первой строке прибавляется вторая строка, умноженная на число q. Это действие эквивалентно умножению почти единичной матрицы на исходную матрицу:

Матрица, задающая элементарную операцию, отличается от единичной матрицы тем, что у нее в первой строке на втором месте стоит число q, а не ноль. Если бы к первой строке прибавлялась не вторая строка, а строка с номером j, то число q стояло бы не на втором месте, а в позиции j. Если строка j прибавляется не к первой строке, а к строке с номером i, то число q появлялось бы в i-ой строке матрицы.

Рассмотрим теперь возможную реализацию алгоритма Гаусса:

public void Gauss(double[,] M)

{

det = 1;

int n = M.GetLength(0);

int m = M.GetLength(1);

double d =0,r=0;

for (int i = 0; i < n; i++)

{

//Приведение столбца i к единичному вектору

d = M[i, i]; det *= d;

//деление на диагональный элемент: M[i,i]теперь = 1;

for (int k = 0; k < m; k++)

M[i, k] /= d;

//Элементарная операция: сложение строк

for (int j=0; j<n; j++)

{

//К строке j прибавляется строка i, умноженная на r

//В результате M[j,i]=0

if(j!=i)

{

r=-M[j,i];

for (int k = 0; k < m; k++)

M[j, k] += r * M[i, k];

}

}

}

Аргументом метода является расширенная матрица M = ||A B||. В результате работы метода матрица M приобретает вид: ||E A^-1*B||. В зависимости от того, как задана матрица B, находится решение одной системы уравнений, нескольких систем или вычисляется значение обратной матрицы. Параллельно в переменной det формируется значение определителя матрицы A.

На рис. 8 показан возможный интерфейс проекта, построенного для работы с линейными уравнениями.

Рис. 8 Интерфейс проекта, предназначенного для решения задач линейной алгебры

Алгоритм Гаусса в том виде, как он выше рассмотрен, не всегда обеспечивает получение результата. Действительно, пусть в матрице А элемент a[1,1] равен нулю. Тогда при выполнении элементарных операций в процессе преобразования матрицы А к единичной матрице Е возникнет ошибка уже на первом шаге при делении первой строки на элемент a[1,1]. Однако равенство нулю диагонального элемента вовсе не означает, что определитель матрицы равен нулю (если речь не идет о диагональной матрице), или что для нее не существует обратной матрицы.

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

В случае, когда а[i, i] равно нулю, алгоритм ищет первую строку ниже i-й, в которой элемент a[i, j] не равен нулю. Эта строка добавляется к строке i, что гарантирует возможность деления на а[i, i].

Алгоритм с выбором главного элемента в столбце

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

Алгоритм с выбором главного элемента во всей матрице

На каждом шаге приведения очередного столбца к диагональному виду в еще не приведенной матрице ищется максимальный элемент и меняются местами не только строки, но и столбцы матрицы, ставя максимальный элемент в позицию а[i, i]. Этот прием гарантирует отсутствие переполнения при выполнении операции деления. Гарантируется также, что при умножениях не будет получено слишком большое число, поскольку деление на максимальный элемент с последующим умножением на один из элементов приводит к тому, что элементы преобразованной матрицы не увеличиваются по модулю. Однако ничто не дается даром. Выбор главного элемента, перестановка строк и столбцов, необходимость обратной перестановки в конце вычислений, - все это усложняет алгоритм. Как правило, страдает и точность вычислений, особенно для плохо обусловленных матриц. Все модификации алгоритма стоит применять тогда, когда в основной схеме возникла исключительная ситуация, требующая корректировки алгоритма. Обработчик исключительной ситуации при делении на ноль, возникновении переполнения, потери значащих цифр, может вызывать модифицированный вариант алгоритма в надежде получить решение, когда отказывается работать основная схема.

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

Интерполяционный полином, определитель Вандермонда и обусловленность матриц

Вернемся к задаче построения интерполяционного полинома, проходящего через заданное множество точек. Напомню, заданы своими координатами (x₀, y₀), …(x_n, y_n) точки, через которые должен пройти полином степени n. Требуется найти коэффициенты a₀, …a_n этого полинома. Ранее был рассмотрен алгоритм Лагранжа, позволяющий построить этот полином. Нетрудно понять, что существует прямое решение этой задачи, состоящее в построении системы линейных уравнений для нахождения неизвестных коэффициентов полинома. В матричном виде эта система уравнений имеет вид:

Матрица X называется матрицей Вандермонда, а ее определитель соответственно определителем Вандермонда. Этот определитель вычисляется достаточно просто:

Он равен произведению разностей координат точек, где произведение берется по всем j, большим i. Очевидно, что определитель будет отличен от нуля, если у точек нет совпадающих координат x. Этот факт отмечался и при рассмотрении полинома Лагранже, когда говорилось, что множество точек «не имеет возвратов».

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

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

Если матрица А плохо обусловлена, то и обратная к ней также является плохо обусловленной матрицей. Во сколько раз могут возрастать ошибки в элементах прямой матрицы при ее обращении? Примерный ответ на это дают «числа обусловленности» матрицы. Предлагаются различные количественные меры обусловленности матриц. Одной из таких мер является М-число обусловленности Тьюринга:

Матрица Вандермонда – потенциальный кандидат на плохую обусловленность. Если посмотреть на ее структуру, то видно, что для ее элементов во многих случаях характерен большой размах – отношение между максимальным и минимальным элементом велико. Действительно, пусть например максимальная по модулю координата x_n имеет значение 100, а степень полинома n равна 6. Это довольно скромные цифры, но уже в этом случае минимальный элемент матрицы равен 1, а максимальный - 10¹² . Примерно такой же размах будет и у элементов обратной матрицы. Ее максимальный элемент будет примерно равен 1, а минимальный - (max a_i,j)^-1, так что число обусловленности M будет примерно равно 10¹². При наличии небольших ошибок в измерении координат, ошибки в определении значений полинома в точках, отличных от измеряемых, могут многократно возрастать. По этой причине интерполяционный полином еще можно применять внутри интервала измерений, но не рекомендуется использовать в задачах экстраполяции.

Последовательность умножений матриц. Динамическое программирование

На практике зачастую возникает необходимость найти произведение нескольких матриц:

Умножать матрицы можно в произвольном порядке, - результат не зависит от порядка умножения, а вот эффективность выполнения этих умножений может существенно зависеть от порядка. От порядка зависит как время выполнения операций, так и память, требуемая для хранения промежуточных результатов.

Умножение двух матриц размерности [m * p] * [p * n], выполняемое в соответствии с определением (7), требует m*p*n операций умножения элементов (операций сложения, которых требуется примерно столько же, учитывать не будем, поскольку они менее трудоемкие). Рассмотрим в качестве примера умножение трех матриц:

Если выполнять умножение в естественном порядке, то для умножения A₁ на A₂ понадобится выполнить 10⁶ операций умножения. Полученная в результате матрица из 10⁶ элементов будет требовать значительных ресурсов оперативной памяти. Последующее умножение потребует выполнения еще 10⁶ операций умножения, так что суммарно потребуется выполнить 2*10⁶ операций умножения. Если же вначале выполнить умножение A₂ на A₃, а затем умножить A₁ на полученный результат, то суммарно потребуется выполнить 2*10³ операций умножения – в тысячу раз меньше! Заметьте, дополнительной памяти вообще не понадобится, так как результатом умножения A₂ на A₃ будет число (матрица размером 1*1).

Найти оптимальный порядок в случае трех матриц нетрудно. Но уже при 10 сомножителях задача становится нетривиальной. Дело в том, что число возможных расстановок скобок в последовательности A₁* A₂*… A_n растет экспоненциально с ростом n, так что простой перебор становится практически невыполнимым при сравнительно небольших значениях n. Рассмотрим сейчас классический алгоритм решения этой задачи, основанный на идеях динамического программирования.

Нашу задачу сформулируем следующим образом. Дана последовательность перемножаемых матриц - A₀* A₁*… A_n-1. Дана последовательность m₀, m₁, …m_n, определяющая размеры перемножаемых матриц, так что матрица A₀ имеет размеры m₀ * m₁, матрица A₁ имеет размеры m₁ * m₂, матрица A_n-1 имеет размеры m_n-1 * m_n. Матрица, полученная в результате этих перемножений, имеет размеры m₀ * m_n. Требуется найти расстановку скобок в последовательности A₀* A₁*… A_n-1 или, что тоже, найти порядок выполнения умножений матриц, минимизирующий суммарное количество операций умножения элементов матриц.

Идея построения алгоритма состоит в следующем. Вначале применяется типичный прием решения сложных задач «Разделяй и властвуй» – исходная задача разбивается на подзадачи меньшей размерности с той же структурой, затем строится общее решение на основе решения подзадач. В таких случаях часто применяются рекурсивные алгоритмы. Но в данной ситуации обычный рекурсивный алгоритм не работает, поскольку для него число подзадач растет экспоненциально. Одна важная особенность спасает ситуацию– одни и те же подзадачи появляются на разных ветвях рекурсивного спуска, так что число новых задач сравнительно невелико. Это и позволяет применить технику динамического программирования – запоминать решение новой задачи при первом ее появлении и использовать это решение, когда в нем появляется необходимость. Для динамического программирования характерно и то, что решение строится снизу вверх вместо рекурсивного спуска сверху вниз, – вначале решаются простые задачи малой размерности, затем из них строится решение подзадач следующего яруса, подъем вверх продолжается, пока не будет достигнуто решение основной задачи.

Покажем теперь, как эта идея воплощается в жизнь. Пусть mop(i, j) – функция, задающая минимальное число операций умножения элементов на отрезке последовательности A_i* A_i+1*… A_j. Если число операций умножения равно mop(i, j), то это означает, что скобки в последовательности расставлены оптимально. Пусть nop(i, j) – функция, задающая номер, определяющий последнее умножение в оптимальной расстановке скобок. Если nop(i, j) = k, то это означает, что в оптимальной расстановке скобки поставлены так: (A_i* A_i+1*… A_k)*( A_k+1* A_k+2*… A_J). Конечно же, внутри каждой скобки существует своя оптимальная расстановка скобок.

Теперь нетрудно выразить mop(i, j) через значения этой функции для двух подзадач меньшей размерности

(11)

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

Все было бы хорошо, но беда в том, что в соотношении (11) неизвестно значение k. К счастью, у k не так много возможных значений. Поэтому можно исключить k из соотношения (11) ценой увеличения числа подзадач. Вот окончательное выражение для mop(i, j):

(12)

Соотношение (12) позволяет построить конкретный алгоритм решения задачи. Напомню, что наша цель (одна из наших целей) вычислить mop(0, n-1) – минимальное число операций для последовательности матриц A₀* A₁*… A_n-1. Для ее достижения будем постепенно заполнять квадратную матрицу mop ее элементами – mop(i, j), пока не получим искомый элемент - mop(0, n-1). Фактически достаточно построить треугольную матрицу, поскольку элементы ниже диагонали нас не интересуют.

Матрица будет заполняться не по строкам и не по столбцам, а по диагоналям. Заполнение этой матрицы начинается с заполнения главной диагонали, элементы которой представляют простейшие задачи размерности 1, решения которых заранее известны и равны 0 в соответствии с соотношением (12). На следующем шаге будут вычисляться по формуле (12) элементы диагонали, примыкающей к главной, - mop(0, 1), mop(1, 2) и так далее. Эти элементы дают решение задач размерности 2. Так, восходя по диагоналям, доходим до правого угла матрицы, дающего решение общей задачи. Каждый новый элемент диагонали использует при своем вычислении ранее построенные элементы левее и ниже его.

Параллельно с заполнением матрицы mop формируется и матрица nop с элементами nop(i, j). Значением nop(i, j) становится значение k, обеспечивающее минимум mop(i, j).

Поскольку обе матрицы являются треугольными, то можно экономит память и хранить их значения в одной матрице, у которой диагональные элементы равны нулю, выше диагонали хранятся значения mop(i, j), а ниже диагонали nop(i, j). Вот как выглядит эта матрица:

Приведу теперь процедуру, выполняющую вычисление этой матрицы:

public void MopAndNop()

{

maxmop = MaxM();

//заполнение матрицы mop

//главная диагональ заполняется нулями

for (int i = 0; i < n; i++)

mop[i, i] = 0;

//цикл по диагоналям

for (int i = 1; i < n; i++)

//диагональ i заполняется элементами mop

//симметричная диагональ - элементами nop

//цикл по элементам диагонали

for (int j = 0; j + i < n; j++)

EvalMopAndNop(j, j + i, ref mop[j, j + i],

ref mop[j + i, j]);

}

Вначале заполняются нулями элементы главной диагонали – mop(i, i). Затем начинается цикл по боковым диагоналям. Процедура EvalMopAndNop вычисляет два симметричных элемента матрицы – mop(j, j+i) и mop(j+i, j), равное nop(j, j+i), где j – номер элемента на диагонали, а i – номер диагонали.

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