Учебное пособие по курсу лекций (1164067), страница 55
Текст из файла (страница 55)
кусочно-линейную аппроксимацию.Такой элемент называют линейным.При кусочно-линейной аппроксимации функции, распределенной водномерной области, коэффициенты полинома оказываются зависящими отзначений неизвестной функции в узлах fi, fi+1 и геометрических размеровэлемента xi, xi+1 (i=1…M, M - число узлов):f − fϕ j = a0, j + a1, j x ; ϕ j = i +1 i ( x − xi ) + f ixi +1 − xi357fffiϕjϕjfi+1j-1j+1jxixi+1xxjИз рисунка видно, что размеры элементов необходимо согласовывать синтенсивностью изменения функции. Там, где функция имеет наибольшийградиент, для более точной аппроксимации размер элементов должен бытьменьше.Более точной аппроксимации можно добиться, используя полиномвторого порядка.
Для однозначного определения коэффициентов полиномавторого порядка необходимо знать значение функции в трех узлах. Такойэлемент называют квадратичным. Один квадратичный элемент обычноаппроксимирует функцию точнее, чем два линейных. В этом случае элементимеет 3 узла, два из которых граничат с другими элементами, а один –внутренний.Аналогичные рассуждения можно провести для двумерных итрехмерных областей. Таким образом, коэффициенты полинома,аппроксимирующего распределение неизвестной функции внутри элемента,однозначно определяются значением функции в узлах и геометрическимиразмерами элементов.
Таким образом, если нам известны геометрическиепараметры разбиения области на элементы и значения функции в узлахэлементов, мы можем определить значение функции в произвольной точкеобласти.Вопрос состоит в том, каким образом по известным граничнымусловиям определить значения неизвестной функции в узлах. В методеконечных элементов существуют два подхода к решению этой задачи.Рассмотрим эти методы на примере задач теории упругости.Первый подход – т.н. «прямой метод». В этом методе разрешающаясистема уравнений строится на основе применения условий равновесия вузлах сетки КЭ.Второй подход – «вариационный». Здесь разрешающая системауравнений строится на основе минимизации функционала дополнительнойработы (энергии).
Согласно вариационному принципу Лагранжа, «среди всехкинематически возможных перемещений точное решение соответствуетабсолютному минимуму функционала дополнительной энергии».358Π = Λ −W ,где П – дополнительная энергия, Λ - энергия деформации, W - работавнешних сил.Для задач теории упругости определяющие уравнения, полученные наоснове обеих подходов полностью совпадают. Прямой метод в задачахтеории упругости для инженера является более наглядным. Но вариационныйметод является более общим. Можно сказать, что развитие МКЭ получилосущественный толчок, когда в 60-х годах 20 века было показано, чтоопределяющие соотношения МКЭ в большинстве задач математическойфизики могут быть получены на основе минимизации функционалов,связанных с физической сущностью задачи.Рассмотрим простейший пример решения одномерных задач теорииупругости методом конечных элементов.
Конечно, решать такие задачиметодом конечных элементов не имеет практического смысла, однако такойпростейший пример позволяет проследить всю процедуру метода конечныхэлементов вручную.Необходимо составить определить перемещения точек стержняпеременного сечения под действием сосредоточенной силыБудем считать, что металл не переходит в пластическое состояние, т.е.остается упругим. В упругом состоянии металл работает как пружина. Силасжатия пружины пропорциональна ее жесткости и деформации P = kδ .Если сечение постоянно, то жесткость пружины также постоянна.ES, здесь E – модуль упругости, S – площадь, H – длина стержня.k=H359Если бы стержень был постоянного сечения, то мы бы легко решилизадачу.Q QHu= =k ESНо опора имеет переменное сечение, поэтому жесткость ее зависит отфункции, описывающей внешнюю поверхность.
Конечно, такая задача можетбыть решена довольно просто аналитически, но мы попытаемся рассмотретьее в конечноэлементной постановке.Аппроксимируем функцию u(x) совокупностью кусочно-линейныхфункций ϕi(x)Lu ( x) = ∑ ϕi ( x)i =1ϕi = a0,i + a1,i xКоэффициентыкусочно-линейнойфункцииопределяютсякоординатами узлов xj и значениями перемещений в узлах Uj⎛U j +1 − U jU j +1 − U j ⎞ U j +1 − U jϕi =x − x j + U j = ⎜U j −xj ⎟+x⎜⎟x j +1 − x jxxxx−−j +1jj +1j⎝⎠Таким образом выполнены следующие действия по составлениюматематической модели:стержень разделен на L КЭ, имеющих по 2 узла каждый;в каждом КЭ неизвестное перемещение аппроксимировано полиномом1-й степени;коэффициенты полиномов зависят от координат и перемещений узлов.В общем виде можно сказать, что выполнена замена неизвестнойфункции перемещения u(x) на совокупность неизвестных значений этойфункции – перемещений узлов сетки конечных элементов Uj⎧ U1 ⎫⎪U ⎪⎪⎪u ( x) ⇒ U j = ⎨ 2 ⎬ = {U } , где {U} - вектор узловых неизвестных⎪… ⎪⎪⎩U N ⎪⎭Для определения узловых неизвестных, используя физическиеуравнения решаемой задачи, составляют основное разрешающее уравнение(уравнение жесткости) в виде СЛАУ{R} = [ K ]{U }[K] - матрица жесткости; {R} - вектор столбец внешних воздействий.Процесс составления основного разрешающего уравнения системыможно подразделить на два этапа:создание уравнения для конечного элементасборка уравнения для системы из уравнений элементов(ансамблирование)В учебных целях воспользуемся прямым методом.()360В прямом методе :уравнения для элемента строят путем решения исходной системыдифференциальных уравнений, используя принятую аппроксимациюнеизвестной функции.сборка уравнений элементов в единую систему уравнений для областиосуществляют путем применения условий равновесия в узлах сеткиКЭ.Выполненнаянамиранеекусочно-линейнаяаппроксимацияперемещения эквивалентна разделению стержня переменного поперечногосечения на совокупность нескольких линейно-упругих конечных элементов,каждый из которых характеризуется длиной L и площадью поперечногосечения S.QL1L2LjHQL1,S1Lj,SjL2,S2Линейный упругий КЭ – это стержень постоянного сечения,подвергающийся упругим деформациям вдоль оси.
Узлами таких конечныхэлементов являются концы отрезков, которыми можно схематичноотобразить этот элемент.L, Su1P1ku2P2Конечный элемент имеет 2 степени свободы, соответствующиеперемещениям узлов вдоль оси стержня. Выполним построение системыуравнений для одномерного линейно-упругого элемента прямым методом.Дифференциальное уравнение среды, описывающее поведениеэлемента можно получить из закона Гука:Pσ = = EεS∂u, откудаДеформация ε =∂x∂u P=∂x SE361Аппроксимируем перемещение внутри стержня полиномом первойстепени:u −uu −uu ( x) = 2 1 ( x − x1 ) + u1 = 2 1 ( x − x1 ) + u1x2 − x1LЗдесь u1, u2 - перемещения узлов, x1, x2 - координаты узлов.Тогда∂u u2 − u1=∂xLОкончательно внутренняя сила в элементеES(u2 − u1 ) = k (u2 − u1 ) , где k = ESP=LLВнешние силы, воздействующие на узлы:P1 = k (u1 − u 2 ),P2 = k (u 2 − u1 )Или в матричном виде:⎧ P1 ⎫ ⎡ k11 k12 ⎤ ⎧ u1 ⎫ ⎡ k − k ⎤ ⎧ u1 ⎫⎨ ⎬=⎢⎥⎨ ⎬ = ⎢⎥⎨ ⎬⎩ P2 ⎭ ⎣k 21 k 22 ⎦ ⎩u 2 ⎭ ⎣− k k ⎦ ⎩u 2 ⎭Общая запись называется уравнением жесткости конечного элемента.{P( )} = ⎡⎣⎢k ( ) ⎤⎦⎥ {u( )} ,здесь{u( )} eeeeматрица-векторузловыхнеизвестных,{ }eразмерностью N (e) × 1, где N (e) – число степеней свободы элемента; P( ) -матрица-вектор узловых внешних воздействий, размерностью N (e) × 1;⎡ k ( e ) ⎤ - матрица жесткости элемента, размерностью N (e) × N (e) , зависящая⎢⎣⎥⎦от свойств и геометрических размеров элемента.Это уравнение является математической моделью поведения конечногоэлемента под действием внешних нагрузок, приложенных к узлам.Составляющие этого уравнения – вектор узловых сил КЭ, матрица жесткостиКЭ, узловые перемещения КЭ.
Поскольку метод конечных элементовпервоначально был применен к задачам теории упругости, то терминологиюраспространили и на все остальные задачи математической физики. В любыхзадачах – тепловых, задачах пластичности, гидравлики и т.д. матрица,связывающая вектор узловых неизвестных с вектором узловых внешнихвоздействий называется матрицей жесткости, хотя в физической основезадачи никакой жесткости может и не быть.Таким образом, матрица жесткости для упругого стержневогоконечного элемента имеет вид:362ES ⎤⎡ ES−⎢L ⎥ = ES ⎡ 1 −1⎤⎡ k (e) ⎤ = ⎢ L⎥⎣⎦ ⎢ ESES ⎥ L ⎢⎣ −1 1⎥⎦−⎢⎣ LL ⎥⎦Для простоты разделим опору на два конечных элемента.
Общностирешения задачи это не нарушит, но позволит легко получить определяющиесоотношения.Тогда система будет состоять из 2-х конечных элементов и иметь 3глобальных степени свободы. Каждой глобальной степени свободы U mсоответствует узловая внешняя сила Rm .Вектора глобальных степеней свободы и внешних узловых воздействийв данном случае имеют следующий вид:⎧U1 ⎫⎧ R1 ⎫ ⎧ Q ⎫⎪ ⎪ ⎪ ⎪⎪ ⎪{U } = ⎨U 2 ⎬ ; {R} = ⎨ R2 ⎬ = ⎨ 0 ⎬⎪ ⎪ ⎪ ⎪⎪U ⎪⎩ R3 ⎭ ⎩−Q ⎭⎩ 3⎭Каждый из двух конечных элементов, на которые мы разбили систему,- линейный упругий элемент, имеет по 2 локальных степени свободы:u1(i ) , u2(i ) . Здесь (i) – номер конечного элемента.
Соответственно узлахкаждого элемента действуют узловые силы P1(i ) , P2(i ) .Свойства каждого элемента описываются уравнениями жесткостиэлемента в локальной системе координат (ЛСК) – системе координат,связанной с элементом:363⎧ P1(1) ⎫ ⎡ k1 − k1 ⎤ ⎧u1(1) ⎫ ⎧ P1( 2) ⎫ ⎡ k2 − k2 ⎤ ⎧u1( 2) ⎫⎨ (1) ⎬ = ⎢⎥ ⎨ (1) ⎬ , ⎨ ( 2) ⎬ = ⎢⎥ ⎨ ( 2) ⎬⎩ P2 ⎭ ⎣− k1 k1 ⎦ ⎩u2 ⎭ ⎩ P2 ⎭ ⎣− k2 k2 ⎦ ⎩u2 ⎭Cоставим систему уравнений для всего деформирующегося тела наосновании уравнений для каждого конечного элемента.
Этот этап носитназвание ансамблирование. Используем для этих целей уравненияравновесия: сумма внутренних сил в каждом узле равна внешней силе,действующей на узел. Иными словами используем «прямой метод»составления разрешающей системы уравнений.Условие равновесия в узлах КЭ в ЛСК:()R1 = P1(1) = k1 u1(1) − u2(1) = Q() ()R2 = P2(1) + P1(2) = − k1 u1(1) − u2(1) + k2 u1(2) − u2(2) = 0()R3 = P2(2) = − k2 u1(2) − u2(2) = −QДля составления системы уравнений в глобальной системе координатвоспользуемся соответствием локальных и глобальных координат, котороеследует из рисунка:U1 = u1(1) , U 2 = u1(2) = u2(1) , U 3 = u2(2) .Окончательно условия системы равновесия в глобальной системекоординат имеют вид:k1 (U1 − U 2 ) = Q−k1U1 + ( k1 + k2 )U 2 − k2U 3 = 0k2 (U 3 − U 2 ) = −QПолученная система может быть записана в матричном виде:{R} = [ K ]{U } ,где:⎧ R1 ⎫ ⎧ Q ⎫{R} = ⎪⎨R2 ⎪⎬ = ⎪⎨ 0 ⎪⎬ - вектор глобальных узловых нагрузок (известен)⎪ R ⎪ ⎪− Q ⎪⎭⎩ 3⎭ ⎩⎧U1 ⎫{U } = ⎪⎨U 2 ⎪⎬ - вектор глобальных узловых перемещений (неизвестен)⎪U ⎪⎩ 3⎭⎡ (1)0 ⎤ ⎢ k11(1)− k2 ⎥ = ⎢ k21⎥ ⎢k2 ⎥⎦ ⎢ 0⎣жесткости (известна).−k1⎡ k1[ K ] = ⎢⎢ −k1 k1 + k2⎢⎣ 0− k2(1)k12(1)(2)k22+ k11(2)k210 ⎤⎥(2) ⎥k12⎥(2) ⎥k22⎦- глобальная матрица364Следует обратить внимание на структуру глобальной матрицыжесткости.