Шабров Н.Н. - Метод конечных элементов в расчётах деталей тепловых двигателей (1061803), страница 10
Текст из файла (страница 10)
РО 20 К=1,1ЧРМ 20  — — В+ВЕ(К,11.) з СК(К) з ВЕ(К,) 1.) з РЕТ К= 1,+) 1. ЯЕ(К) = ВЕ(К)+ Э 10 СОЧТ1ЫБЕ Параметрам 1ЧРЕ, 1ЧВМ, 0ЕТ присваиваются значения соответственно 3, 2 и А('), 1ЧРŠ— число узлов конечного элемента, двухмерный массив ВЕ был определен ранее в (2.70).
Очевидно, что эта часть программы должна быть окаймлена внешним циклом по конечным элементам. Здесь также используется суммирование с накоплением результата в массиве ЯЕ, необходимость которого будет видна позже, при рассмотрении конечных элементов высокого порядка. Матрица поглощения Р)(е) в соответствии с (3.15) будет иметь вид Ю; Л~,Ю, Л,Ж, Л,Л» Рй =- Л'у Р(( Л(1 Л(») Й "() = И;)(7( 1(71Ы1 7()1Л(» Р до).
л (е) Л'» л(е) У»Л( 7(7»Ц И»Л» (3.31) Заменим в (3.31) функции формы Л-координатами и воспользуемся при интегрировании формулой (2.55). Значение коэффициента Я примем постоянным. Получим 2 1 1 Р"= 1 2 1 12 1 1 2 (3.32) Для матрицы конвекции Р»(' в соответствии с (3.14) и для функций форм, выраженных Ь-координатами, можно записать: л Р» = ~ 1.з (У.) й У.з1й Й. (е) ~(') С, Для определенности предположим, что конвекция осуществляется на стороне, для которой Ь, =- О и значение 6 постоянно. Тогда (3.33) будет иметь вид О О О )( у(е) ~-з~-з ~з~з Й Й = .2'(е) О г.з~.з ~-з(-з 1» О О О О 2 1 .
(3.34) О ! 2 Для векторов узловых тепловых сил в предположении постоянства в элементе коэффициентов /т, Т, Ч, 9 можно получить в соответствии с (3.1б) — (3.18) следующие выражения: Е.т Г /) =- 7.2 л(~) / з 1 Я([и) = 1 ~)А (') 1 (3.35) О О ~(е) Ез Ч [1=,,/' 1 /-з 1 (3.3б) у(Е) /'л О йТ .У('„) О Ез ЬТ й — 1 Е.з 1 (3.37) у(Е) /'/~ Выражения (3.36) и (3,37) могут быть обобщены в предположении линейной зависимости для д, /т, Т на стороне элемента.
Пусть для теплового потока () на стороне элемента, для которой у 1., = О, имеет место интерполяционное соотношение (рис. 3.2) [ ~'2 1'з1 Тогда для вектора узловых теп- Р . З.~. Интериолннии граничных (е) условий теилообмена на контуре ловых сил Г~( получим двухмерного конечного элемента О О О О О О Г~~,'~ = Ез [О У.з Е.з1Ж д; = '~ О 2 1 д/ . (3.38) / (~) д, О 1 2 д, /е Далее пусть для определенности на той же стороне элемента задана зависимость для Ь и Т ' следующими интерполяционпыми соотношениями (рис. 3.2): Тогда для вектора узловых тепловых сил Рз можно записать (е) 1.' 1.г1.з 1 (-з(-г )-з 1 Затем получим 1/(; Л),) с у(е) ))е ()), ))з) ,~ (е) О Й,-(ЗТ,+Т,)-;6,(Т,+Т„) !2 Ь (Т + Т„) +)г„(Т + 3Т„) (3.39) Э.З.
СООТНОШЕНИЯ МКЭ ДЛЯ ТЕТРАЭДАЛЬНОГО КОНЕЧНОГО ЭЛЕМЕНТА В ТРЕХМЕРНОЙ ЗАДАЧЕ ТЕПЛОПРОВОДНОСТИ Интерполяционное соотношение для температуры Т(') в трехмерном симплекс-элементе будет иметь вид Т(') =- [Ж( И, Жз Ж(1ч(е), (3.40) где Т; Т; ч(') —.
== аьо11. Т„ Т Л) а=( ' ~з й Функции формы „(, 1,, даются соотношениями (2.43). Матрица градиентов в данном случае образуется в соответствии с (3.11) следующим образом: д дх (ж; М) 1Ч .Ч,)=(В, В В В,). В(е) = введем обозначение По аналогии с (3.25) (е) в„ (3.41) е г еге-3 1 г с йг ~гЕ,з 1.зЕг д ду д дх И 'И Подставим в (3.41) выражение для Уа нз (2.43) и получим Ва (е) с . (3.42) Матрица теплопроводности трехмерного симплекс-элемента имеет размерность (4 х 4) и в соответствии с (3.14) и (3,42) запишется так: р В; в' Н [В( В; Ви В()(1п, с р(е) ) (3.43) где Н представлено в (3.7). Рассмотрим подробнее блок матрицы Р~'иа, который по понятным причинам также имеет размерность (1 х 1). Имеем Рис.
З.З. Организация хранения матрицы теплопроводности трехмерного симплекс-элемента Рси'а = ~ В,",) НВ~Р'~ г)о р (е) (а, р = (, ), Й, 1). Далее на основании (3.7) и (3.42) получим (3.44) Выражение (3.44) легко программируется. Б соответствии с принятой стратегией программирования будем вычислять только симметричную часть матрицы и хранить ее в виде одномерного л(ассива (рис. 3.3). Если определить двухмерный массив ВЕ так же, как в (2.97), а одномерный массив С дополнить одним элементом, содержащим коэффициент Х„то для вычисления симметричной части матрицы теплопроводности и сортировки ее элементов в массив ЯЕ годится без изменения фрагмент программы, приведенной на стр.
Б2. Параметры ЯРЕ, ХЭМ и 1:)ЕТ должны соответствовать рассматриваемой задаче и конечному элементу, т. е. Х1)М = 3; БЕТ =- Ге). МРŠ—. 4; 3 шабров и. н. д дх д ду д дг (ети+ ЬаХ+ СаУ+ б(ир) (е) (е) 1 1 Рси)з = ()вхЬиЬа +. ) мсис() — )- ЛАиг1))) (е) (е) Для матрицы поглощения Р)(() в предположении постоянства коэффициента 1~ получим в соответствии с (3.15) 2 1 1 1 1 2 1 1 1 1 2 1 1 1 1 2 Е1 г(е) Е2 р1е(е) ! 1 2 3 411~ (3.45) Ез Для матрицы конвекции Р~3' и векторов узловых тепловых сил Го'~, Г ), Г») в предположении постоянства коэффициентов Я, д, Ь, Т па стороне, для которой Е, = О, получим на основании (3.14) и (3.16) — (3.18) и формул интегрирования (2.55) и (2.59) О О О О О 2 1 1 О 1 2 1 О 1 1 2 еА(е) 1О Е2 ЕЗ Е4) е" (1(е) = 12 (е) Рь (3.46) Ез А„ 1Е, Е2 1, 1, (~1е(е) Я сЬ = (3.47) к(е) О Ч 4(е) () (10) = з (3.48) Е3 А (е) 1 О )~у я(е) !)Т (1а = (3.49) Ез А( е) 1 Формулы (3.46) — (3.49) могут быть обобщены на случай более сложной зависимости подынтегральных функций.
Однако соответствующие выражения получаются громоздкими. Для примера рассмотрим вычисление вектора Г(' в предположении линейной зависимости д на грани А(е), т. е. 66 Тогда -0 -0- ~г 0 00--0- 2 1 1 ()) 1 2 1 1 1 2 д( (3.50) 1(е) ! 12 (е) Рд [О У,г Ха У41 й) ~з 4 '.(. 0 ЗА. СООТНОШЕНИЯ МКЗ ДЛЯ ТОРОИДАЛЬНОГО КОНЕЧНОГО ЗЛЕМЕНТА В ОСЕСИММЕТРИЧНОЙ ЗАДАЧЕ ТЕПЛОПРОВОДНОСТИ Рассмотрим функционал (3.8) в цилиндрической системе координат (г, (р, г) для осесимметричных задач теории теплопроводности. Тогда элемент объема (1с и элемент поверхности й можно представить (1() = г (1(р (1г (1г; Ж = г (1(() (11. (3.51) Подставляя (3.51) в (3.8) и интегрируя по угловой координате, получим т =- 2л ( — „~ ата(1 Те Н дга(1 Тг (1г йг — — 1 ЯТ'г (1г (1г— ( 2 3 2 А А — 1~тгд дг+ 1дт дн+ 1дт, ( — т — т )гй).
(352) А с с„ Функционал (3.52) можно рассматривать как функционал для двухмерной плоской задачи теплопроводности -с модифицированными значениями коэффициентов теплопроводности и теплоотдачи и заданных функций источников теплоты, теплового потока и температуры среды. Я=в; (*( .= Ьг. (3,53) 61 Таким образом, преобразования (3.53) в функционале позволяют использовать уже рассмотренную ранее технику МКЭ для решения двухмерной плоской задачи теплопроводности.
С физической точки зрения преобразования типа (3.53) соответствуют переходу от осесимметричной конструкции к пластине переменной толщины, которая изменяется по закону 2пг (рис. 3.4). Поскольку интерполяционные соотношения для температуры и градиента температуры для осесимметричной задачи формально так ие же, как и для плоской задачи, сразу перейдем к рассмотрению выражения элемента матрицы теплопроводности.
Имеем Р(.)„= 2п ~ В()'НВ()г,(где. (3.54) А(е) Интеграл от радиуса в (3.55) по треугольнику в предположении зависимости для радиуса (2.127) будет иметь вид ~ г(1» й = »(',)А" = — (г(+г). +гл) Аьо (3.56) А(е) Окончательно для (3.55) получим Р(')р = 2пг") (Х,У~Ьа + "л,с ср) 7П 1(а) (а, 'р = (, (, й, 1). (3.57) Сопоставление выражений (3.29) и (3.57) позволяет заключить, что фрагмент программы вычисления матрицы элементов теплопроводности и их формирования в массив БЕ, приведенной па стр.
62, пригоден совершенно без изменения и в данном случае. Значение параметра ВЕТ для осесимметричной задачи теплопроводности будет равно 2пг('А(', т. е. объему тороидального симплекс-элемента. Выражение для матрицы поглощения Р("') при условии постоянства по элементу значения К и зависимости (2.127) для радиуса будет следующим: Рис. 3.4. Интерпретация осесим метринной области пластиной пе ременной толщины Подставляя в (3.54) выражения для блоков матрицы градиентов В('), Вр(') из (3.26) с соответствующей заменой декартовых координат цилиндрическими, получим Рсар = 2Й (~г()а()() + ~асаср) — ( ° ~ г (1» ((а. (3.55) А(е) Р)( = 2)т (а (Е) Е~ Лз]К»с]»с1г = (е) А(е) Е,а 6г, + 2», + 2»а 2»; + 2г; + гл 2»; + г; + 2»а ~д(Е) = 2л — 2»;+2»;+гл 2г; (-6»,+2»л г;+ 2»;+2»а 2г; + г~+ 2»л г; — 2г;+ 2»л 2»,:+ 2»~+ 6»а (3.58) При выводе выражения для матрицы конвекции и векторов тепловых сил Р(," и Ра(' предположим опять для определенности, что конвекция и тепловой поток д действуют на стороне элемента, для которой Е, = О.
Получим в соответствии с (3.14), (3.16)— (3.18) и при условии зависимости (2.127) для радиуса и постоян- 68 ства в пределах элемента функций (;), (), Ь, Т следующие формулы: 0 2 ейск'(') 1» (О Л~ Е,з) (11= лз 0 0 0 3«)+ «л 0 «,+«л 0 «)+«л «« —, 3«л (3.59) ее)е(Е) )з 2 1 1 1 2 1 1 1 2 О 0 0 0 2 1 0 1 2 2~~е» 4(е) «(1«((г = 12 0 2„ч Г(.) «й= е.з "о = 2з( (3.60) А(е) Е,"' = 2» (3.61) «) «л у(е) !'з Вектор тепловых сил получается из (3.61) заменой д на ЙТ Глава 4 ИЗОПАРАМЕТРИЧЕСКИЕ КОНЕЧНЫЕ ЭЛЕМЕНТЫ ВЫСОКОГО ПОРЯДКА В ТЕОРИИ УПРУГОСТИ 4.1. ОДНОМЕРНАЯ ЗАДАЧА В предыдущих главах подробно изложена процедура МКЭ на примере симплекс-элементов, т.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
