Шабров Н.Н. - Метод конечных элементов в расчётах деталей тепловых двигателей, страница 2

В матричной форме уравнение (1.2) запишется так Результатом перемножения матриц в (1.8) будет Ооозначим И~е1 ~х1" ~ х1 '1<'. Ъ'(е1 гх х")1~1 м (1.9) Тогда (1.б) запшпется так: д (х) у(~)д(~] + уапц(о (1. 1О) Точки интерполяции х~ ~ и х~ ~ называются узловыми точкома или узлами. Матрица И~о называется матрицей функций фарлафы. Функции И," и Л'," называются функциями формы. Важно понять, что функции формы предназначены для интерполяции по выбранной форме (в нашем примере линейной) решения при помощи узловых значений. Аналогичные рассуждения можно провести применительно к подынтервалу Лх„~, для которого узловые точки, узловые значения и функции формы имеют индекс (е — 1) (рнс.

1.1, в). Таким образом, узловыми точками х~' и функциями формы Ж„" (а — ~', 1) определен одномерный линейный конечный элемент. Линейность конечного элемента означает линейную зависимость функций формы относительно х. Формула (1.10) дает решение задачи в пределах конечного элемента при условии, что узловые значения у ~ (а = 1, 1) известны. 1.2.

ВЫВОД РАЗРЕШАЮЩИХ УРАВНЕНИЙ МКЭ и, (1.11) и векторы узловых неизвестных элементов: Рассмотрим задачу на осевое растяжение стержня постоянного сечения под действием сил собственного веса. Точное решение этой задачи представлено на рис. 1.2. Будем решать эту задачу методом конечных элементов. Разобьем стержень на два конечных элемента одинаковой длины и пронумеруем узловые точки так, как это показано на рис. 1.3. В качестве неизвестных узловых значений выберем перемещения и„(а = 1, 2, 3) центров тяжести соответствующих сечений стержня по оси х. Введем глобальный вектор узловых неизвестных всей систсмы Введем матрицу кинематических связей конечных элементов таким образом, чтобы имело место равенство В нашем примере равенство (1.13) имеет вид и, 1 0 0 из 0 1 0 из 0 1 О иа 0 0 1 (1.14) нз Рис.

