Курс лекций - Математическое моделирование технических объектов (1075784), страница 16
Текст из файла (страница 16)
Само интегрирование в (15.4) может быть выполнено черезL – координаты. Вспоминая выражение (13.39) и интегральную формулу (13.38), имеем для первого элемента:N1(1)N2(1)R(1) =0dVZY (e)(15.6)(1)N4V00Заменяя ФФ L – координатами:L1=N1 (1),L2=N2 (1), L3=N4 (1),имеем:(15.7)А(1) ZY (e)[1 1 0 1 0 0 ]T3Подставляя численные значения, имеем: R(1) = А(1) [78 78 0 78 0 0]TАналогичные вычисления для остальных конечных элементов приводят к следующим результатам:R(2) = А(2) [0 213 213 0 213 0]TR(3) = А(3) [0 165 0 165 165 0]TR(4) = А(4) [0 0 0 165 165 165]TОбъединяя матрицы по методу прямой жесткости и принимая А = А (1) = А(2) =А(3) = А(4) =получим выражение для столбца свободных членов в (15.2):R = А [78 456 213 408 543 165]TДля получения матрицы [C(1)] выразим ФФ первого конечного элемента черезL – координаты и воспользуемся выражением (15.7):N(1) = [N1(1) N2(1) 0 N4(1) 0 0] = [L1 L2 0 L3 0 0]Тогда произведение матриц в выражении (15.3) примет вид:R(1) =83[N(1)]T[N(1)] = [L1 L2L12 L1 L2L1 L2 L22=00L1 L3 L2 L300000 L3 00 L1 L30 L2 L3000 L3200000]T [L1 L2 0 L3 0 0]=0 00 00 00 00 00 0Запишем интегральную формулу для элемента а11 полученной матрицы:2A(1) 2!(2+2)!2SL1 dS ==A(1)12[2]Аналогично вычислим интеграл для элемента а12:SL1L2 dS =2A(1) 1! 1! 0!(1+1+2)!=A(1)12[1]Вычисляя аналогично остальные интегралы, получаем матрицу [C(1)]:A(1)[C ] =12(1)210100120100000000110200000000000000Таким образом, приходим к следующей системе уравнений:210100781120100456210000002133=12110200408400000054350000001656TРешение дает следующий результат: { } ZY = [71 437 724 354 671 476]T.Рассмотрим еще один пример, иллюстрирующий применение метода сопряженной аппроксимации.Пример 15.1.
Стальной конический стержень, показанный на рисунке 13.3-а,имеет длину 90 см и площадь поперечного сечения: у широкого закрепленного в стенеоснования S=12см2, у свободного узкого торца 6 см2. Стержень подвержен двум видам нагружения: а) в осевом направлении силой F=42 кН, приложенной с узкому торцу, и б) воздействию температуры t=40oC, действующей равномерно по всей длине.Требуется определить значение напряжений в основаниях стержня и в точках Аи В, делящих стержень на 3 равные части, если ТОС=200C, а температурный коэффициент расширения материала стержня равен: ТКС=710–6 [0C]–184Решение.Разбиваем стержень на три конечных элемента по L=30 см.
Задачу решаем в дваэтапа: на первом вычисляем перемещения в узлах, а на втором – вычисляем напряжения в элементах и методом сопряженной аппроксимации находим значения напряжений в заданных точках.I. Рассчитываем перемещения в узлах. С этой целью рассчитываем матрицы элементов.
Поскольку стержень ориентирован только в направлении оси абсцисс, имеет место только одна компонента тензора напряжений ХХ, связанная с компонентой тензора деформаций формулой:ХХ = Е(ХХ – ХХ0) = Е(ХХ – ТКТ)Здесь член ХХ0 учитывает начальную деформацию, связанную с тепловым расширением стержня. Для одномерного симплекс – элемента функция перемещенийимеет вид: u(е) = Ni(e)Ui+ Nj(e)Uj = [N]{U}.Выразим деформацию через функцию перемещений: ХХ=du/dx и, вычислимматрицу градиентов [В(1)], входящую в выражение для матрицы жесткости элемента:UiNi(1)Nj(1)1(1)(1)=[– 1 1] {U}ХХ =[В ] {U} =UjxxL[]{}Вспоминая формулу (12.23) и обозначая среднюю площадь элемента=(Ai+Aj)/2, записываем матрицу жесткости 1-го элемента, имеем:(1)K=Â(1) Е[1-1через Â=]-1 1LПодставляя численные значения и, вычисляя среднее значение площади сечения(1)1-го элемента Â = (12+10)/2 = 11 см2, имеем:(1)K =116,7106301-1-11= 1062,46-2,46-2,462,46Вектор нагрузки для первого конечного элемента связан только с вычислением интеграла, определяющего изменение объема от теплового расширения, который по определению равен: dVf(1) =V{B(1)}T EХХ0dV = [(EÂ(1)ТКТ)/L] [-1 1]T LdxПодставляя численные значения и учитывая, что интеграл равен L, имеем:85f(1) = 710-66,710611(40 – 20) [-1 1]T = [-10318 10318]TТаким образом, получаем систему уравнений для первого элемента:2,46-2,46U1-10318=106-2,462,46U210318[]Система уравнений для 2-го и 3-го элемента получается аналогично:106[2,01-2,01-2,012,01]U2106[1,56-1,56-1,561,56]U3-8442=U38442-6566=U46566К столбцу свободных членов последней системы необходимо добавить интеграл f(3),суммирующий нагрузку по поверхности тонкого (свободного) торца стержня.
Чтобы получить выражение для вычисления f(3), запишем дифференциальное уравнение, описывающееперемещения стержня при его осевом нагружении: du/dx – F/AE = 0 или, что то же:E(du/dx) – F/A = 0(15.8)Сравнивая выражения (15.8) и (13.28), можно по аналогии с (13.3) записать:f(3) =(F/A) [N(3)]T dV =F/A(3)V(3)[0 1] T dS(15.9)S(3)Так как нагрузка приложена к правому узлу 3-го конечного элемента, где х=L, то выражение для N(3)принимает вид: [N(3)]=[(1-x/L) (x/L)]=[0 1].
Далее, площадьS(3)= A(3), поэтому окончательно можно записать:f(3) = [ 0 42000]TОбъединяя матрицы по методу прямой жесткости, приходим к следующей системе уравнений:1062,46-2,460-2,464,47-2,01000-2,013,57-1,56U1U2U3-10318=1876187600-1,56 1,56U448566Преобразование (так как U1=0) и решение системы дает следующий результат:{ U }T =[0 2,07 4,5 7,53]T{ U }Tтеоретическое =[0 2,10 4,6 7,80]TДанные приведены в мм.
Теоретические значения получены интегрированием деформации подлине.II. Рассчитываем напряжения в узлах по методу сопряженной аппроксимации.1. По найденным узловым перемещениям находим деформацию элементов:ХХ(1) = (– U1 + U2)/L = (-0,0000 + 0,0207)/30 = 0,6910-386ХХ(2) = (– U2 + U3)/L = (-0,0207 + 0,0450)/30 = 0,8110-3ХХ(3) = (– U3 + U4)/L = (-0,0450 + 0,0753)/30 = 1,0110-32. Напряжение в е-м элементе равно: ХХ(е) = Е(ХХ(е) – ХХ0) = Е(ХХ – ТКТ).Подробно для 1-го элемента имеем:ХХ(1) = 6,7106 {0,6910-3 – 710-6(40-20)} = 3685 [H/см2]Аналогично вычисляем напряжения в остальных конечных элементах:ХХ(2) = 4480 [H/см2];ХХ(3) = 5820 [H/см2]3. Составляем уравнения сопряженной аппроксимации: [C] {} = {R}, для чего предварительно вычислим произведение матриц [N (e)]T [N (e)]:L1L12L1L2[L1 L2] =L2L1L2L22Учитывая, что L L12 dx = 2!L/(2+1)! = L/3 и L L1L2 dx = 1!1!L/(1+1+1)! = L/6, вычисляем левую часть искомой системы уравнений:L2 11[C] {} =[]61 22Составление матрицы {R} требует вычисления интеграла вида:LL1 dx = 1!L/(1+1)! = L/2Приходим к системе уравнений для 1-го конечного элемента:13[2112]12=36853685Аналогичные вычисления для 2-го и 3-го элементов дают системы:13[2112]23=4480448058203121=31258204Объединение полученных матриц (по методу прямой жесткости) и решение системыдает следующий результат: {} T = [3558 3935 5222 6132]T (Н/см2).
Теоретическое значениевектора получим делением величины приложенной нагрузки на площадь поперечного сечения в соответствующем узле. Данные сведем в таблицу:[N п/п1234Теоретическое3500420052507000]Рассчетное3558393552226132Внутри элементов36854489582916.1 Элементы высокого порядкаОсновные причины, по которым прибегают к использованию элементов высокого порядка – комплекс- и мультиплекс- элементов:87-более высокая точность решения при таком же количестве элементов (или достижение той же точности при меньшем числе элементов);невозможность аппроксимации с помощью симплекс – элементов градиентов искомых величин кусочно-линейными функциями;при использовании элементов высокого порядка отпадает необходимость в применении теории сопряженной аппроксимации.Определение квадратичного элемента. Рассмотрим порядок вычисления ФФ для одномерного квадратичного комплекс – элемента и методику его использования для решения конкретной задачи.Одномерный квадратичный комплекс – элемент (второго порядка) представлен на рисунке 15.1.
Его аппроксимирующий полином имеет вид: = 1 +2 x +3 x2(16.1)Число узлов квадратичного элемента равно трем (i, j, k). Коэффициенты 1, 2 и 3 определяются из условий: = Ф при x = X ( = i, j или k). Если поместить узел i в начало координат, указанные условия примут вид: = Фi при x = 0; = Фj при x = L/2 ; = Фk при x = L;Эти узловые условия приводят к системе уравнений, решив которую получим:1= Фi; 2 =(4Фi-3Фj-3Фk)/L; 3 =2(Фi-2Фj+Фk)/L2.Рис. 16.1Подставляя выраженные через узловые значения коэффициенты в (16.1) и перегруппировывая члены, получим интерполяционный полином для квадратичного комплекс – элемента в матричном виде: = [Ni Nj Nk] { Фi Фj Фk }T(16.2)здесь: Ni = (1-2x/L) (1-x/L); Nj = (4x/L) (1-x/L); Nk = (x/L) (1-2x/L).Существуют и широко применяются на практике формальные способы вычисленияфункций формы, использующие их известные свойства.Применение квадратичного элемента.
Элементы высокого порядка применяются также, как симплекс – элементы, поскольку выбор интерполяционного полинома не связан с исходными дифференциальными уравнениями. В качестве примера рассмотрим одномернуюзадачу переноса тепла. Ее решение было рассмотрено в разделе 12 с использованием симплекс – элементов. Задача касалась определения распределения температуры по длинестержня, подверженного конвективному теплообмену. Исходные уравнения для произвольного элемента, выведенные ранее, имеют вид:88[K(e)] {T} = {F(e)},где:[K(e)] = q [N(e)]T dS + [B(e)]T [B(e)]dV +V[F(e)] = -s1S2h [N(e)]T [N(e)] dSS2h TOC[N(e)]T dSВсе интегралы должны быть вычислены заново, если мы используем квадратичный (вместолинейного) элемент.