Шабров Н.Н. - Метод конечных элементов в расчётах деталей тепловых двигателей (1061803), страница 13
Текст из файла (страница 13)
Тогда алгоритм ничной нормали к криволинейному определения внешней нормали контуру двухмерного конечного эле- на любои стороне элемента будет один и тот же. Радиус- вектор точек, принадлежащих одной и той же стороне элемента, будет иметь вид Кв = хе, + уе„ (4.54) где х и у даются соотношениями (4.52) и (4.53), р означает номер стороны элемента, р = 1, 2, 3, 4. Нумерация сторон элемента показана на рис.
4.3. Рассмотрим для определенности первую сторону конечного элемента второго порядка. Соблюдение условия о совпадении направления обхода узлов при одномерной интерполяции кривой в соответствии с (4.52) и (4.53) требует, чтобы в качестве координат х, и у, были взяты последовательно координаты точек с но- т 7 мерами узлов 1, 5 и 2. Выражение касательной к контуру для первой стороны будет К,,~ — — х,~е, + у,~е„ (4.55) где $ обозначает дифференцирование по $. Для орта внешней нормали к первой стороне элемента получим й~,~Х (е,Хе,) е,у,~ — е,х,~ (4.56) ~К1Х(едХеа)) у„а 1уа (4.59) р= Х л ($)р (4.62) где р — значения функции р в узловой точке с номером у (у = =с, (,/г). Очевидно, что порядок интерполяциопного полинома в (4.62) может отличаться от порядка соответствующего полинома в (4.52) или (4.53).
Важным моментом при использовании метода численного интегрирования Гаусса является правильный выбор числа точек интегрирования. Этому вопросу уделено достаточно внимания в рабоче (19), где приводятся таблицы для различных вариантов интерполяционпых полиномов в интегралах для одномерных и двухмерных конечных элементов. Так, например, если интеграл (4.60) вычисляется на контуре, представляющем собой кривую 83 Далее на основании (4.56) имеем для первой стороны элемента р, = рп ~ е1 = ру / ~/ х' ~+ у' ~, (4.57) р„= рп1.
е2 = р( — х ~) /~/ х' ~ + у' ~, (4.58) где р — скалярная функция нормального давления. Йля элемента длины дуги кривой, заданной в параметрической форме, можно написать й = ф х-' ~ - ~- у- '~ д.'. Тогда для интеграла (4.51) получим выражение 1 Г~~',~ = рЕ~Л'т ' ' с)Е (1 — с, у, А). (4.60) — 1 В случае простой зависимости нормального давления р на контуре интеграл (4.60) довольно легко может быть вычислен в явном виде.
Однако при сложных зависимостях для р интеграл (4.60) удобно вычислять методом численного интегрирования Гаусса. Применяя к интегралу (4.60) одномерный метод численного интегрирования Гаусса, получим 1= — 1 где $~ — точки интегрирования; Н, — весовые коэффициенты Гаусса. При полиномиальных или близких к ним законах изменения нормального давления р на контуре элемента можно воспользоваться техникой функций формы для интерполяции давления по его значениям в узловых точках. По аналогии с (4.49) можно записать второго порядка при квадратичной зависимости давления р, то наивысшая степень полинома подынтегральной функции будет равна пяти. Для точного вычисления интеграла от полинома пятой степени требуется взять три точки интегрирования Гаусса.
Обсуждение вопросов программирования выражения (4.61) будет сделано позднее, в разделе, посвященном организации данных для решения задачи методом конечных элементов. 4Д. ОСЕСИММЕТРИЧНАЯ ЗАДАЧА Для решения осесимметричных задач с использованием изопараметрических конечных элементов высоких порядков вводится тороидальный элемент, образованный вращением плоского изопараметрического четырехугольного конечного элемента первого или второго порядка вокруг оси симметрии.
Интерполяционные соотношения для таких элементов строятся в плоскости, проходящей через ось вращения, и, очевидно, совпадают с рассмотренными ранее подобными соотношениями (4.32) и (4.33) для плоских четырехугольных элементов с тем лишь отличием, что в качестве глобальных координат теперь выступают цилиндрические координаты г и г. По аналогии с (2.120) напишем выражение для блока матрицы градиентов В тороидального элемента ов д (4.63) д д Ла д где Л'„— функции формы, которые определены в локальной системе координат соотношениями (4.29) — (4.31) или (4.28); а =- = — 1, 2, ..., т; т — число узлов конечного элемента. Связь между производными функций форм по глобальным координатам в (4.63) и производными тех же функций по локальным координатам устанавливается точно так же, как при рассмотрении плоского элемента в предыдущем разделе.
Следовательно, в данном случае справедливы соотношения (4.38) — (4.40), которые нужно отнести к цилиндрической системе координат. Введем обозначения д — Л' дт (4. 64) д 84 и далее Ь„о 0 с„ 0 0 с„д„ (е) Ваа (4.65) Тогда выражение для блока ма-рицы жесткости тороидального изопараметрического конечного элемента будет иметь вид (е) Каа =1л+1в (4.66) где 1я — — 2я ~ ) Во'ч'НоВОг с(е1.1 о$ г(т); †! — 1 (4.67) ! ! 1в = 2я ) ~ (В~оо Н)К~+ В))!,'„) Н)Вф + В!!"„! Н)В~)а)) г г(е1 Я !1$ д)1. — 1 — ! (4.68) Матрицы ВЙ, Но, Н! даются соотношениями (2.123) и (2.124).
Поскольку форма конечного элемента интерполируется по значениям координат узловых точек, то для радиуса г, входящего в выражения (4.67) и (4.68), можно написать гЯ, )1) =- ~ У,,г„, а (4.69) где ㄠ— радиус узловой точки с номером с!. Видно, что подынтегральные функции в (4.67) и (4.68) являются сложными полиномиальными функциями локальных координат, поэтому для вычисления интегралов (4.67) и (4.68) удобно воспользоваться методом численного интегрирования Гаусса: (4.70) 1в =- 2я Е Е (Во"Н)В$ + В)~)' Н)ВО~' + ! .! !.=! + В)!„') Н)В)ач).де1,1 -,,,Л)о,. (4.71) где ~! и )1; — координаты точек интегрирования Гаусса.
Подьштегральные выражения в (4.67) и (4.68) с точностью до скалярной функции совпадают с аналогичными выражениями в (2.126). Это обстоятельство позволяет использовать без изменений подпрограмму ЯТ1ГЛ для вычисления выражений (4.70) и (4.71) в одной точке интегрирования Гаусса. При этом двухмерный массив ВЕ определяется точно так же, как в плоской за- 85 даче в соответствии с (4.45). Для одномерного массива КЫ и параметра ЭЕТ справедливы следующие представления: (4.72) (4.73) )ЭЕТ = где1 3 ~1„.Н;Н;, САРЕ ТАВЕЕ(ХО) РО 400 1О=1,ХО СА(.1.
ОА1)ЗЗ(1О,Х,Н1) РО 400 ЗО=1,ХЙ СА1.1. ОА()ЗЗ(ЗО,У,НЗ) СЛЕ1. РХРХЦХРЕ,Х,У,Х1,У1,ЕХ1) СЛЕЕ РХРХ1(ХРЕ,Х,У,Х1,У1,РХ1) САРЕ РХДЛС(ХРЕ,ХРЕ,УРЕ,РХ1,КЗС) СЛ1Л. М1ХУ (РДС,ХРГьРЕТ,1.1,12) СА1.1 ХХВЕЕ(ХРГьХРЕ РХ1,КЛС) ЙЕ=ЧЛ1БЕ(ХРЕ,ХРЕ,РХ1) ТЕ=ЧА11)Е(ХРЕ,ТРЕ,РХ1) РО 410 1=1,ХРЕ 410 КХ(1)=ГХ1(1)/ЙЕ С1.= РИАЗ э У1 1ХО(ТЕ) !((1+РУЛЯ) е (1 — 2 е Р(ЗЛБ)) СМ=УТЗХО(ТЕ)/(2 э (1+Р()АЬ)) РЕТ=РЕТэ Н1э НЗ АКЕА-АКВА-';РЕТ РЕТ=РЕТ э КЕ Я=(3 ж С1.+2 ж СМ) е А1.ЕА(ТЕ) э ТЕ РО 420 11 =1,ХРЕ Х=- (11. — 1) э ХРР РЕ(Х+1)=ЕЕ(Х+1)+(РХ1(1,11)+ЦХ(11.)) ж Я ж РЕТ РЕ(Х+2) =-.РЕ(Х+2)+ РХ1(2,! Ь) э Я+ РЕТ РО 420 )1=1,11. СА1.1. БТ1РЛ(ХРР,Н.,З(.,С1.,СМ,РХ1,ЙХ,ЗТ,РЕТ) СА1Л.
РАЗЕ(Х РР, Н.,Л 1.,ЗТ,ЗЕ) 420 СОХТ1ХУЕ 400 СОХТ1ХНЕ 86 где значения функций форм, радиуса и определителя матрицы Якоби вычисляются в точке интегрирования Гаусса. Поскольку двухмерный массив ВЕ в осесимметричной задаче имеет точно такую же структуру, как и аналогичный массив в плоской задаче, то программа формирования массива ВЕ для тороидального элемента с четырехугольным сечением будет такой же, как и для соответствующего плоского элемента. Ниже приведен фрагмент программного модуля для вычисления ~ имметричной части матрицы жесткости К(') тороидального элемента с четырехугольным поперечным сечением и размещение ее в одномерном массиве 5Е. Вычисление блока вектора узловых сил, эквивалентных действию распределенного в пределах тороидального элемента температурного поля, сводится к определению интеграла ! 1 Г !3, = 2п ~ ) В,',! Нв!!г г(е1 .1 г$ дт1, (4.74) — ! — ! где е, — матрица-столбец начальных деформаций — дается соотношением (2.131), а для радиуса г справедливо соотношение (4.69).
Матрица-столбец е, содержит в качестве скалярного множителя функцию температуры, для которой, как и прежде, можно использовать интерполяционное соотношение (4.49). Таким образом, подынтегральная функция в (4.74) представляет двухмерный полипом, наивысшая степень которого определяется порядком интерполяционных полипомов для поля перемещений, температуры и радиуса. Выполняя матричные преобразования в (4.74) и применяя одномерный метод численного интегрирования Гаусса, получим Ь„+— (ЗХ+ 2р) ОТ!'!г' де( Л Н;'г1;.