Диссертация (1151675), страница 23
Текст из файла (страница 23)
Компоненты тензора деформаций вычисляют относительно неизменного глобального базиса. Наибольшее распространение находят простейшие восьмиузловые, шестиузловые и четырехузловые изопараметрические КЭ [139, 260], требующие для формирования матрицы жёсткостиансамбля элементов относительно небольших затрат машинного времени, а такжепозволяющие автоматизировать процесс ансамблирования (построения конечноэлементной сетки).
Необходимая точность достигается за счёт изменения шага168сетки. Практический интерес представляет распространение существующей стандартной схемы МКЭ на случай, когда напряжённо-деформированное состояниеопределяется в локальных осях элемента, совпадающих с главными направлениями анизотропии материала.Ниже приведён матричный алгоритм построения матрицы жёсткости ивектора обобщённых сил в осях ортотропии для объёмных изопараметрическихКЭ [28, 29].Матрицы жёсткости объёмных КЭ. Функции формы.– 8-ми узловой КЭ⎡ − 1 1 − 1 1 − 1 1 − 1 1⎤[ p ] = ⎢⎢ − 1 − 1 1 1 − 1 − 1 1 1⎥⎥( 3×8 )⎢⎣ − 1 − 1 − 1 − 1 111 1⎥⎦– 6-ти узловой КЭ:1ϕ 1 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 11 ) (1 + x 2 p 21 ) (1 + x 3 p 31 ) + (1 + x p ) (1 + x p ) (1 + x p ) ]1 13223333818;18;ϕ 2 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 1 2 ) (1 + x 2 p 2 2 ) (1 + x 3 p 3 2 ) ]ϕ 3 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 1 4 ) (1 + x 2 p 2 4 ) (1 + x 3 p 3 4 )];1ϕ 4 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 15 ) (1 + x 2 p 2 5 ) (1 + x 3 p 3 5 ) + (1 + x p ) (1 + x p ) (1 + x p ) ]1 17227337818;18.ϕ 5 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 16 ) (1 + x 2 p 2 6 ) (1 + x 3 p 3 6 )]ϕ 6 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 18 ) (1 + x 2 p 2 8 ) (1 + x 3 p 3 8 )]– 5-ти узловой КЭ:18ϕ 1 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 11 ) (1 + x 2 p 21 ) (1 + x 3 p 31 )];18;18;ϕ 2 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 1 2 ) (1 + x 2 p 2 2 ) (1 + x 3 p 3 2 )]ϕ 3 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 13 ) (1 + x 2 p 2 3 ) (1 + x 3 p 3 3 )];16918ϕ 4 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 1 4 ) (1 + x 2 p 2 4 ) (1 + x 3 p 3 4 )];1ϕ 5 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 15 ) (1 + x 2 p 2 5 ) (1 + x 3 p 3 5 ) + (1 + x p ) (1 + x p ) (1 + x p ) +1 162263368+ (1 + x 1 p 17 ) (1 + x 2 p 2 7 ) (1 + x 3 p 3 7 ) + (1 + x 1 p 18 ) (1 + x 2 p 2 8 ) (1 + x 3 p 38 ) ].– 4-х узловой КЭ:1ϕ 1 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 11 ) (1 + x 2 p 21 ) (1 + x 3 p 31 ) + (1 + x p ) (1 + x p ) (1 + x p ) ]1 13223333818;18;ϕ 2 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 1 2 ) (1 + x 2 p 2 2 ) (1 + x 3 p 3 2 )]ϕ 3 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 1 4 ) (1 + x 2 p 2 4 ) (1 + x 3 p 3 4 )];1ϕ 4 ( x 1 , x 2 , x 3 ) = [ (1 + x 1 p 15 ) (1 + x 2 p 2 5 ) (1 + x 3 p 3 5 ) + (1 + x p ) (1 + x p ) (1 + x p ) +1 162263368+ (1 + x 1 p 17 ) (1 + x 2 p 2 7 ) (1 + x 3 p 3 7 ) + (1 + x 1 p 18 ) (1 + x 2 p 2 8 ) (1 + x 3 p 38 ) ]Ковариантные составляющие тензора деформацийствующих” осях{x1 , x 2 , x 3}εij.в местных „сопут-в случае малых перемещений определяют по формуле[139]:12ε i j = ( z m , j ∇ i u m + z m , i ∇ j u m ) , i , j = 1, 2 , 3 ,(4.1)где {z 1, z 2, z 3} - в общем случае криволинейные глобальные координаты;z m , j= ∂ z m / ∂ xj- компоненты вектора перемещений в осях;umz- компоненты тензора преобразования координат;mковариантная производная от ковариантных компонент вектора переме∇iu m = u m , i − z m , i Γ nm k u n ;щений:u m, i = ∂ u m / ∂ x i;Γ nm k(4.2)- символы Кристоффеля второго рода.
В выражениях(4.1), (4.2) принято правило суммирования по повторяющемуся („немому”) индексу. Индекс m последовательно принимает значения 1, 2.Рассмотрим случай, когда глобальные оси являются прямолинейными де-170картовыми, а местные криволинейными. Как известно [260], для декартовых координатΓ nm k = 0и ковариантная производная совпадают с обычной производнойкомпонент вектора по координате. Тогда выражение (4.1) принимает вид:ε i j=1(z m , jum , i + z m , iu m , j ) .2(4.3)Согласно „моментной схеме”, деформации элемента представляем в видеотрезка ряда Тейлора в окрестности начала местной системы координат:8-ми узловой КЭε11 = ε~11 + ε~11,2 x 2 + ε~11,3 x 3 + ε~11,23 x 2 x 3 ε 22 = ε~22 + ε~22,1 x 1+ ε~22,3 x 3 + ε~22,13 x 1 x 3,,ε 33 = ε~33 + ε~33,1 x 1+ ε~33,2 x 2 + ε~33,12 x 1 x 2 ε 12 = ε~12 + ε~12,3 x 3,ε13 = ε~13 + ε~13,2 x 2,ε 23 = ε~23 + ε~23,1 x 1,;– 6-ти узловой КЭε11 = ε~11 + ε~11, 3 x3 ε 22 = ε~22 + ε~22 , 3 x3 ε 33 = ε~33 + ε~33,1 x1 + ε~33, 2 x2,,,ε12 = ε~12 + ε~12, 3 x3 ε13 = ε~13 ε 23 = ε~23,,;– 5-ти узловой КЭε 11 = ε~11 + ε~11, 23 x 2 x 3 ε 22 = ε~22 + ε~22,13 x 1 x 3 ε 33 = ε~33 + ε~33,12 x 1 x 2,,,ε12 = ε~12, ε13 = ε~13 , ε23 = ε~23 ;– 4-х узловой КЭε11 = ε~11, ε22 = ε~22 , ε33 = ε~33 ,ε12 = ε~12, ε13 = ε~13 , ε 23 = ε~23 .Здесьε~ i j , ε~ i j ,α , ε~i j ,α β- „моменты деформаций” в точкеИндексы α , β = 1 , 2 , 3 .Введём векторы – столбцы деформаций{ε } = {ε11 ε 22 ε 332ε122ε132ε 23 }Τи узловых перемещений{w} = {{um(1) }{u }( 2)m{ }} , {u } = {u(n )K um rΤ(k )m(k )1}u 2( k ) u3( k ) ,x0x 1 =0, x 2 =0, x 3 =0.171nr– число узлов КЭ (nr =4, 5, 6, 8).Связь между {ε } и {w} с учётом выражения (4.3) представим в матричной{ε} = [D] {w},формегде(1)( 2)( 3)[ D] = [ [ D] 1 [ D] 2 K[ D] n r ] ; [ D] k = [{D } k {D } k {D } k ] ,( 6×3 n r )( 6×3)k =1, 2, …8;8-ми узловой КЭz m ,1 + ( ~z m ,12 + ~z m ,1 p 2 k ) x 2 + ( ~z m ,13 + ~z m ,1 p 3 k ) x 3 + ⎫⎧ p1k [ ~⎪⎪ ~~~~⎪⎪ ( z m ,123 + z m ,12 p 3 k + z m ,13 p 2 k + z m ,1 p 2 k p 3 k ) x 2 x 3 ]zm,2 + (~z m ,12 + ~z m , 2 p1k ) x1 + ( ~z m , 23 + ~z m , 2 p3 k ) x3 + ⎪⎪ p2k [~⎪⎪ ~~~~⎪⎪ ( z m ,123 + z m ,12 p 3 k + z m , 23 p1k + z m , 2 p1k p 3 k ) x1 x3 ]z m ,3 + ( ~z m ,13 + ~z m , 3 p1k ) x1 + ( ~z m , 23 + ~z m ,3 p 2 k ) x 2 + ⎪⎪ p3k [ ~⎪⎪z m ,123 + ~z m ,13 p 2 k + ~z m , 23 p1k + ~z m , 3 p1k p 2 k ) x1 x 2 ]1 ⎪ (~⎪(m){D } k = ⎨ ~⎬ ;~~~~8 ⎪ z m ,1 p 2 k + z m , 2 p1k + ( z m ,13 p 2 k + z m ,1 p 2 k p 3 k + z m , 23 p1k + ⎪⎪⎪+ ~z m , 2 p1k p 3 k ) x3⎪~~~~~z p + z m , 3 p1k + ( z m ,12 p 3 k + z m ,1 p 2 k p 3 k + z m , 23 p1k + ⎪⎪⎪ m ,1 3 k⎪⎪+ ~z m , 3 p1k p 2 k ) x 2⎪~~~~~z p + z m , 3 p 2 k + ( z m ,12 p 3 k + z m , 2 p1k p 3 k + z m ,13 p 2 k + ⎪⎪⎪ m ,2 3k⎪⎭⎪⎩ + ~z m , 3 p1k p 2 k ) x1– 6-ти узловой КЭ{D }(m)~z ϕ + ( ~z~⎧⎫m ,1k ,1m , 13 ϕ k , 1 + z m , 1 ϕ k , 13 ) x 3⎪⎪~z~~m , 2 ϕ k , 2 + ( z m , 23 ϕ k , 2 + z m , 2 ϕ k , 23 ) x 3⎪~⎪~~~⎪ z m , 3 ϕ k , 3 + ( z m , 13 ϕ k , 3 + z m , 3 ϕ k , 13 ) x1 + ( z m , 23 ϕ k , 3 + ⎪⎪⎪1 ⎪ + ~z m , 3 ϕ k , 23 ) x 2⎪⎨ ~z ϕ⎬~~~k=8⎪m ,1k , 2 + z m , 2 ϕ k , 1 + ( z m , 13 ϕ k , 2 + z m , 1 ϕ k , 23 +⎪⎪ + ~z m , 23 ϕ k , 1 + ~z m , 2 ϕ k , 13 ) x 3⎪⎪⎪~z ϕ~m ,1k , 3 + z m , 3 ϕ k ,1⎪⎪~z~⎪⎩⎪⎭m ,2 ϕ k ,3+ zm ,3 ϕ k ,2;– 5-ти узловой КЭ{D }(m)⎧ ~z m ,1 ϕ k ,1 + ( ~z m ,123 ϕ⎪~~⎪ z m , 2 ϕ k , 2 + ( z m ,123 ϕ~~1 ⎪⎪ z m , 3 ϕ k , 3 + ( z m ,123 ϕ⎨k=8⎪⎪⎪⎪⎩z m ,12 ϕ k ,13 + ~z m ,13 ϕ k ,12 + ~z m ,1 ϕ k ,123 ) x 2 x 3 ⎫+~⎪~~~k , 2 + z m , 12 ϕ k , 23 + z m , 23 ϕ k , 12 + z m , 2 ϕ k , 123 ) x 1 x 3 ⎪~~~⎪k , 3 + z m , 13 ϕ k , 23 + z m , 23 ϕ k , 13 + z m , 3 ϕ k , 123 ) x 2 x 2 ⎪⎬~z ϕ + ~z ϕm ,1k ,2m,2k ,1⎪~z ϕ + ~z ϕ⎪m ,1k ,3m,3k ,1⎪~z ϕ + ~z ϕ⎪⎭m,2k ,3m,3k ,2;k ,1172– 4-х узловой КЭ{D(m )}~z ϕ⎧m ,1k ,1⎪~zm,2 ϕ k ,2⎪~z m ,3 ϕ k ,31 ⎪⎪⎨~k=8 ⎪ z m , 1 ϕ k , 2 + ~z m , 2 ϕ⎪ ~z m , 1 ϕ k , 3 + ~z m , 3 ϕ⎪~⎪⎩ z m , 2 ϕ k , 3 + ~z m , 3 ϕ⎫⎪⎪⎪⎪⎬k ,1 ⎪⎪k ,1⎪k ,2 ⎪⎭.Производные компонент тензора преобразования координат и функцийформы в центре КЭ вычисляются по формулам [164]:∂z m∂2zm~z m, j =z m ,α β =, ~, ~z∂ x α x = x = x =0∂xα ∂x βx1 = x2 = x3 =0123m ,1 2 3 =∂2zm∂ x 1 ∂ x 2 ∂ x 3 x = x = x =0123;2ϕ k ,α = ∂ ϕ k / ∂ x α x = x = x =0 ϕ k ,α β = ∂ ϕ k / ∂ x α ∂ x β x = x = x =0123,.1В выражениях23{D ( m ) } kдля 6-ти, 5-ти и 4-х узловых КЭ функции формыϕkвычисляются, соответственно, по формуле (4.3).Физические соотношения записываем в тензорном видеσ i j= E i jkl εгде σijkl,(4.4)- ковариантные компоненты тензора напряжений в местных осях;E i j k l - контравариантные компоненты тензора модулей упругости.
В дальi jklнейшем изменением Eв пределах элемента пренебрегаем, считая его равным~значениям в точке x m .Представим выражение (4.4) в матричной форме {σ } = [E] ({ε} − {ε 0 }) + {σ 0 } ,112233121323где вектор – столбец напряжений {σ } = {σ σ σ σ σ σ } ; матрица упруΤгости в общем случае имеет вид⎡ E 1111 E 1122⎢E 2222⎢⎢[E ] = ⎢⎢⎢⎢симметрично⎢⎣E 1133E 1112E 1113E 2233E 3333E 2212E 3312E 1212E 2213E 3313E 1213E 1313E 1123 ⎤⎥E 2223 ⎥E 3323 ⎥⎥;E 1223 ⎥E 1323 ⎥⎥E 2323 ⎥⎦173{ε 0 } и {σ 0 } - векторы - столбцы начальных деформаций и напряжений, задаваемые в местных осях.
Для термоупругой задачи вектор-столбец{ε 0 } = α Δ t {g 11g 22 g 33 0 0 0}T,где g ii , i = 1, 2, 3 - диагональные элементы метрической матрицы.EВыразим элементы матрицы [E ] через физические компоненты ( i j k l ) тензоE( i j k l )E i jkl=g ii g j j g k k g l lра модулей упругости(4.5)(по индексам в круглых скобках суммирование не производится). Здесьg ii , gjj, ...-ковариантные компоненты метрического тензора.Для изотропного и ортотропного материалов матрица упругости имеет следующую структуру:⎡ E 1111⎢⎢⎢[E ] = ⎢⎢⎢⎢⎢⎣E 1122E 113300222222330000E 12120EEE 3333симметричноE 13130 ⎤⎥0 ⎥0 ⎥⎥0 ⎥0 ⎥⎥E 2323 ⎥⎦.i jklДля изотропного материала компоненты Eопределяем через коэффици-енты Ляме μ и λ [260]:Ei j k l = μ (gi k g j l + gi l g j k ) + λ gi j g k l ,ijгде g - контравариантные компоненты метрического тензора КЭ.
Величины μ иλ связаны с техническими постоянными - модулем Юнга E и коэффициентомПуассона ν , следующими зависимостями:μ=EνE; λ=.2 (1 + ν )(1 − 2ν ) (1 + ν )Для ортотропного материала дискретизацию континуума осуществляем та~ким образом, чтобы оси x m КЭ совпадали с главными направлениями упругости.На основании принципа возможных перемещений (4.1) получаем выраже-174ния для матрицы жёсткости КЭ [139][h] = ∫ [D ] Τ [E ][D ] d v(4.6)vи вектора-столбца узловых сил{ f } = ∫ ( [F ] Τ{q} + [D ] Τ [E ]{ε 0 } − [D ] Τ{σ 0 } )d v + ∫ [F ] Τ{p}d s .v(4.7)sЗдесь обозначено:[ F ] = [ [ F ] 1 [ F ] 2 K[ F ] n r( 3×3 n r )– блочная матрица, субматрицы кото-рой имеют структуру [F] k = diag[ϕ k M ϕ k M ϕ k ] ;{q} = {q(1)q ( 2 ) q (3) }T,{p} = { p (1)p ( 2) p (3) }T– соответственно векторы-столбцыобъёмных и поверхностных сил, заданные в глобальном базисе (q (i ) , p (i )- физиче-ские компоненты);v – объём, занимаемый КЭ;s – площадь поверхности КЭ, к которой приложена распределённая нагрузка.Вычисление интегралов, входящих в выражения (4.6) и (4.7), осуществляемчисленно с использованием квадратурных формул Гаусса [260].Учитывая, что физические компоненты вектора напряжений наиболее точно вычисляются в точках интегрирования элементавые напряженияσ (i j ))σ(i j ), соответствующие узло-определяем с помощью процедуры “локального” сглажива-ния и последующего усреднения в узлах конечно - элементного ансамбля.Порядок вычислений: формируем вектор эквивалентной узловой нагрузкидля КЭ:){ p } = ∫ [ P ]{σ ( i j ) } d v( n r ×1)где матрицаv( n r ×8 )(8×1),[ P ] =[{ P } 1 {P} 2 ...{P } 8 ],{ P } k = {ϕ 1 ( x , x , x ) ϕ 2 ( x 1k , x 2k , x3k ) K ϕk1kk2k3nr( x 1k , x 2k , x3k ) } T= 1, 2,…, 8; - вектор физических компонент напряженийтегрирования.)σ (i j ),в точках ин-1754.3 Численные расчёты напряжённо-деформированного состоянияводопроводящих сооружений4.3.1 Результаты численных расчётов железобетонной облицовкиводопроводящих каналовРасчёты напряжённо-деформированного состояния водопроводящих магистральных и межхозяйственных каналов выполнялись на модели их характерногоэлемента − железобетонной облицовки (рисунок 4.2), рассматривающейся кактонкостенная оболочка, с учётом условия её опирания на грунтовое основание.Всё теоретическое обоснование математической модели характеризуются наоснове экспериментальных исследований по оценке надёжности сооружения приразличных сочетаниях разрушающих воздействий, выполненных на основе математического моделирования для установления степени опасности и допустимыхразмеров участков, существующих и возможных разрушений объектов, например,железобетонная облицовка рассматривается как тонкостенная оболочка, где учтены условия её опирания на грунтовое основание.34334266324158315030282718933510111914Y14145522591366306021134Z5112338612915773123168322425174133574973896588878685848317940645648728180797877767516839635547717372717069686715466237522012353284354634562645354364419922452765615544709082746669565246353657604743372658684851383449593942295067402625242322212019189XРисунок 4.2 – Конечно-элементная модель железобетонной облицовки действующегомагистрального канала176Конечно-элементная модель железобетонной облицовки действующего магистрального канала приведена на рисунке 4.2.