Дульнев Г.Н., Парфенов В.Г., Сигалов А.В. Применение ЭВМ для решения задач теплообмена (1185899), страница 13
Текст из файла (страница 13)
(2.! 6! Значение рм+'> соответствуег точке, в которой касательная к кривой в точке рИ> пересекает ось р (рис. 2.3). В результате корень уравнения находится как предел последовательности, вычисляемой по формуле (2.16). Можно доказать, что если 1(р) = О, ~' (р) чь О, а ~п непрерывна, то около корня р суЩествует интервал, при попадании в который начального приближения р<"> метод сходится.
Поэтому основная трудность реализации метода Ньютона состоит в выборе начального приближения р<'>. Обычно этот выбор производится с помощью какого-нибудь безусловно сходящегося алгоритма, который может быть основан, например, на методе половинного деления. После определения р<~> выполняется переход на ньютоновские итерации, имеющие более быструю сходимость. На рис. 2.4 показано, что для нашей функции )'(р) — (2.15) ите. рации по методу Ньютона всегда сходятся при выборе начальной точки левее корня р„, но могут привести к выходу из интервала при вы.
боре р"„' правее корня. Если нахождение производной!' (р) в методе Ньютона затруднено, то можно воспользоваться его модификацией — мелюдом секущих. В этом алгоритме начинают с двух приближений к корню — р(»> и р<". Функция 1(р) заменяется прямой, проходящей через точки Ряс. 2.5 Рис. 2Л 1(реп) и((р<')), т. е.
вместо касательной в методе Ньютона строится секущая (рис. 2.5). Тогда по двум приближениям р<"-"> и р<мследующее приближение к корню находится по формуле р(»+ < > — р(») (2.! 7) ( (»)) 11(Р" ") — 1(и(">)1!(Р(» "— Р(»>) в которой по сравнению с формулой(2.16) производная заменена конечной разностью на отрезке (р(» '), р(»'!. Метод простой итерации. В заключение остановимся на самом простом методе решения нелинейных уравнений, который так и называется — метод простой итерации.
Для применения этого метода уравнение (2.!4) представляется в виде Р =<2(1(). (2.18) Очевидно, что это всегда можно сделать, если ввести, например, р(р) =1(р)+р. Итерационная формула метода имеет вид р(»+ < > ~ (р(»>) (2.19) Ясно, что если такой алгоритм сходится, то он сходится к корню уравнения. Для сходимости требуется, чтобы в окрестности корня, в которой выполняются итерации, соблюдалось неравенство (<р' (р) ~ ( Рис.
2.7 Рис. 2.6 ( 1. На рис. 2.6 показан пример сходящегося, а на рис. 2.7 — расходящегося процесса простой итерации. Отметим, что переход от уравнения (2.14) к уравнению (2.18) может производиться различными способами и соответственно будут получаться разные схемы простой итерации.
Например, в нашем случае уравнение (2.15) можно представить в виде (2.18) следующими способами: !) р= В]с1др, <р' (р) =- — В1/з!п'р; 2) р=р — р/В!+с!ар, ср'(р)=1 — —.— —, 1 ' 1 В! мни р ' 3) р=р+р/В] — с!яр, <р'(р) = — 1+ — +— 1 1 В1 и1п' Р Итерационный процесс для варианта 3 будет всегда расходиться, так как <р' (р) ) 1. Варианты 1 и 2 при некоторых числах В1 в некоторой окрестности точки (2п — 1)и/2 могут давать сходящееся решение.
$ 2.3. численнОе интегвивование Задачи вычисления определенных интегралов возникают в расчетах по многим аналитическим решениям. Приведем ряд примеровПри использовании метода конечных интегральных преобразова. ний ]3] приходится вычислять интегралы вида где д, (х), о, (х) — пространственные распределения объемной и поверхностной плотностей теплового потока; 1„(х) — собственные функции краевой задачи. Решения задач с переменными воздействиями д, (т), д, (т) с по. мощью теоремы Дюамеля [13) представляются в виде т [(х, т) = — ( 7(х, д,(т — г'), д,(т — т'), т')от', т д о где 1 — решение задачи с постоянными воздействиями.
Решения задач в полуограниченном теле с помощью интеграль. ного преобразования Фурье представляются в виде несобственных интегралов 1(х, у)= — [ 0(у, в)созвхдв, Г2 Г о где 0 (у, в) — косинус-преобразование Фурье (изображение) функции ) (х, у). Часто в приведенных интегралах аналитическое выражение для первообразной найти не удается, даже если подынтегральная функция содержит элементарные функции, а во многих решениях под интеграломвстречаются специальные функции (например, функции Бесселя). В этих случаях приходится прибегать к численному интегрированию.
В данном параграфе рассмотрим методы вычисления одномерных интегралов, аметоды вычисления кратных интегралов будут изложены в главе 6 на примере задач нахождения угловых коэффициентов излучения. Итак, ставится задача приближенного расчета определенного интеграла ь (=- ) 7(х)дх. а Эти расчеты проводят по так называемым квадрагпурным формулам, в которые входят значения подынтегральной функции 1 (х) в конечном числе узлов, расположенных на отрезке [а, Ь[. Простейшая из квадратурных формул — формула прямоугольников. При ее построении отрезок интегрирования [а, Ь[ разбивается на У элементарных интервалов [х;, х;+,[, ( = 1, ..., У, х, =- = а, хи+, — — Ь, и на каждом из них функция ) (х) заменяется постоянным значением (рис.
2.8) 1(х) ж ) [(х~+ х;+,)!2[. Тогда на 1-м интервале разбиения можно определить приближенное значение интеграла 1;л+, по формуле к~+ 11 у+о= ) ~(х)дхж~( )(хьы ха). 2 58 Суммируя по всем элементарным интервалам при равномерном разбиении с шагом )! = х,„, — к,,получаем следующую приближенную формулу: При более общем подходе на отрезке (а, Ь) выбирается (У + 1) узел (х!),".+,', и приближенную формулу для вычисления интеграла (2.20) записывают в виде «!+ ! !л = ~ с!1(х!). (2.22) с= ! Задача построения квадратурной формулы состоит в выборе расположения узлов х; и значений весовых коэффициентов со Довольно часто функция 1 (х) бывает задана в табличном виде и поэтому узлы х! оказываются фиксированными. Рассмотрим общий подход к определению коэффициентов с! при заданных узлах.
Общий интервал интегрирования (а, Ь) разбивается на совокупность интервалов, ограниченных узловыми точками (рис. 2.9). Хотя в принципе интервалы разбиения могут содержать т'!'к! Ч г:!'«21 к. кэч «,=а".» к„...»;- Ь=»„„ Рис. 2.9 Рис. 2.8 "/+Эй !'+ «$ )'г — — ) !Р,„(х)дх= ч~; ад<р (х,), »=! 59 разное количество узлов, на практике их число делают одинаковым.
В /-м интервале разбиения !х;, х;+ ), содержащем т + 1 точку, подынтегральная функция 1 (х) заменяется интерполяционным поли- номом степени и — сэ (х), значения которого в узлах х;, х;+„... ..., хэ+ совпадают со значениями подынтегральной функции. Далее вычисляют интеграл от интерполяционного полинома по 1-му интервалу разбиения: а затем для расчета интеграла (2.20) суммируют 11 по всем интервалам: М+! 1н =~р,11 — — ~~~ с1)(х1). 1 к=и Перейдем к рассмотрению наиболее употребительных квадратурных формул. Формула трапеций. При ее получении функция 1(х) интерполируется на каждом элементарном отрезке [х1, х1+,] линейной функцией (рис.
2.[0). Тогда 1+~ 1 11= ~ 1(х) с[х ж — [)(х1) +) (х1+,)] (х1+,— х1), 2 к1 и после суммирования по всем интервалам [х1, х1+,] придем к формуле 1л — — — ](а)+ ~~~~1(х1) + — [(Ь))! х,=а, хя.ь~ =Ь. (2.23) 2 2 формула Симпсона. Она основана на использовании отрезка разбиения с тремя узлами и интерполяции подынтегральной функции полиномом второй степени по трем равноотстоящим узлам х;, х1+,, б!х 21 ' 1+2 2;+1 Х Х1,1 Х+а Рис.
2.11 х 1+1 Рис. 2.10 х1+ . На отрезке [х1, х1.+ 2Ь] функция 1 (х) заменяется параболой, проходящей через точки 11, ~1+„~1+2 (рис. 2. [!). Интеграл от интерполяциоиного полинома на отрезке [х1, х1 + 2Л] равен 1+2 а 11= ~ ~(х)дхж — Д1+4~1.и+11„.2).
3 к Считая, что полное число интервалов й! = (Ь вЂ” а)/Ь четное, и суммируя по всем парам интервалов, получим формулу Симпсона В/2 х12 1я == 12+4 ~~' 122+2 ~ ~22-2+)я+~ й13 (224) Увеличивая число узлов л! в интервалах разбиения и применяя для интерполяции подынтегральной функции полиномы более высоких степеней, можно получить еще ряд других квадратурных фор. мул. Однако на практике они применяются редко.
Теперь рассмотрим случай, когда узлы х, не фиксированы, и таким образом в квадратурной формуле (2.22) можно выбирать не только весовые коэффициенты с!, но и расположение узлов хь, в которых вычисляется подынтегральная функция. Использование этих дополнительных «степеней свободы» позволяет повысить точность квадратурных формул. Формулы Гаусса. Если на отрезке интегрирования в качестве искомых параметров квадратурной формулы рассматривать (А/+ !) коэффициент с, и (А/ + 1) узел х!, то получим (2А/ + 2) неизвестных. Эти параметры можно выбрать так, чтобы квадратурная формула бы.ла точна для любого многочлена степени не выше (2А/ + 1).
Решение такой задачи известно: расположение узлов х; вычисляется с помощью корней полиномов Лежандра. Узлы х, и весовые коэффициенты с! для различных У приведены в [2). Построенные таким образом квадратурные формулы называют формулами Гаусса или формулами наивысшей алгебраической точности. Для гладких функций эти формулы дают очень высокую точность. Оценка погрешности численного интегрирования. Различают два вида оценок: априорные и апостериорные. Априорную оценку получают заранее, до проведения расчетов, на основе теоретического анализа квадратурной формулы. Апостериорную оценку определяют после вычислений на основе сопоставления результатов расчетов, проведенных при разных числах отрезков разбиения.