Fletcher-1-rus (1185917), страница 40
Текст из файла (страница 40)
Т(к,з)+у(к,зр)) - ест*(т(км,з3 1 - 2.*У(К,З) + У(КР,З)) ист ист + 2 1 соиттиве 2 СОИТ1)ПТЕ Рис. 6.)2. Распечатка подпрограммы )<ЕЕВ<). лировкой той же задачи (см. $ 6.4). Решение данной задачи с помощью псевдонестапиоиарного метода Ньютона излагается в и. 6.4.1. б.!,4. Квазиньютоновсний метод Другой вариант подхода к решению систем нелинейных алгебраических уравнений, часто называемый квазиньютоновским методом, сводится к замене формулы (6.4) соотношением Ч<а+') = Чко — <о<а)Н(а)Н(а), (6.
20)~ 2 зс 4С 5 С 6 7 В ВС 1О 11 12 1З 14 15 С 16 17 1В 19 20 21 22 С 23 24 25 26 27 2В 29 зо Зз 32 зз Зб 35 Зб 37 зв ЗВ 40 41 бг 4З 44 45 Лб 47 44 49 50 51 52 53 54 С 56 ВОВЙООТ1ИЕ дйсво(о,т,ад) ЕУЙЬОЙТЕВ ТНЕ дйСОВ1ЙИ От ТНЕ 20 ВТЕЙОУ ВОЙОЕЙВ' ЕЯОЙТ10ИВ ЙЕЙЬ*В О,Ч,СХ,СУ,ССХ,ССУ 01НЕИВ10И О(21,21),Ч(21,21),хд(50,50) СОННОМ ОХ,ОУ,ЙЕ,НХ,ИУ,ОТ1 И (ИХ вЂ” 2)*(ИТ - 2)*2 ВО 2 д и 1 И 00 1 Е 1,И 1 Йд(х,д) О. 2 СОИТ1ИОЕ ИХР ИХ вЂ” 1 МТР МУ - 1 СХ 0.5/ОХ СТ 0.510У СХХ 1./ЕЕ/ОХ/Вх СТУ 1.ЯЕ/От/Рт ВО 7 д 2,ИХР ОН 3-1 дР 3+1 00 б Е 2,МУР ЕН Е-1 ХР К+1 НВ = 2*(К-2) + 2*(ИУ-2) а(7-2) ЬО НВ + 1 Ь7 НВ + 2 НЬ НВ " 2*(ИУ - 2) 1Р(НЬ .ЬТ.
0)СОТО 3 Йд(ЬО,НЬ+1) - СХ*О(к,д) — СХХ Йд(ЬЧ,НЬ+г) - Скао(х,д) - СХХ 3 Йд(ьо,ив+1) ~ 2.*(схх+стт) + сх*(О(к,дР)-О(к,дн!) + Втх Йд(ЬО,ИВ+2) СТ*(О(КР,д) - О(хк,д)) Йд(ЬЧ,НВ+1) Сха(Ч(х,дР! — Ч(х,дН)) Йд(ЬЧ,НВ+2) 2.*(СХХ+СУТ) + Ст~(Ч(КР,д) - Ч(КН,д!) + ВТ1 НЙ НВ Ф 2и(ИТ-2) 1Р(НЙ .СТ.
И-2)ООТО 4 Йд(ЬО,НЕ+1! СЕ*О(х,д) - СХХ Йд(ЬЧ,НЕ+2) = СХ*О(х,д) — СХХ 4 НК НВ " 2 1Р(НК .ЬТ. 0)СОТО 5 Йд(ЬО,НК+1) = " Ст*Ч(К,д) — Стт Йд(ЬЧ,НЕ+2! = - Ст*Ч(х,д) — Стт 5 НТ НВ + 2 1Р(НТ .СТ. И-2)СОТО б Йд(ЬО,НУ+1) " СУкт(х,д) - Стт Йд(ЬЧ,НУ+2) Ст*Ч(х,д! - Стт 6 СОИТ1ИОЕ 7 СОИТПП)Е ЙЕТОЙН ЕМВ Рис. б.!3. Распечатка подпрограммы ПАСВ(). Гл. б.
Стационарные задачи где матрица Н<й) представляет собой аппроксимацию величинья (1(й))-<, систематически модифицируемую на каждом итерационном шаге, так что она приближается к (8(й))-< по мере приближения к сходящемуся решению. Модификация матрицы Н(йй на каждом шаге оказывается значительно более экономичным иеатоиз нетеОВ гоя н = 1$ 1тмх 400 180 1 1РИ 20 Вт лоооок+оз кРЗ» лоок-оа як лоок+ог он .15ок+ав ИОЛЭ ИОЛЗ .Оа .ОО 1ЛО .4888 ".0559 .0991 .4863 †.0565 ".0991 .4786 ".0584 ".0991 .аь55 -.оь14 -.о952 .аааг -.0657 -.О992 .оооо .оооо .оооо .0753 .0090 .0012 .1515 .0179 .0023 .2297 .0266 .0034 .иоэ .озаэ .ооа5 Я -.21690+00 -.28500-01 -.ИЗВВ+Оа я- -.Ваотп-ог †.ио50-ог -.в4$зп-ог Я -.3259Р-03 ".4288пе04 -.32890-03 я -.12640-04 -.16620-05 -.12770-оа я -.491ап-оь †.ьатоп-от †.49770-06 я -.19эоп-от -.25590-08 -.19810-07 и5 гткялтгонз тнк .9990 .9586 .9990 .9605 .9990 .9603 .9989 .9577 лэзв л52а .оооо .оооо .1317 .1282 .2679 .2605 .4142 .4017 .5774 .5577 Рис.
6.14. Типовая выдача результатов по программе ХЕФТВ1). процессом, чем факторизация матрицы 9(й). Уравнение (6.26) может быть переписано в форме 7<<(й+и — ()г(й) 03(й)В»(й) (6.21) где вектор 6(й) = Н(й)К(й) может рассматриваться как определяющий направление поиска. Скаляр О)(й) выбирают так, чтобы. величина гсг+," оказывалась минимальной в направлении поиска й(й). При больших <т' последнее свойство обеспечивает зна- пе пк пк= Вк пе те те те тк те гт- о 1Т» 20 1Т 40 ХТ 60 гт во гт-го о лттке и и и и и т т» т .Т» т» .9990 .9990 .9990 .9989 .998$ .оооо Лээт .2679 .4142 .5774 ямв= янз янз янз янз ЯН5 .9586 .9583 Л572 .9553 .9524 .оооо .1277 .2598 .аозо .5577 .14960+00 .28540-01 .55040-02 .10240-02 .18920"03 .34990-04 ЯНЗ ЯБ51ВПЛЬ 15 .48$$ -.0559 .4835 †.0567 .474$ -.О5$5 .4625 ".0615 .4462 -.0657 .оооо .оооо .074$ .0090 .1507 .0179 Л290 .0266 .3109 .0349 .98726Р-05 -.0991 ".0991 -.0991 -.0992 -.0992 .оооо .оогг л$02Э .0034 .0045 $6.2, Прнмые методы длн линейных систем 233 чительио больший радиус сходимости, чем при применении обычного метода Ньютона.
Характерный алгоритм, известный как формула ВРО5 [Р)е1сЬег, 1980), позволяет строить сле- дующее приближение Ч<'е» с помощью формулы (6.21). В ре- зультате имеем у<м х((л+» я(м Р<" = - (о(мЫ(' В(">, й(м = 1/Р(м гу<е<, П(~+» — П( < Я(м(П( <у(мр( < + Р(му(мгПав)+ еи(м (1 + [)(му(ю гЯ(му(м) Р(мр(м г (6.22) 5 6.2. Прямые методы для линейных систем Для эллиптических задач, описываемых линейными уравнениями типа уравнений потенциального течения или тех, которые соответствуют промежуточным этапам метода Ньютона, необходимо строить решение линейной системы алгебраических уравнений (6.23) АЧ=В относительно компонентов вектора Ч, где все элементы матрицы А и вектора В известны. Метод исключения по Гауссу [1)аЬ1(1»1з1, В)огсз 1974) является предпочтительным средством Решения уравнений (6.23), однако эффективность реализации Несмотря на то что матрица Л является, как правило, разреженной, структура выражений (6.22) свидетельствует о том, что матрица П оказывается плотной.
При й = 0 матрица Я(о< полагается равной единичной матрице 1. Эффективность квазиньютоновских методов зачастую зависит от того, обладает ли матрица 1 некоторыми специальными свойствами, такими, как положительная определенность [Левш!пцз, 1977а). Следовательно, применимость этих методов к задачам вычислительной гидроаэродинамики должна исследо<ваться для каждого случая отдельно. Значительная часть литературы по квазиньютоновским методам связывает их приме(нение с безусловной минимизацией (см., например, [Роч<е!1, 1976; Р!е1сЬег, 1980) ).
Обсуждение квазиньютоновских методов <в применении к специальным задачам можно найти в работах [Вгоудеп, 1965; Ог1еца, КЬе(пЬ»161, 1970; ЯЬаппо, 1983]. Характерный вариант применения квазнньютоновского метода к задаче гидроаэродинамнки дается в работе [Епце1шап е1 а1., 3981) . Гл. 6. Стационарные задачи этого метода зависит от структуры матрицы А. Заслуживают выделения три категории подобного рода матриц: (а) Матрица А содержит лишь несколько или не содержит совсем нулевых элементов.
Тогда говорят, что матрица А полная или плотнал. (Ь) Матрица А содержит много нулевых элементов. Такую матрицу называют разреженной. (с) Матрица А содержит много нулевых элементов, причем ее ненулевые элементы группируются вблизи главной диагонали. Тогда матрицу А называют разреженной и ленточной. Компьютерные программы, основанные на предположении о том, что матрица А является полной, оказываются более универсальными, но и более дорогостоящими в вычислительном смысле. Такие программы могут использоваться в применении к разреженным и/или ленточным матрицам, однако для последних категорий существуют специальные и притом более экономичные процедуры.
Какой бы ни была структура матрицы А, метод исключения по Гауссу должен применяться в два этапа. Сначала матрица А факторизуется, превращаясь в нижнетреугольную матрицу Е и верхнетреугольную матрицу 11 в качестве двух сомножителей. После этого факторизованная форма уравнения (6.23) решается как (П/= Е 'В. (6.24г Вследствие структуры матрицы 11 процесс решения (6.24) сводится к повторным подстановкам, осуществляемым после формирования комплекса 1 В.
6.2.1. ЕАСТ/$0(.ЧЕ: решение систем с плотными матрицами Если матрица А плотная, то все ее элементы участвуют в образовании матриц Е и 11. Подпрограмма РАСТ (рис. 6.15) реализует требуемую факторизацию по формуле Ш = А. Последующая серия подстановок согласно (6.24) осуществляется с помощью подпрограммы 801 ЧЕ (рис. 6.16). На каждом шаге факторизации множество элементов из текущего ряда й вычитается из элементов последующего ряда 1 с целью исключения оь При этом вводится множитель аь з/аз, ы а величина ам з называется ведущим элементом. Если ведущий элемент равен нулю или очень мал, то исключение по Гауссу либо не проходит вообще, либо дает очень неточные результаты. Во избежание этого вводится частичный отбор ведущего элемента.
Для исключения оз проводится поиск в й-м столбце в интервале 1= й+ 1, й/, ставящий целью найти 1 2 С зс 4С 5 С б С 7 в 9 С 1О с 11 С 12 13 14 С 35 С 16 С 17 1В 19 го гз 22 гз гб 25 С 26 С гт с гв 29 С зо с Зз С зг зз Зб Зб С 36 С ЗЧ С зв зз 4О 41 42 43 44 45 46 47 4$ 49 5О 51 6 6.2.
Прямые методы для линейных систем зпвлоптзне Глст(И,А,$РУТ) ГАСТОА1$Е$ А 1МТО РЕАНОТЕВ Ь.о $0 ТИАТ РЕКИ~А Ь"О ЗРЧТРК) ОТЧЕВ 1МВЕХ ОГ КТН Р1ЧОТ АОН $ЕТ5 БЕРУТ(н) = "1 1Г КЕАО Р1ЧОТ ОССОЛВ 01ИЕН$10И А(50,50),ВРУТ(50) ИН=М" 1 слп$$ еынгиАт10м м1ти Рлнтзль Р1чот1ио 00 5 Х 1.ИИ ГР К+1 $ЕЬЕСТ Р1ЧОТ Ь Х ВО 1 1 КР,М 1Г(АВ$(А(1,К)) .ОТ. Анв(А(Ь,Х]))Ь ~ 1 1 СОИТТМПЕ АРЧУ(Х] $ = А(Ь,Х] А(Ь,Х! Л(Х,Х] А(Х,К) = $ СНЕСК ГОЛ ЕЕЛО Р1ЧОТ 1Т(АВВ($] .ЬТ. 1.0Е-15)ООТО 6 САЬСПЬЛТЕ )Н)ЬТ1РЫЕЯВ 00 2 1 = ГР.И А(1,К] = -А(1,К)/$ 2 сомтзное 1МТЕАСНЛМПЕ АНВ ЕЬТН1МАТЕ ВТ СОЬШ(ИЕ 00 4 д КР,И $ = А().!) А(ь,г) А(к,а] А(н,д) $ 1ГЬАВ$($) .ЬТ. 1.0Е-15]ООТО 4 0031 ХР,М Л!1,3) Л(1,2) + А(1,Х)*$ 3 СОИТ1ИПЕ 4 СОМТ1И1)Е 5 СОИТ1ИПЕ АЕТПАМ 6 $РЧЧ(н) -1 НЕТПАИ ЕИВ Рис.