Сегерлинд Л. Дж. - Применение метода конечных элементов (1051193), страница 32
Текст из файла (страница 32)
Однако эти дополнительные затраты машинного времени компенсируются, вероятно, экономией времени в процессе обработки исходных данных. В тех прикладных областях, где градиенты искомых величин не могут быть надлежащим образом аппрокснмированы системой кусочно-постоянных функций, использование элементов высоких порядков позволяет получать более точные результаты. Сокращение времени обработки исходных данных в сочетании с большей точностью результатов расчета служит достаточным основанием для изучения возможностей применения указанных элементов. В этой главе основное внимание будет сосредоточено на обсуждении одномерного элемента.
В первом разделе рассматриваются функции формы, второй раздел иллюстрирует применение квадратичного элемента. Использование естественной системы координат и методы численного интегрирования обсуждаются в третьем и четвертом разделах, причем этн вопросы включены в данную главу потому, что при рассмотрении одномерного элемента упрощается иллюстрация нх реализации. Действительная необходимость в методах численного интегрирования будет пояснена в Элементы высокого нередко. Одномерный элемент следующих двух главах, где рассматриваются треугольный и че- тырехугольный элементы.
Данная глава заканчивается введением в теорию изопараметрических элементов. 13.1. Квадратичные и кубичные элементы Аппроксимирующий полипом в общей форме для одномерного элемента имеет вид ф=а,+аах+аасса+... +а,х'-'е (13.1) где г — число узлов элемента. Рассмотрим элемент с тремя узлами, причем один из его узлов расположен посередине между край- Фиг. 13,1, Одномерный квадратичный элемент. ними точками (фиг. 13.1). Ему соответствует интерполяционный полипом ср=а,+а,х+а,ха, (13.2) где ас, аг и аг определяются из условий (13.3) при х= 1..
1'ак же как в гл. 3, эти условия приводят к алгебраической систе- ме уравнений, решив которую, находим ат=Ф„ 4Фс — ЗФС вЂ” Фа а,= 1 (13.4) а,= —, (Фс — 2Фс+ Ф,). ф=Фс ф Ф1 ф=Ф, при х=О, с, при х= —, Глава 18 Подставив выражение (13.4) в (13.2) и выполнив образования, получим ф=У~Ф~+ ИгФ~+ ЫлФ„, где некоторые пре- (13.5) ( 2х)( х) 4х (1 х) х( 2х) и У =а,+а,х+а,х'+...+а,х'-'. (13.8) Первое свойство уже обсуждалось в предыдущих главах. Второе свойство, хотя и не рассматривалось ранее, очевидным образом следует из того факта, что все аз в (13.1) могут быть записаны в виде линейных комбинаций Фа.
Например, формулы а,= а,Ф, + а~Ф~+ а,Ф,, а =б~Ф~+бгФ~+блФм а =с,Ф~+с~Ф,+свФ„ можно использовать для элемента с тремя узлами. Сравнивая формулы (13.9) н (13.4), можно заметить, что а,=1, а;=ах — О. б,= —, 4 3 = б' б= —— Р= (13.9) ! б = —— 2 с, =-~-, 4 сг — — — — ь-, 2 с = —. Й бв Легко видеть, что эти функции формы удовлетворяют критерию, установленному в гл. 3: l ~, "У„=1, р = 1, 1, К..., г. Только что описанная процедура определения коэффициентов а становится утомительной прн большом числе узлов. Другой способ определения функций формы состоит в использовании формул (13.6) и еще двух дополнительных свойств функции формы: )1 в узле с номером р, (О во всех других узлах, отличных от (), Элементы высокого порядка.
Одномерный элемент узким образом, если в интерполяционном полиноме (13.2) заменить а линейными комбинациями узловых значений Фа, то после надлежащей перегруппировки получим, что функции формы должны быть такого же типа, как аппроксимирующая функция. Напри- р=Х т т М'=1-~~ р-I-'~ Фиг.
13.2. Узловые функции лля квадратичного элемента. мер, функции формы (13.5) могут быть записаны в форме Мр — — ар+бах+срхз, совпадающей по виду с многочленом (13.2). Выражение (13.8) может быть представлено в виде произведения сомножителей типа $1 + ч'за (13.10) или гт т = (тра+ эрах) (фа+ трах) ° ° ° (три-а+ фнх) (13 11) Найти произвольные константы в (13.11) просто, если заметить, что й( =0 во всех узлах, отличных от узла р. Это требование может быть удовлетворено за счет подбора сомножителей в (13.11) таким образом, чтобы каждый из них обращался в нуль в одном из узлов. Совокупность функций ~ь каждая из которых обращается в нуль в определенном узле, показана для квадратичного элемента на фиг.
13.2. Положим х 6=— д ке 1=1 — —, Ь 1=1— х и определим функцию г" а следующим образом: 1а, если бч" р для гт'а, Р = 1, если б=() для У (р фиксировано, б=(, /, й,..., г). Представим 1ч а в виде (13.13р ~а 1 П~а» 6-~ Глава 18 где е л — число узлов элемента, à — произвольная константа, определяемая из условия, что Уз = 1 в узле )). Функция формы У дается формулой У 1Ур= П (13.14) р Еыа Описанная здесь процедура может быть обобщена на двумерные и трехмерные элементы.
Соотношения, аналогичные (13.14), для четырехугольного и треугольного элементов будут даны ниже. а) Лонвйныо 61 Квадрамочныо в) кубооныа хриг. 13.3. Линейный, квадратичный и кубичный одномерные элементы. х х а)у 1 — —. и= —; с ь г=ь ( ех ) ( х ) Чх ( .х ) х ( тх ) ( ь)( эь)( ь)' У у ( эо)( у.)' и„-- ф(ю — ф) (г — — „'), н,-ф(1 — +) (~ — ~ ). 247 Элементы высокого норядко. Одномерный элемент На фиг. 13.3 приведены выражения функций формы для линейного, квадратичного и кубичного элементов в одномерном случае. Пример 126. С помощью формул (13.14) и (13.12) требуется опреде- дить функции формы квадратичного одномерного элемента для узлов ! и /. Начнем с узла Е Имеем р=( и Ха=Х;=О. Функциям Ре 8=й /, й соответствуют выражения Р,=1, так как б=(=р, Р/=/г= (1 — 2х/Е), Ре =/е — — (1 — х/Е).
Вычисляя значения этих функций в точке х=Ха=О, получаем Р/1 в=1 и Ре) о=1. Подстановка в формулу (13.14) дает 11 — (2х/Ц! 11 — (х/Ц1 / 2х 1 / х 1 1 1 что совпадает с (13.5). Повторяя выкладки для центрального узла /, имеем () /, х=Хг=Е/2, х Р= —, 11 — (х/Е)1 4х /1 х ~ (1/2) Е (, Е )' (х/Ь) ! (1/2) 1З 2. Применение квадратичного элемента Элементы высокого порядка применяются так же, как симплекс-элементы, поскольку выбор интерполяционного полинома не связан с исходными дифференциальными уравнениями.
Однако есть смысл рассмотреть применение квадратичного элемента, который обсуждался в предыдущем разделе„с тем, чтобы закрепить "ащи знания, связанные с реализацией метода конечных элементов. В качестве примера рассмотрим одномерную задачу переноса ~сила из гл. 8. Задача состоит в том, чтобы определить распределение температуры по длине стержня, подверженного конвективному теплообмену. Глава 18 Уравнения для элемента, выведенные в гл. 5, имеют внд (й >1(т)=() ), (13.15) где (й< >1= ~(В1г(О)(В)(Ч+ ) 5РП'РУЗ т, т=()у) (т) =(ж, й(, М„) 'Г т (13.15) .Матрица градиентов имеет вид (13.17) С помощью функций формы в (13.5) получаем ~( 4х 3) (4 ах) (4х 1)1(л) Учитывая, что в данном случае [01=[Кхх) и ИР=Адх, интеграл )' (В)'(.О1 (В) (У г запишем в виде ( 4 зх)2' (4х 1 )2 Симметрично Так как мы пользуемся квадратичным элементом вместо линейного элемента, который применялся в гл.
8, все интегралы в (13.15) должны быть вычислены заново. Во внутренних точках элемента температура определяется с помощью матрицы функций формы [й11: Элементы вэиеокого порядка. Одномерный элемент Отсюда имеем 14 — 16 2 — 16 32 — 16 2 — 16 14 (Я(е)) ~кк 6Ь (! 3. 18) где А — площадь поперечного сечения элемента. Конвективная часть [й<'>) дается формулой И,Ут'и'И,У, У,У, У,У, У,У, У,У„ УяУт УнИт УяУя (13.19) где Р— периметр элемента.
Конвективная составляющая вектор- столбца (~<и>) имеет вид Если конвективный теплообмен наблюдается на конце элемента, например в узле 1, то У;=1, У;=Ил=О и поверхностный интег- рал принимает вид 1 ьТ 1 (У)где ьТ А, Бэ О (13.21) где А; — площадь поверхности в узле й Наличие теплообмена в узле 1 сказывается и на матрице теплопроводности 19и>] благода- ря поверхностному интегралу У~У~ УЗУТ И1Ия И,У, И,У, У,Иэ И~Ии ИтИи ИиИи 1 "ни~'~иии-и) ~ Ю.
(13,22) Интегрируя по поверхности, содержащей узел 1, получаем )и~и~'щиэ-л,и ~ (13.23) ~ Ь(И)'(И~ЙБ= Рй ~ Бэ о 1 Ых = а 4 (13.20) 1 Глава 18 Интеграл от теплового потока д идентичен уже вычисленному интегралу (13.21), поэтому можно сразу записать 1 ~~(й()г„3 А, О з1 О где д — заданный поток в узле 1. Объемный интеграл, включающий источник тепла Я, вычисляется также легко: ь Ф, 1 ~ (Щ ЯгМ= АЯ ~ М1 Ых= — 4 .
У е У» 1 (13.25) Пример 127. Нужно определить распределение температуры в стержне кругового сечения, изображенном ниже. Радиус 1см 150'С К задаче 127. Запишем матрицы теплопроводности И) ЗП> 1Ш 30 14 — 16 2 — 16 32 — 16 2 — 16 14 п1 К, А 10 АН1 61.<П и фз>1=ф >)+А~РЬ<з> Применение полученных соотношений иллюстрируется ниже на примере того же стержня кругового сечения, который был рассмотрен в гл. 8 (стр. 139). Элемента( емеокоео порядка Одномернад элемент 2Я )(4атрица [н(а)1 содержит дополнительное слагаемое в связи с тем, что на свободном конце второго элемента тоже происходит тепло- обмеи.
Для вектор столбцов 11())1 и 11(а)1 имеем 1 у(() ~ 2птп~ 5(()ЬТ, 4 6 ,о (У(а)!=у(1>1+А57 0 . 1 После подстановки числовых значений исходных данных получаем 54,8 — 46,2 3,9 — 46,2 142,4 — 46,2 3,9 — 46,2 54,8 [й(а)1= и (1п)! = 1ь(а)1 500 2000 500 54,8 — 46,2 — 46,2 142,4 — 3,9 — 46,2 — 46,2 500 (~(а))=к 2000 . 900 Соотношения, включающие этн элементы в единую дискретную мо- дель, имеют вид первый элемент: (=1, )=2, н=3, второй элемент: 1=3, )=4, й=б. Объединение уравнений, определяющих отдельные элементы, произведем методом прямой жесткости.
В результате придем к системе уравнений 3,9 — 46,2 64,8 54,8 — 46,2 3,9 0 0 — 46,2 142,4 — 46,2 0 0 3,9 — 46,2 109,6 — 46,2 3,9 0 0 — 46,2 142,4 — 46,2 7) 7, 7а 7 7( 500 2000 1000 2000 900 Так как Т,=150'С задана, эта система уравнений должна быть предварительно модифицирована. Модифицированная система уравнений имеет вид Решая эту систему„получаем следующие значения температуры в узловых точках: (Т[~=[150, 80,8 55,8, 46,3 43,51, 'С.
Эти значения хорошо согласуются с аналитическим решением исходной задачи: (Тмм „, „„,[~=[150, 80,9, 55,4, 46,2, 43,3[, 'С. Сделаем несколько замечаний, касающихся матриц квадратичного элемента. Во-первых, поверхностный интеграл, который используется при составлении матрицы теплопроводности ~ 6[У[~[У[Ю, з~ содержит несколько отрицательных коэффициентов, чего не было в случае линейного элемента с двумя узлами. Отрицательные чле- ны встречаются при рассмотрении всех элементов высокого поряд- ка. Во-вторых, результат вычисления поверхностного интеграла 1 ЪТ [/т'[ а= 4 1 перестает быть интуитивно очевидным.