Вычислительная математика лекции. Конспект лекций. Санкт петербург 2011 2 Оглавление
Скачать 3.86 Mb.
|
3.7.2.B сплайны. Идею B сплайнов (базисных сплайнов) поясним на примере простейшего линейного сплайна. Используя набор B сплайнов можно сформировать кусочно-линейный интерполяционный сплайн согласно формуле S 1 (x)= 1 0 n i i i y B . Здесь B 1i базисный сплайн в узле i. B 1i =0, если x≤x i-1 и x≥x i+1 B 1i = 1 1 i i i x x x x , если x i-1 ≤ x < x i B 1i = 1 1 i i i x x x x , если x i ≤ x x i+1 33 B 1i =1, если x=x i B 1i X 1 X i X i-1 X i+1 Разработаны способы построения многомерных базисных сплайнов произвольной степени. Базисные сплайны находят широкое применение при решении краевых задач (в том числе дифференциальных уравнений в частных производных). Примером может служить метод конечных элементов. 4. Среднеквадратичное приближение 4.1. Постановка задачи для f(x). Требуется найти приближающий полином P m (x) , минимизирующий функцию расстояния между f(x) и P m (x): ρ(f(x),P m (x)) → min; 0 ( ) ( ) m m i i i P x x φ i (x) – линейно независимые базисные функции. В качестве функции расстояния для среднеквадратического приближения используется следующий функционал 34 2 0 1 ( , ,...., ) ( ( ) ( )) b m m a J P x f x dx Функционал – отображение из некоторого пространства функций или векторов на числовую ось. 4.2. Постановка задачи для приближения векторов Для векторов задача среднеквадратичного приближения ставится следующим образом. Пусть имеется вектор у R (n+1) и задана линейная комбинация 0 m i i i линейно независимых векторов φ i с неопределенными коэффициентами α i . Каждый из φ i имеет размер n+1. Нужно заменить вектор y на 0 n i i j . Параметры α i нужно найти таким образом, чтобы достичь минимума функционала: 0 2 1 0 0 ( , ,..., ) ( ( ) ) n m m k i i k k i J y Задача среднеквадратичного приближения векторов возникает, если таблично заданную функцию y k (x k ) (k=0,1,…,n) нужно аппроксимировать алгебраическим полиномом P m (x) = α 0 +α 1 x+…+α m x m степени m ≤ n. В этом случае имеем СЛАУ: α 0 x 0 0 +α 1 x 0 1 +…+α i x 0 i +…+α m x 0 m = y 0 α 0 x 1 0 +α 1 x 1 1 +…+ α i x 1 i +…+α m x 1 m = y 1 …………………………… α 0 x i 0 +α 1 x i 1 +…+ α i x i i +…+α m x i m = y i …………………………… α 0 x n 0 +α 1 x n 1 +…+ α i x n i +…+α m x n m = y n Эквивалентная форма записи этой системы: 35 α 0 φ 0 +α 1 φ 1 +…+ α i φ i +…+α m φ m = y , где φ i = (x 0 i x 1 i …x n i ) T , i = 0…m. Таким образом, в данном случае мы аппроксимируем заданный вектор y = (y 0 y 1 …y n ) T размера n+1 набором из m+1 базисных векторов φ i размера n+1 φ i = (x 0 i x 1 i …x n i ) T , i = 0…m. В частном случае m = n, получаем задачу интерполирования. При m < n имеем дело с переопределенной системой. 4.3.Определения. Определение1. Будем называть скалярным произведением функций или векторов следующие числа: для функций v(x), y(x) x [a,b] 1 0 , , ( , ) n n i i i v y R v y v y ( , ) ( ) ( ) b a y v y x v x dx для векторов v,y R n+1 0 ( , ) n i i i v y v y Воспользовавшись ведёнными обозначениями, соответствующие функции расстояния можно записать следующим образом ((y-v), (y-v)), где y -приближаемый вектор или функция, а v – приближающий вектор или функция. Пространство вещественных функций (векторов), для которых определено скалярное произведение, называется евклидовым. Определение2. Пусть в евклидовом пространстве функций (векторов) дана система базисных элементов φ 0 , φ 1 , …, φ m . Эта система называется линейно независимой, если равенство 36 α 0 φ 0 + α 1 φ 1 + … + α m φ m = 0 имеет место тогда и только тогда, когда все α i = 0 для i = 0, … , m. Определение3. Определитель матрицы вида 0 0 0 0 ( , ) ( , ) ( , ) ( , ) m m m m называется определителем Грама для системы функций или векторов. Определитель Грама равен нулю тогда и только тогда, когда система базовых элементов φ 0 , φ 1 , …, φ m линейно зависима. Определение 4. Будем называть функцию ||x|| нормой вектора или функции, отображающей вектор или функцию на числовую ось и удовлетворяющую следующим условиям: 1) Неотрицательность ||x|| ≥ 0; ||x|| = 0 ⟺ x= 0. 2) Однородность ||αx|| = |α| ||x||, где α скаляр. 3) Неравенство треугольника ||x+y|| ||x|| + ||y||. Примеры норм. Для векторов: 2 1 0 0 | |; ||x|| ; ||x|| max | | . n n i i i i i i x x x x Для функций: 2 1 2 [ , ] || ( ) || | ( ) | ; || ( ) || ( ( )) ; || ( ) || max | ( ) | . b b a b a a f x f x dx f x f x dx f x f x 37 4.4.Решение задачи. Для того, чтобы найти min(J(α 0 , α 1 , …, α m )) требуется решить систему m+1 линейных алгебраических уравнений для m+1 неизвестных α 0 , …, α m : 0,( 0, ) i J i m Такая система называется нормальной. Можно показать, что она имеет вид: 0 0 0 0 0 0 ( , ) ( , ) ( , ) ( , ) ( , ) ( , ) m m m m m m f f Здесь f приближаемая функция или вектор. Действительно, пусть осуществляют среднеквадратическое приближение функции f(x). Тогда 2 [ ( ( ) ( )) ] 2 ( ( ) ( )) ( ( ) ( )) 2 ( ( ) ( )) ( ) 0. b b m m m i i i a a b m i a J P x f x dx P x f x P x f x dx P x f x x dx Отсюда следует (φ i , φ 0 )α 0 +(φ i , φ 1 )α 1 +…+(φ i , φ m )α m -(f, φ i )=0, i=0,1,…,m. Приведенное доказательство легко модифицируется на случай среднеквадратичного приближения заданного вектора. 38 5. Ортогонализация 5.1. Постановка задачи При отсутствии специального способа выбора базисных элементов φ i , (i = 0…m ) даже в случае их линейной независимости обычно уже при m > 5 нормальная система оказывается настолько плохо обусловленной, что её решение на ЭВМ практически бесполезно. В частности, такими свойствами при больших m обладает система базисных функций 1, x, … x m . Дело в том, что, будучи формально линейно независимой, с ростом m система базисных элементов может оказаться очень близкой к линейно зависимой. Максимально линейно независимой оказывается система базисных элементов φ 0 , φ 1 , …, φ m , ортогональных на множестве точек x 0 , x 1 , …, x n для векторов или на отрезке [a,b] для функций. Система базисных элементов называется ортогональной, если (φ i ,φ j ) = 0 при i ≠ j и (φ i ,φ j ) ≠ 0 при i = j. 5.2. Процедура построения системы ортогональных базисных элементов. Пусть задана система базисных элементов φ i i = 0, 1, …, m. Нужно построить систему соответствующих ортогональных базисных элементов ψ i . Последовательно рассчитывают: ψ 0 = φ 0 ; ψ 1 = φ 1 + а 10 ψ 0 ; ψ 2 = φ 2 + а 20 ψ 0 + а 21 ψ 1 ; ψ j = φ j + a j0 ψ 0 + a j1 ψ 1 + ... + a jj-1 ψ j-1 ; ψ m = φ m + a m0 ψ 0 + a m1 ψ 1 + … +a mm-1 ψ m-1 39 Коэффициенты а ji i < j, i = 0,…,m отыскиваются из условия ортогональности ( ψ i ,ψ j ) = 0. Можно показать, что a ji = -(ψ i ,φ j )/ (ψ i ,ψ i ). Например, из условия ортогональности (ψ 0 ,ψ 1 ) = 0 следует (ψ 0 ,ψ 1 ) = (ψ 0 ,φ 1 +а 10 ψ 0 ,ψ 0 ) = (ψ 0 ,φ 1 ) + (а 10 ψ 0 ,ψ 0 )= (ψ 0 ,φ 1 ) + а 10 (ψ 0 ,ψ 0 )=0. Отсюда а 10 = -(ψ 0 ,φ 1 ) / (ψ 0 ,ψ 0 ). При использовании ортогональной системы базисных элементов матрица нормальной системы становится диагональной Дополнительные замечания. 1.Если в качестве базисных функций задают алгебраические полиномы, для построения ортогональной системы можно использовать полиномы Лежандра как удобную альтернативу методу Грама-Шмидта. Полиномы Лежандра ортогональны на отрезке [-1,1] . Для их определения используется рекуррентная формула (n+1) ℒ n+1 (t)-(2n ℒ n (t)t+n ℒ n-1 (t , где ℒ 0 (t , ℒ 1 (t)=t. Переход от отрезку x [a,b] к отрезку t [-1,1] осуществляется заменой переменных 2 ; t= 2 2 2 a b b a a b x t x b a 2 Нормальная система имеет единственное решение, если базовые функции линейно независимы При таком условии можно говорить о существовании и единственности обобщенного полинома, удовлетворяющего критерию среднеквадратического приближения Разработаны приемы построения обобщенного полинома в случае линейной зависимости базисных функций В этом случае решение не единственно 3 Нелинейная задача метода наименьших квадратов предполагает, что приближающая функция g(x,α 0 ,…, α m не линейна относительно искомых коэффициентов В типовых 40 программах встречаются функции следующего вида: ae bx +ce dx экспоненциальная аппроксимация, Q m (x)/P n (x рациональная дробь, ax b потенциальная функция 6. Равномерное приближение функций. Постановка задачи. Пусть f(x) заданная на отрезке [a,b] непрерывная функция. Будем говорить, что многочлен P n (x) приближает f(x) равномерно на отрезке [a,b] с точностью ε, если ||f(x) - P n (x)|| ∞ = [ , ] | ( ) ( ) | max n x a b f x P x Естественно поставить задачу следующим образом: среди всех многочленов фиксированной степени n найти многочлен Q n (x), для которого ε минимальна. Q n (x) называют многочленом наилучшего равномерного приближения. Теорема 1. Для любой непрерывной на отрезке [a,b] функции f многочлен наилучшего равномерного приближения Q n (x) существует и единственен. Теорема 2.(Теорема Чебышева): для того, чтобы многочлен Q n (x) был многочленом наилучшего равномерного приближения непрерывной на отрезке [a,b] функции f(x), необходимо и достаточно, чтобы на [a,b] нашлись по крайней мере n+2 точки x 0 < x 1 <… x n+1 такие, что [ , ] ( ) ( ) ( 1) | ( ) ( )|, 0,1,..., 1 max i i n i a b f x Q x n f x Q x i n . Здесь α – постоянная, равная 1 или -1 для всех i одновременно. Точки x 0 , x 1 ,…, x n+1 , удовлетворяющие условиям теоремы, называют точками чебышевского альтернанса (alternare – 41 чередоваться). Теорема накладывает на точки чебышевского альтернанса два требования: 1) В точках x i погрешность |f(x i )-Q n (x i )|= [ , ] | ( ) ( )| max a b n f x Q x 2) Для всех i=0,1,…,n+1 погрешность f(x i ) - Q n (x i ) меняет знак при переходе от точки x i к следующей точке x i+1 Следовательно, внутри диапазона функция погрешности имеет n+1 нулей. Следующий пример подтверждает справедливость последней теоремы. Рассмотрим частную задачу наилучшего равномерного приближения функции f(x) = x n+1 многочленом n-ого порядка Q n (x) на отрезке [-1.1] . Заметим, что погрешность R n+1 =x n+1 –Q n (x) есть многочлен n+1-ой степени. Задачу можно переформулировать так: среди всех многочленов степени n+1 с фиксированным старшим коэффициентом, равным 1, найти многочлен R n+1 (х) такой, для которого величина max|R n+1 (x)|, x [-1,1] была минимальной. Таким многочленом является нормированный многочлен Чебышева T n+1 (x)/2 n , а решением поставленной задачи многочлен Q n (x) = x n+1 - T n+1 (x)/2 n . С другой стороны, легко доказать, что R n+1 (х) удовлетворяет условиям теоремы Чебышева. Действительно, T n+1 (x m ) = (-1) m , где x m = cos(mπ)/(n+1). Здесь x m – точки максимумов многочлена Чебышева, m = 0,1,…,n+1. Максимум модуля погрешности max|R n+1 (x)| на отрезке [-1,1] достигнут в n+2 точках и равен 1/2 n . При этом знаки погрешности чередуются при переходе к каждой следующей точке. Таким образом, Q n (x) удовлетворяет условиям теоремы Чебышева. Дополнительные замечания. В большинстве случаев задача о наилучшем равномерном приближении является очень трудной. Для ряда частных случаев разработаны специальные численные методы, реализованные в типовых программах. 42 Часто достаточно ограничится нахождением многочлена, близкого к наилучшему, либо просто найти многочлен, равномерно приближающий f(x) с заданной точностью ε. 7. ЧИСЛЕННОЕ ДИФФЕРЕНЦИРОВАНИЕ Постановка задачи. Пусть заданы значения функции в узлах ( ) i i f f x , где 0 , 0, 1, 2, , i x x ih i h-постоянный шаг изменения аргумента. Будем полагать [ , ] ( ) n a b f x C . Требуется определить 0 ( ) f x Решение задачи. Используем соотношение 0 0 0 0 ( ) ( ) ( ) lim x f x x f x f x x . Для приближенного решения можно положить 1 0 0 0 ( ) f f f f x h . Для оценки погрешности рассмотрим степенное разложение 2 0 0 0 ( ) ( ) ( ) ( ) 2 h f x h f x h f x f , из которого получим соотношения 1 0 0 f f f r h , (1) где 0 0 ( ), [ , ] 2 h r f x x h − погрешность (1). Соотношение (1) имеет первый порядок точности, определяемый степенью h в остаточном члене r. Формула (1) не единственна. Используя разложения 2 3 1 0 0 0 0 ( ) ( ) 2 6 h h f f h f x f f , 2 3 1 0 0 0 0 ( ) ( ) 2 6 h h f f h f x f f , можно получить 43 2 (3) 1 1 0 ( ) 2 6 f f h f f h (2) Формула (2) имеет второй порядок точности. Складывая приведенные выше степенные разложения, получим формулу для вычисления производной второго порядка 2 (2) (4) 1 0 1 0 2 2 ( ); [ , ]. 12 f f f h f f h (3) Для нахождения производных любого порядка можно строить формулы любого порядка точности. Основной прием заключается в замене функции, заданной в (n+1)-ом узлах интерполяционным многочленом p n (x) и принятии ( ) ( ) ( ) ( ) m m n f x p x (4) |