Диссертация (1025882), страница 19
Текст из файла (страница 19)
Тем самым осуществляется переход от непрерывной системы к системе сконечных числом степеней свободы. Аппроксимация проводится при помощифункций формы . Для восьмиузлового КЭ она имеет вид ∼=8∑︁(4.65) .=1Поскольку геометрия тела может быть произвольной, функции формы КЭудобно задавать в локальной системе координат. На Рисунке 4.7 представленэлемент в локальной системе 0 в виде кубика со сторонами, равными 2, ипронумерованными вершинам (узлами КЭ). Для этого элемента функции формы представимы в виде =1(1 + ) (1 + ) (1 + ) ,8(4.66)где ( , , ) – локальные координаты -го узла КЭ.
Для изопараметрического23416785Рисунок 4.7. Восьмиузловой конечный элементКЭ взаимно-однозначное соответствие между локальными (, , ) и глобаль-130ными координатами КЭ устанавливается соотношением =8∑︁ () ,(4.67)=1где - глобальные координаты узлов КЭ. Выражение (4.67) является точным заданием координат произвольной точки КЭ через координаты узлов . Используя представления (4.65), (4.67) можно вычислить производные отперемещений по глобальным координатам недеформированной системы, необходимые для определения деформаций тела и скоростей движения его точек.Для этого необходимо сначала определить матрицу Якоби :8 () ∑︁= ⊗.=(4.68)=1Градиент вектора перемещения вычислим путем дифференцирования перемещения как сложной функции8∑︁ ()grad ( ) ==·= ⊗· −1 .
(4.69)=1Воспользовавшись соотношениями Коши, получим выражение для деформации=)︀1 (︀grad ( ) + grad ( ) =2)︂8 (︂∑︁ () ()−1−1= ⊗· +· ⊗ .(4.70)=1Равенства (4.65) – (4.70) вводятся в сформулированный ранее принцип виртуальных работ (4.36). Вычисляемые при этом интегралы по объему Ω и поповерхности Ω тела приближенно заменяются интегрированием по объему ℬи по поверхности ℬ конечно-элементной системы. Эти интегралы представляют собой суммы интегралов, вычисляемых по каждому из конечных элементов131в отдельности∫︁(...) ∼=∫︁(...) =Ωℬ∫︁∫︁Ωc(...) ∼=∑︁ ∫︁(...) =ℬ(4.71)∑︁ ∫︁ℬc(...) ,(...) .c(ℬ )В полученных суммах осуществляется переход от интегралов, записанных в глобальной системе координат, к интегралам по объему и поверхности конечныхэлементов в локальных координатах, которые вычисляются методом квадратурГаусса-Лежандра.
После чего собираются все слагаемые при одинаковых вариациях узловых перемещений. Такую операцию ассемблирования будем обозна⋃︀чать символом . В результате принцип виртуальных работ (4.36) принимаетвид⎞⎛∫︁⋃︁ ⎜∫︁⎟⎟ = 0.⎜ ·· −f·⎠⎝cℬ(ℬ )(4.72)Поскольку вариации узловых перемещений независимы друг от друга, то множители перед ними должны обращаться в ноль, чтобы удовлетворить уравнению (4.72). Выполнение этого условия приводит к системе нелинейных алгебраических уравнений относительно узловых перемещений системы. Посколькуформирование и решение этой системы производится с использованием ЭВМ,ее принципиальный вид должен быть получен в векторно-матричной форме.Введем матрицу функций форм[︁]︁[ ] = [ ]1 .
. . [ ]8 ,(4.73)состоящую из восьми блоков, отвечающих каждому из узлов КЭ⎡⎤ 0 0⎢⎥[ ] = ⎣ 0 0 ⎦ .0 0 (4.74)132Координаты узлов КЭ и узловые перемещения представим в виде векторов{︁}︁{ } = { }1 . . . { }8 ,{︁}︁{} = {}1 . . . {}8 ,где координаты и перемещения -го узла выделены в векторы{︁}︁{ } = 1 2 3 ,{︁}︁{} = 1 2 3 .Используя эти обозначения, перепишем соотношение (4.65) и (4.67):{ } ∼= [ ] {} ,{ } ∼= [ ] { } .(4.75)(4.76)Применяя соотношения Коши к функциямперемещений (4.75),}︁выразим ком{︁поненты тензора деформации {} = 11 22 33 12 23 13 :{} = [] {} ,(4.77)где [] – матрица форм деформаций, имеющая вид⎡⎤/100⎢⎥⎢ 0⎥/02⎢⎥⎢ 00/3 ⎥⎢⎥[] = ⎢⎥ [ ].⎢/2 /1⎥0⎢⎥⎢ 0/3 /2 ⎥⎣⎦/30/1(4.78)По известному полю деформаций тела из определяющих соотношений (4.5)–(4.8) вычисляются компоненты тензора напряжений.
Выбранное представление поля перемещения в виде функции гладкости 0 приводит к разрыву полей деформаций и напряжений между соседними конечными элементами. Длявычисления вязких деформаций в виде (4.54) на текущей итерации вводитсяотдельная аппроксимация поля деформаций, полученного на предшествующей133итерации. В качестве иллюстрации используемого приближения на Рисунке 4.8схематично представлена область тела, разбитая на КЭ. Рассмотрим произвольный ряд элементов, выстроенных в окружном направлении, по которомупроизводится численное интегрирование истории деформирования. Между одноименными гауссовыми точками (т.е.
точками с одинаковыми локальными координатами) двух смежных КЭ, в которых деформации равны {}−1 и {} ,вводится линейная аппроксимация(︁)︁{} = {}−1 + {} − {}−1 ,где ∈ [0 , 1] .(4.79)В методе Рунге-Кутта шаг по времени выбирается равным углу Δ междууказанными точками соседних КЭ. В результате интегрирования определяетсявектор вязких деформаций { }. На текущей итерации он является известнойвеличиной.Ряд КЭOГауссова точкаРисунок 4.8. К вычислению напряжений в телеРезультирующий вектор компонент напряжений {} представим в виде{} = [] {} − 2 { } ,(4.80)где [] – матрица упругих постоянных, определяемая объемным модулем134упругости и суммарным касательным модулем Σ = + ,⎡+⎢⎢ −⎢⎢ −⎢[] = ⎢⎢⎢⎢⎣43 (Σ )23 (Σ )23 (Σ )−+−23 (Σ )43 (Σ )23 (Σ )−−+⎤23 (Σ )23 (Σ )43 (Σ )ΣΣΣ⎥⎥⎥⎥⎥⎥.⎥⎥⎥⎦(4.81)Конечноэлементный аналог уравнения (4.51) с учетом выражений (4.55), (4.57),(4.64) запишем в виде(︀[︀ ]︀ )︀{︀ }︀[ ] + [ ] + st ] {} = {f } + {f } + fst ,(4.82)[︀ ]︀где [ ] , [ ] , st – матрицы жесткости всего конечно-элементного ансамбля, обусловленные работой сил упругости, контактных сил давлений, сил трения сцепления; в левой части записаны соответствующие им векторы узловыхсил; {} – искомый вектор перемещения узлов конструкции.
Ниже приведенокончательный вид используемых в расчете матриц и векторов.Из выражения (4.55) получим следующее представление для матрицыжесткости упругости [ ] и вектора нагрузки {f }, порождаемого вязкостьютела:⋃︁ ∫︁[ ] =[] [] [] ,(4.83)ℬ{f } = 2⋃︁ ∫︁[] { } .(4.84)ℬОбозначив за { } вектор нормали к поверхности барабана, получим матрицужесткости и вектор сил нормального контакта[ ] =⋃︁ ∫︁ [ ] { } { } [ ] ,c(ℬ )(4.85)135{f } =⋃︁ ∫︁ g− [ ] { } .(4.86)c(ℬ )Для представления матрицы жесткости и вектора нагрузки, обусловленных силами трения сцепления в контакте, введем обозначения для вектор-столбцовединичных ортов {2 } = {0, 1, 0} и {3 } = {0, 0, 1} .
В результате будем иметь[︀ st ]︀ ⋃︁ =∫︁ [ ](︁{2 } {2 } + {3 } {3 })︁[ ] ,(4.87))︁(4.88)st{︀ st }︀f =(ℬ )⋃︁ ∫︁ 2 [ ] {2 } +st+(ℬ )⋃︁ ∫︁ [ ](︁{2 } {2 } + {3 } {3 }{(0 )} .st(ℬ )Алгоритм и основные процедуры решения вязкоупругой задачи каченияприведены в приложении П.1.4.7Исследование поля температур при стационарном качениимассивной шиныЭкспериментально [10] и теоретически [12] установлено, что при стационарном режиме качения температурное поле внутри массивной шины не меняется со временем. Во всех ее осевых сечениях распределения температур одинаковые. Эти утверждения позволяют сформулировать температурную задачукак осесимметричную.На Рисунке 4.9 изображено осевое сечение массивной шины.
Посколькусечение симметрично относительно центральной оси, температурное поле определяется только для половины сечения, которая показана сплошными линиями.Запишем уравнение стационарной теплопроводностиdiv ( grad ) + w = 0 ,(4.89)где – температура, являющаяся функцией координат точек области ; –коэффициент теплопроводности; w – мощность внутренних источников теплоты, которая положительна в случае выделения теплоты. При записи уравне-136ния (4.89) использовано определяющее соотношение, известное как закон Фурье, связывающее тепловой поток с градиентом температуры(4.90) = − grad .металлРисунок 4.9.
Осевое сечение резинового массиваМощность внутренних источников теплоты w определяется напряженнодеформированным состоянием резинового массива по формулеw=∫︀2 · ˙ 02,(4.91)где точкой ( ˙ ) показано дифференцирование по угловой координате . Причисленном моделировании задачи качения напряжения и деформации точнеевсего определяются в гауссовых точках конечных элементов. Для этих точекосуществлялся подсчет интеграла, стоящего в числителе выражения (4.91).Сформулируем краевые условия задачи.