!.2. Зпюры продольной силы, напряжений и перемещений в задаче о растяжении стержня под действием собственного веса; у = 6Гг; в = Ог; и = (4(е — ге) 0Г1е Матрицу кинематических связей а можно разделить на подматрицы а(') таким образом, чтобы т) (') — а(е) 11. (1. 15) Рис. 1.3. Схематизация стержня конечными эле- ментами Полная потенциальная энергия рассматриваемой системы запишется так: 2( з( П = — ~ ео(1х — з ~ би(х)(1х, 2 (!.1б) где е — продольная деформация; о — нормальное напряжение; и (х) — продольное перемещение точек стержня; б — объемная сила; з — площадь поперечного сечения стержня. Выразим е и а, входящие в (1.16), через перемещения. На основе (1,10) имеем (е) (Х1 е) (е) (е) ~ е) (е)П(е) ( -Г- ) ) (1.17) 9 Далее з(е) .

(1и(е) (х)((с1х), или Йх (1.18) В (1.18) дифференцированию подлежат только функции формы, поскольку и; ) и и)'~ являются числами. Дифференцирование в (1.18) в соответствии с (1.9) дает (е) 1 и (е) [ 1 и(е) 1(Е) ' 1(О ! или в матричной форме в(е) ),(е) , е *~(е) 1(е) (1.19) Введем в рассмотрение матрицу градиентов В(') конечного элемента В(') = — — Х(')(х) =- [ — 1 1) — „, = [В,") 8)(').

(1.20) Тогда деформация з(') и напряжение а(') в элементе запишутся так: е(е) В(е) у(е) . о(е) — Ез(') = ЕВ(е) ч(е) (1.21) где П() = — ~ е(е)а() с1х — з ~ 6и(е) (х)([х (г = 1, 2). 2 ((е) ((е) С учетом (1.15) и (1.21) выражение для П(') запишется так: Пво = ~3'а(')' — [ В(')'ЕВ( ) (1ха(') $/ — $3'а(')'з [ бй(')'([х. 2 )(е) ) (е) (1.22) Введем обозначение: К(е) =з ~ В()'ЕВ(~) ([х.

(1.23) (() 10 где Š— модуль Юнга. Соотношения (1.21) показывают, что дсформации и напряжения в элементе определяются через узловые значения перемещений при помощи матрицы градиентов. Вернемся к выражению полной потенциальной энергии П. Разобьем величину П па две части, каждая из которых соответствует полной потенциальной энергии одного конечного элемента.

Итак, П= П')+ П"', Матрица К(е) называется матриуей жесткости рассматриваемого конечного элемента. С учетом соотношения (1.20) получим для К(') следующее выражение: К(е) з ~ Е( — 1 1! (,)е ((х г- — 1 ! 1 ((е) или в итоге к (1.24) (1.25) Наша цель состоит в том, чтобы найти вектор узловых сил конечного элемента Г('), который статически эквивалентен объемным силам 6, т. е. П(е) = ч(')'Г(е) =з ~ (7и(') (х)(1х. 6 (() Подставим в (1.28) выражение для и" (х) из (1.17) и получим ч(е)'Г(е) = )Га(е)'з ~ дЧ(е) ()х. (1.27) ( (е! Из (1.27) следует выражение для вектора узловых сил Г(') Г(') = з ~ Х(')' (х) 6(1х. (1.28) ((е) Подставляя в (1.28) выражения для функций форм (1.9), получим х х,".

Г (е) ~е (е) ( (е) (1.29) (1х. !! Матрица жесткости К(') конечного элемента является симметричной матрицей. Этот факт не вызывает удивления, поскольку потенциальная энергия деформации представляет собой симметричную, положительно определенную квадратичную форму. Выясним физический смысл элемента 1(,'~) матрицы К' '.

Для этого закрепим узел 1, а узлу 1 придадим единичное перемещение. Тогда элемент К;; матрицы жесткости численно равен (е) реакции в направлении 1-й обобщенной координаты от единичного перемещения в направлении 1пй оообщенной координаты. Введем в рассмотрение вектор внешних сил Г(') конечного элемента. Второе слагаемое в (1.22) представляет работу объемных сил. Рассмотрим подробнее это слагаемое, Выведем обозначение П(е) = — я ) 6и(е) (х)(1х.

((е) Выполняя интегрирование в (1.29) по длине каждого конечного элемента, получим г( ) Я01 (1.3О) Теперь выражение для П(') можно записать следующим образом; П(') =. — 11'а("')'К(')а(е) $1 — $)'а(')'Г('). 2 (1.31) Подставим выражения для П(') (е = 1, 2) из (1.31) в вариационное уравнение Лагранжа (1,1) и последовательно выполним необходимые преобразования. Будем иметь 6П = 6 ~ П (') = О. Далее (1.32) где о1) означает произвольные вариации узловых неизвестных. В силу производительности последних следует — "УП() =О.

д13 (1.33) Дифференцирование н сложение в (1.33) дает ~ а(е)'К(е)а(е) 11 ~~~ а(е) 'г (е) () (1.34) Введем следующие обозначения: — ~~~~ а(е)'К(оа(е) ° е (1.35) а(е) ' г'(е) е (1.36) где К вЂ” глобальная матрица жесткости, а à — глобальный век- тор узловых внешних сил вссй системы. Тогда с учетом введенных обозначений уравнение (1.34) запишется (1.37) Уравнение (1.37) представляет собой систему линейных алгебраических уравнений относительно неизвестных узловых значений. Решение 11 системы уравнений (1.37) используется для получения вектора узловых значений элемента ~(') и затем по вектору ч(') с помощью матрицы функции форм получаем решение и(') (х) в пределах каждого конечного элемента. 12 1.3.

ФОРМИРОВАНИЕ МАТРИЦЫ~(ЖЕСТКОСТИ И ВЕКТОРА НАГРУЭКИ СИСТЕМЫ УРАВНЕНИЙ МКЭ Проиллюстрируем процесс формирования глобальных матрицы К и вектора Г на основе формул (1.35) и (1.36), который сводится к следующим матричным преобразованиям: Г 1 2 а ())' Г (() а(2) ' 0 0 +бз1 1 О О 1 Г (2) 1 2 1 2 1 0 = 6з1 0 О 1 — 2 1 2 (1.38) К(() а()) К (() а(') а(')' К (') а (') 1 — 1 0 ьЕ 0 О О а(')' К"'а" 1 0 0 1 — '; (1.40) 0 0 К (2) а(2) К (2) а (2) 0 — 1 1 1 — 1 1 0 0 1 а(2) К(2)а(2) 0 0 0 — 0 1 — 1 0 — 1 1 (1.42) 1 — 1 0 Х а(е)'К(е) а(е) 1 1+ 1 е Π— 1 1 (1.43) Результатом преобразований (1.38) является получение глобального вектора внешних сил Г, а результатом преобразований (1.39) — (1.42) является получение глобальной матрицы жесткости К всей конструкции. При подробном анализе рассмотренных выше преобразований можно заметить, что матричные операции умножения удается заменить логическими операциями пересылки компонент векторов Г(') и матриц К(') в глобальные вектор Г и матрицу К всей конструкции.

Остановимся на этом подробнее. !3 Существует две системы нумерации узловых значений — локальная и глобальная. Локальная нумерация устанавливает локальные номера узловых значений в пределах каждого конечного элемента, а глобальная нумерация устанавливает глобальные номера узловых значений в пределах всей конструкции. Для примера рассмотрим второй элемент стержня ца рис. 1.3. Локальными номерами его:, узлов является последовательность чисел 1, 2, а глобальными номерами — последовательность чисел 2, 3. В соответствии с двумя системами нумерации узловых значений существуют две системы адресации компонент векторов Г<'~ и матриц К~') конечных элементов — это локальная и глобальная адресации.

Локальный адрес устанавливает адрес компонент вектора Г(') или матрицы Кич в пределах тех же вектора или матрицы. Глобальный адрес устанавливает адрес компонент вектора Г<'> или матрицы Кип в пределах глобального вектора Г или глобальной матрицы К. Соответствие локального адреса глобальному и наоборот устанавливается матрицами кинематических связей конечных элементов а~'>. Каждая из матриц кинематических связей конечного элемента а'4 содержит в строке единственный ненулевой элемент, который равен единице.

Это неудивительно, поскольку при помощи а('~ каждому локальному узловому значению элемента ставится в соответствие единственное глобальное узловое значение вссй системы. Легко заметить, что номера столбцов ненулевых элементов матрицы аич совпадают с глобальными номерами узловых значений. Например, для матрицы аиз номерами столбцов ненулевых элементов являются два и три, что совпадает с глобальными номерами узловых значений (и, и и,).

В связи с этим можно отказаться от хранения матриц а('~ и вместо них ввести в рассмотрение двухмерный массив ХЧ списков глобальных номеров узловых значений. Для примера на рис. 1.3 этот массив будет иметь вид (1.44) Каждый столбец этого массива содержит глобальные номера узловых значений конечных элементов. В то же время номер строки г1Ч является локальным номером узлового значения конечного элемента. Введем обозначения для локальных и глобальных адресов элементов матрицы жесткости и компонент вектора нагрузки конечного элемента.