XIII Власова Б.А., Зарубин B.C., Кувыркин Г.Н. Приближенные методы математической физики (1081417), страница 80
Текст из файла (страница 80)
11.2. Задачи теплопроводиооти в твердом теле 581 Отсюда в случае симпдежспых КЭ, используя (10.58), получаем з з )ч. л,„=~~~„. 6„,. ~Ае)м )! е~'~)и)не. )11.22) )=1 1п! ео=! Интеграл в (11.22) не зависит от номера т узла. Поэтому для трехмерной задачи при Х, = 4 при помощи (10.31) находим з з Л" ее ~ С-С-Ь"5')Л" ел л П еб (11.23) в=! 1=1 где через 1', обозначен объем области конечного элемента с номером е, а М, то=! з Л1е1 = Л1 1Ъ' Ч 51е151е) )п е~ 11 ив~ рв! где Л вЂ” среднеарифметическое узловых значений Л в спмв (е) плексном КЭ. Пусть в трехмерном КЭ с номером е одна из граней имеет площадь Я) и аппроксимнрует участок поверхности 5з, на котором задано граничное условие (11.3). Зависимость коэффициента теплообмена ЯР) от координат точки Р Е о", аппроксимируем соотношением Х' является среднеарифметическим узловых значений компоненты Л), тензора Л в симплексном КЭ.
В случае двумерной задачи (Х, = 3) в (11.23) $', следует заменить на площадь Р, КЭ, а для одномерной (Х, = 2) — на его длину 1,. Для изотропного материала вместо (11.23) получим 582 (!. ПРИКЛАДНЫЕ ЗАДА ЧИ где Л!,' — число узлов КЭ, принадлежащих рассматриваемой грани; (1~ ) — значение ЯР) в узле с номером т. Тогда второй интеграл в (11.21) примет вид Так как грань трехмерного симплексного КЭ является двумерным симплексным КЭ, для которого М,' = 3, то, используя (10.28), находим Л(') 311! +(уг' + )3з ~ Л(') 2А' + 2)уг' + дз у =() а остальные значения Л!„получаем соответствующей перестановкой нижних индексов. Для двумерной задачи Л(,' = 2. =(е) ! (е) (е) г — (е) Используя (10,19), вычисляем Лы — — — (3)3 + (1 )1'„Л, = —,(9,' +()г' )1,', где 1,' — длина стороны треугольного симплексного КЭ, на которой заданы граничные условия вида (11.3).
Аналогично можно вычислить интегралы в (11.12). Учет нелинейности, Зависимость Л и 1® от координат точки М Е 1~ в (11.1) и Л, (1 и (г от координат точки Р Е ог в (11.3).позволяет путем итераций учесть и зависимость этих функций от температуры, т.е. приближенно решить нелинейную задачу стационарной теплопроводности. Уточняя узловые значения этих функций по найденным в первом приближении из решения СДАУ (11.20) значениям Т„, Л( = 1, Ж*, темпера; (!) туры в узлах сетки КЭ, можно вычислить уточненные значения элементов матриц Л и („), а затем, решая снова СЛАУ (11.20), получить узловые значения Т(ч во втором приближении и т.д.
(г) Вместе с тем благодаря известному экстремальному свойству функционала (! 1.4) для решения нелинейной задачи применимы и различные модификации методов оптимизации (Х1Ч). Кратко остановимся на одном из наиболее простых вариантов с точки 11дн Эадачи теплопроводноети в твердом теле 583 зрения построения алгоритма решения нелинейных задач, называемом методом локальных вариаций.
Вариация бТ)ч температуры лишь в одном Ю-м узле сетки КЭ относительно принятого в качестве нулевого приближения значения Т, вызывает изменение И(Т~Д лишь той части (о) функционала (11.18), которая непосредственно зависит от узлового значения Тв1. Варьируя только это значение, можно добиться минимума 3(Тм ~1, а затем перейти к варьированию температуры в соседнем узле и т.д. После обхода всех узлов подобную процедуру приходится повторять, так как каждое последуюшее варьирование температуры в соседнем узле может отклонить значение 1(Т))Д от минимального. Поэтому сначала рационально вести варьирование во всех узлах с достаточно грубым шагом ЬТ, а затем, получив совокупность значений Т,, Ж = 1, Ж', первого приближения, для которых (1) переход к температурам Тм = Те1 ~.ЬТ уже не понижает зна- (1) чения 1(Т)ч Д, следует уменьшить шаг варьирования (например, вдвое) и повторять описанную процедуру до тех пор, пока шаг варьирования не станет менее заданной точности вычисления температуры.
Нестацнонарнан теплопроводность. Перейдем к рассмотрению задачи нестационарной теплопроводности в твердом теле, описываемой дифференциальным ураннением (прн 1 > О, М б )') вида (2.53) с(М)~™~) = ЧЯМ)ЧТ((,М))+1к(е)(1,М), (11.24) где с(М) — объемная теплоемкость материала тела, с граничными условиями вида (11.2), (11.3) Т(1, б') = ~,(1, б'), Р Е,Э„ (11.25) )((Р)~УТ(1,Г') (~')+а(1,~ )Т(1, ) =)з(1,~ ), ~'б 51, (11.2б) и начальным условием Т(0, М) = ~о(М), М Е т' ()оз, в момент нремени 1 = О, принимаемый за начало отсчета.
Сформулиро- 584 ГЬ ПРИКЛАДНЫЕ ЗАДА ЧИ ванной задаче не удается поставить в соответствие функционал, для которого ее решение было бы его точкой экстремума*. Поэтому в качестве интегральной формулировки этой задачи п ри мем условия вида (6. 76) (и) ДТ('М) — ~7 (А(М) С Т(1, Л4)) — 7И) (1, М),(М) бУ— д1 ~-~~М~Р~~'ть,Р) ~Р~<-д(Й,Р~т~~,Р~ — яЙ,Рд ~Р)нее=О ортогональности проекций невязки, возникающей при подстановке в (11.24) искомого приближенного решения, на элементы ис, Ь = 1, Жв, базиса Мг.-мерного подпространства Я,у гильбертова пространстпва Я непрерывных в У функций, имеющих в У кусочно непрерывные производные и обращающихся в нуль на Яы Преобразовывая последнее равенство при помощи первой формулы Грина и учитывая, что ис(Р) = О при Р Е.ЧП получаем (обозначения аргументов функций опущены) (Гнжт>)~~., + (Л вЂ” ф),) в + + (РТ Л)ис ПЯ = О, Е = 1, Мв.
(11.27) Приближенное решение рассматриваемой задачи будем искать в виде, аналогичном (11.16): Ти,(1, и) = ф'(М) й(1) = ~~' Т„(1)~Р„(М), (11.26) Мм! где ТВ(1) — зависящие от времени значения температуры в узлах сетки КЭ, а в качестве базиса подпространства М~~ примем систему (ум(М))ч функций (11.17), т.е. базисные н 'Смс Зарубин В.С., Селиванов В.В. 585 ! !.2, Заяачн теалоаровояноетв в твердом теле проекционные функции возьмем одинаковыми, что характерно для процедуры метода Бубнова — Галеркина. В примере 6.15 получена система обыкновенных дифференциальных уравнений (ОДУ) для нахождения зависящих от времени козффицнентов а„(!) в приближенном решении (6.173) рассматриваемой задачи.
Однако отличие (11.28) от (6.173) состоит в том, что при использовании МКЭ нет необходимости в построении функции, которая удовлетворяла бы главному граничному условию (11.25), так как достаточно в (11.28) положить Тч(!) = 1!(1, Рм) = ~!!ч(С) в Ж! узлах Ре Е 5! конечнозлементной сетки, принадлежащих участкам 5! поверхности 5. Это позволяет нз числа базисных функций исключить те функции !рл!, ко торые соответствуют указанным узлам, но сохранить их в представлении (11.28) приближенного решения. Номера Х оставшихся базисных функций !р!ч(м) упорядочим так, что . .!т' = 1 %* где Х" = )Уп — Ю! (ясно, что в частном ..5 =5 имеем !'т' =О и случае отсутствия участков 5!, т.е.
5з = 5, имеем Ж" = Юп). В силу свойств функций формы КЭ имеем рн(Р) = = О, % = 1, Ж, Р Е 5!. Подставляя (11.28) в (11.27) вместо Т. приходим к системе ОДУ (обозначения аргументов функций опущены) ~~!, — | с!!оь!!ол!ей'+ т1Тк Й !ч=! г +~т ! !л~ е' ~ а «-|Ет,т а) = !т=! Г ее! о2 !чс 1~ т~лт+|ьт нл — ~ — 1 е т~й !ч=л '+! Яе Мс з — Е !!!=!!!'+! !, з=! л2 586 1!.
ПРИКЛАДНЫЕ ЗАДАЧИ или, используя обозначения элементов матриц, ~ '(с,~ „'", +л„,т„(()) = (т„(() л(=1 =Яь — ~~!, (Сия „+Льн~1н(()), 1 =1,№ (1129) е(Лл((() Л(=Л( +1 С!.е( = с(М) 1Р1.(М) еР(у(М) е((е, Ь, Ле = 1, Ле~, (11.30) которые являются элементами симметрической матрицы порядка Л1г., называемой влобальной матприцеб тпеплоемкосяпи. Ее (как и глобальную матрицу теплопроводности) удобнее формировать из вкладов отдельных КЭ. Учитывая (11.17) и (11.30), эти значения можно представить в виде Е № )ее с,„=т' /,1м1(Еа)е~ее1м!) (Ее~>„е)1м1)ее= е=11, 1=1 и=\ — С(.) й(е) й( е=! 1=1 е=! где С(„') = с(М) уе~(') (М) !Р(') (М) е(Р'. (11.31) Начальные условия для этой системы составляют значения тл (О) = 1о(МИ), где Мм — узлы сетки КЭ с номерами Л( = =1,№.
Отметим, что элементы Л!.1е, Ь, Ж = 11, Леп, глобальной матрицы теплопроводности Л могут зависеть в (11.29) от времени в силу возможной зависимости от времени коэффициента теплообмена )3 в граничном условии (11.26). Аналогичное замечание относится и к элементам Яь матрицы-столбца Я. В (11.29) входят также значения 587 11.2. Задачи теплопроводпоети в твердом теле Элементы С(„), 1, и = 1, Х„образуют симметрическую лвапзрид(у пзееьдоеммосеви ХЗ с номером е, имеющую порядок Х,.
Элементы этой матрицы распределяют полную теплоемкость КЭ С(') = с(М) 1(1" с(М) ж ~ с~'~~р~'~(М), где с — значение с(М) в узле с номером т, то вместо равенства (11.31) получим Л(е о,'„'= ~'е':|~еи1и)е(е1м)ем1и1а, 1п.зе) еем1 и Отсюда в случае симплексных КЭ для трехмерной задачи при Ж, = 4, используя (10.31), находим 2с1 +2сз +сз +с4 (е) (е) (е) (е) 12— 120 (е) Зс1 +сз +сз + 4 (е) (е) (е) (е) '11 60 а остальные значения С(„) получаем соответствующей перестановкой нижних индексов.
Для двумерной (11(, = 3) и одномерной задач значения С(„') нетрудно вычислить, используя (10.28) и (10.19) соответственно. Систему ОДУ (11.29) запишем в матричном виде (11.33) по его узлам в виде сосредоточенных теплоемких масс таким образом, что изменение внутренней энергии дискретной системы совпадает с ее изменением в КЭ при непрерывных распределениях с(М) и Т((, М), М б р;. Если зависимость с(М) от координат точки М Е Ъ, аппроксимировать в пределах этого КЭ соотношением 588 РЬ ПРИКЛАДНЫЕ ЗАДА ЧИ где С* и Л* — симметрические матрицы порядка М с элементами Сса и Лс,л, А, Ж = 1, Х*, соответственно; В и Я матрицы-столбцы размера Ж х 1 с элементами Тл~(1) и соответственно.