Анализ обобщенного оператора ламе и

1. Введение Согласно натурным экспериментам (, , , и др.), на расстоянии от дневной поверхности более одного метра для глины и полутора метров для торфа имеются остаточные избыточные поровые давления после окончания процесса консолидации водонасыщенного грунта, то есть в стабилизированном состоянии, не зависящем от времени. Моделями теории упругости и пластичности (, -Посадов, , , и др.) это поровое давление не описывается, потому что грунт рассматривается однофазным. По линейным фильтрационным моделям (К. Терцаги, , М. Био, , , -Марти-росян и др.) остаточные избыточные поровые давления обращаются в нуль. Эти модели описываются системами уравнений параболического типа. В нелинейных фильтрационных моделях (, и др.), в которых учитывается начальный градиент порового давления, вводятся две зоны — активная и пассивная. В первой зоне остаточные избыточные поровые давления обращаются в нуль, а вторая зона считается абсолютно твердым телом. В статье рассмотрена модель водонасыщенного грунта, в которой избыточные остаточные давления описаны с помощью нового варианта закона уплотнения грунта, не зависящего от времени, Pli = Eliвl8ij, где Pli = hiДjаl8ij, 8ij — символ Кронекера. Новые уравнения взаимосвязи между относительными компонентами деформаций по Коши твердой и жидкой фаз гs = -К jеl8ij позволяют использовать только деформации твердой фазы. В результате получается система уравнений эллиптического типа, как это принято в механике деформируемого твердого тела и, в частности, теории упругости (уравнения Ламе). Однако в уравнениях типа Ламе за счет наличия жидкой фазы появились младшие производные, что привело к несимметрии оператора и вызвало необходимость адаптировать метод конечных 181 элементов (МКЭ). Аналоги эллиптической модели, по мнению авторов статьи, в литературе отсутствуют. Модель подтверждена новыми лабораторными (, , ) и натурным () экспериментами. Обобщенный оператор Ламе представлен в виде суммы A +B+C, где A — известный симметричный оператор Ламе. Симметричный оператор B включает в себя новые механические характеристики bi, которые отражают изменение механических свойств скелета грунта за счет наличия поровой воды. Оператор C, образованный производными первого порядка, описывает разгружающий вклад поровой воды и является несимметричным, то есть (Cu,v) Ф (u,Cv). Энергетические произведения типа ((A +B)u,v) введены для симметричных операторов, поэтому произведение (Cu,v) назовем квазиэнергетическим. В теории упругости формулы Бетти и Клапейрона получены на основе энергетического произведения (Au, u). В статье аналоги этих формул получены путем применения квазиэнергетического произведения ((A+B+C)u,u). Третий вариант новой формулы Бетти представляет собой разность двух квазиэнергетических произведений: ((A + B + C )u, v) — (u, (A + B + C)v) и описывает несимметричность обобщенного оператора Ламе. Вклад оператора C в удельную работу деформации отражается в новой формуле типа Клапейрона. Новый конечный элемент введен на базе квазиэнергетического произведения. 2. Постановка краевой задачи Для описания напряженно-деформированного состояния двухфазного тела в стабилизированном состоянии относительно перемещений только твердой фазы u получена система дифференциальных уравнений второго порядка с положительными постоянными коэффициентами G, ?, b, ci (i = 1, 2, 3), описывающими механические свойства (модуль сдвига, постоянная Ламе) твердой (индекс опущен) и жидкой (индекс l) фаз: Dijuj=Fi, bi=Eli/?2 i , ci =Eli/hl ? i, (1) Экспериментальное определение механических характеристик жидкой фазы Eli (МПа), безразмерных положительных коэффициентов К , определяющих долю относительных деформаций скелета грунта от относительных де формаций жидкой фазы, и геометрические параметры hi (м) приводятся в работе . Через Fi обозначены компоненты объемных сил. Дифференциальный оператор представим в виде D = -(A+B+C ), Aij=(G + X)8i8j+G8ij8k8k, Bij =bii^dj, Cij=ci8ij8j. Знак минус перед суммой операторов введен потому, что положительно определенным является отрицательный оператор. Здесь A- оператор Ламе, оператор B симметричен, механические постоянные bi описывают изменение трех диагональных элементов в тензоре четвертого ранга упругих постоянных за счет наличия поровой воды: 182 [Ea+b\ fX + 2G + b1 X X 0 0 0 X X + 2G + b2 X 0 0 0 X X X + 2G + b3 0 0 0 0 0 0 G 0 0 0 0 0 0 G 0 0 0 0 0 0 G Оператор C моделирует разгружающее влияние жидкой фазы и несимметричен. Дополним уравнения смешанными граничными условиями: tij = nij+ (G + bfiij )nfl + G*ijnk k, t — оператор внутренних напряжений в скелете грунта записан на основе тензора , заданная нагрузка q приложена к поверхности тела, п — внешняя нормаль к поверхности S = S1+S2. Касательные напряжения по модели действуют только в скелете грунта. Поскольку на дневной поверхности тела S2 за счет дренирующего покрытия избыточные остаточные поровые давления обращаются в ноль, то остаются только нормальные напряжения в скелете грунта, поэтому статические граничные условия (2) записаны по аналогии с записью, применяемой в теории упругости. В приложениях под областью V понимается часть полупространства, нижняя граница S1 которой есть сфера большого конечного радиуса, дневная граница S2 является плоскостью, поэтому один направляющий косинус равен -1, а два других обращаются в нуль. 3. Анализ обобщенного оператора Ламе Анализ основан на свойствах квазиэнергетического произведения (Du,u ) . Для преобразования интегралов по объему применяется формула Остроградского-Гаусса. Первое слагаемое известно : (-Au, u) = -\u\Aijuj dV = 2\WA (u,u )dV — \u\lijuj dS, V V S WA(u,u ) = WA(u ,u), lij =Xnidj+Gnjdi +G5ij nA. При u = u получаем известный упругий потенциал WA(u) = 1 (XQ2+2Geijeij), 9 = 6ij8ij. Два оставшихся слагаемых имеют вид : B (-Bu , u) = -\u\BjujdV = 2\WB(u,u )dV — \biu^njdiuj dS, V V S (3) (-Cu,u?) = — ? ui ?Ciju dV = ? W C (u,u)dV =- ? ciu ?? iju ni dS, VV 2 S 183 WB (u,u ) = 1 bi, WB (u,u ) = WB (u ,u), WC(u,u ) = -cieiiu;, WC(u, u)*WC(u,u). Форма WC(u,u ) не обладает свойством коммутативности аргументов Операторы A и B можно объединить на основе тензора Гука . При u = u объемный интеграл от WC в (3) положителен, так как на S2 один направляющий косинус внешней нормали отрицателен. Аналог первой формулы Бетти описывается выражением ju Dijuj dV = 2J\ W A (u,u ) + W B (u,u ) + 2 W C (u,u) )dV-ju tijuj dS. (4) Полагая u = u, получим аналог второй формулы Бетти \uiDijuj file:///uiDijuj dV = 2UW A + W B +-W C )dV-\uitijuj dS. (5) Оставляя только упругий потенциал WA , приходим к известным формулировкам двух формул Бетти. Покажем несимметричность квазиэнергетического произведения. Вычтем из (4) взаимное выражение и запишем аналог третьей формулы Бетти \(u\Dijuj — ui Diju j )dV = j(WC(u,u ) — WC(u ,u))dV -\(u itijuj — ui tiju j )dS. (6) VV S При введении объемных сил по уравнению (1) и учете статических граничных условий (2) имеем (u\Fi -ui Fi )dV-(u i qi -uiq )dS=( WC(u,u )-WC(u ,u))dV. (7) V S2 V При отсутствии объемных сил и задании в двух точках поверхности S2 сосредо-точенных сил q(x) = Q 5(x — y0), q (x) = QЪ(x — y 1) получим u (y0)-Q = u(y0)-Q +(WC(u, u)-WC(u ,u)) dV. (8) V По теореме Бетти о взаимности работ для упругого тела объемный интеграл отсутствует, поэтому несимметричность оператора D описывается этим интегралом. На основе уравнений (1) и условий (2) введем объемные и поверхностные силы, тогда на основе выражения (5) будем иметь аналог формулы Клапейрона uiFidV + S uiqidS = 2 V(W A +W B + 1 W C)dV, в котором правая часть описывает энергию деформации двухфазного тела в виде квадратичных и билинейного функционалов. При сохранении только WA получим в правой части потенциальную энергию и известную запись теоремы Клапейрона. 4. Основная идея МКЭ Умноженные скалярно уравнения (1) на вектор возможных перемещений u приводят к выражению (Du, u) = (F, u), в котором по принципу Лагранжа работа 184 внутренних сил Du на возможном перемещении u равна работе внешних сил F на том же перемещении: (Au,u) + (Bu,u) + (Cu,u) = (F,u), 2J»(A,92 + 2Gzijzij)dV + 2|bi8ij8ij5ij dV + (Cu,u) = JFiui dV — qiui dS. V V V Для известного оператора Ламе A имеем : qiui dS, VS поверхностный интеграл употребляется для граничного конечного элемента. Второй объемный интеграл объединим с первым на основе тензора Работа деформации внутренних сил в жидкой фазе в пределах объема треугольного элемента равна (Cu, u). Для случая плоской деформации с использованием коэффициента Пуассона V для твердой фазы тензор Гука имеет вид: 1- v 2G = — + b 1 2G 0 v X V 1-2v 1-V+b2 0 0 0 G (9) Поскольку МКЭ является численным методом, то требуется матричная запись каж-дого слагаемого, которая для оператора A известна . 5. Построение конечного элемента, моделирующего избыточные остаточные поровые давления Согласно , введем стандартный конечный элемент треугольной формы пло-щадью ?. Fm 2 t x u1 m x 1 ]u 2 j Fm 1 Fix 2 u 1 i 1 u 1 j* Внутри треугольного элемента перемещения u = (u u 2) записываются через узловые перемещения (u 1 i u2 i u 1 j u2 j u 1 m um) по формулам : uk = 1 [(pi +dix 1 +nix 2)uki+(pj +djx 1 +njx2)ukj + nmx2 +(pm+dmx1+nmx k )ukm] (10) 185 (11) di=x2 j-x2 m, k = 1,2. pi=x 1 j x 2m-x 1mx2 j, ni=x 1m-x2 j, Узловые перемещения являются искомыми величинами и образуют матрицу-столбец {6}. В дальнейшем фигурные скобки обозначают либо матрицу-столбец, либо матрицу-строку; квадратные скобки — полную матрицу. На основании уравнений Коши: дu 1 1 *11= u 1 = 1i( d iu1 i+dju1 j+dmu1 m ), в22 = 1 еu 2 = 1 (niu i2+nju2 j+nmu2 m), У12= u 2 + u 1 = 1 (niu1+nju 1+nmu 1 +diu2 i+d ju2 j+d mu ), относительные деформации внутри конечного элемента запишем через {6}: i {г} = {8} или О 2 в 12 V 12 У = 2Л di 0dj 0 dm 0 0ni 0nj 0 nm ni di j dj nm dm u1 V2 J где матрица образована на основе выражений (11). Удельная работа внутренних сил для скелета грунта равна 2(WA + WB) = {s} T{o}, {6} T={8} T T, индекс T обозначает операцию транспонирования. От удельной работы перейдем к работе внутренних сил в пределах объема элемента, интегрируя по площади элемента А и умножая на единичную толщину: j{s} T{a}dS = {s} T{a}.A.1, jdS = A = 1x1i 1 x1j 1x1m x2i x2j . Поскольку {е} = T и {о} = {г}, то {s} T{o}-A = {b} T[ Nm-A. Работа внешних сил, приложенных в узлах элемента, равна {d} T. Равенство работ внутренних и внешних сил имеет вид: {8}T[Nm.A + (Cu,u) = {8}T{F}. (12) Ниже сомножитель {8}T будет сокращен, тогда в левой части уравнения первое слагаемое примет вид: A.{S}. Сомножители перед искомым вектором перемещений {6} образуют матрицу жесткости для скелета грунта или 186 А = 4Л di 0 d j 0 dm 0 0 ni 0 nj 0 ni ni di n j d j nm dm к +b 2G у 1-2v 0 = v 2G 1-2v +b2 Я 1-у V 0 0 0 G 0 ni 0 ni di d j 0 nj 0 nj d j dm 0 nm 0 nm dm В матричном виде для плоского случая перепишем квазиэнергетическое произве-дение дu2 сu1 )dS. u1+c 2 u2 (-Cu,u) = 1\ На основании нового варианта закона уплотнения внутренние удельные силы Pl =(n1?11, n 2822), вызванные избыточными остаточными поровыми давлениями в жидкой фазе, описываются произведениями c 1^; c 2^2. Выражение c 1z11u1 + c2e22u2 представляет собой удельную работу деформации жидкой фазы, возникшую за счет избыточных остаточных поровых давлений. При переходе от удельной работы к работе внутренних сил в пределах объема элемента возникают определенные интегралы f x 1 dS и f x2dS. Рассмотрим взятие одного из интегралов на основе известной теоремы о статическом моменте площади, согласно которой статический момент равен произведению площади треугольника А на расстояние xc до центра тяжести треугольника, совпадающего с пересечением его медиан: + x2 j+x2 m \x 1dS = x1c-K x x 1c . = x 1 i+x 1 j+xm , x2c = 3 2c 3 Вместо координат x 1 и x2 используем координаты центра тяжести x1c и x2c. Перемещения внутри элемента в новых координатах описываются формулами uk = 1 [(pi +dix 1c+nix2c ) uki +(p +dx 1+njx2c)ukj + 2Л~ + (pm+dmx 1 c+nmx2c)ukm ]. Введем новую матрицу с элементами fi=pi+dix 1c+nix2 c , fj=pj+djx 1c+njx 2 c , fm=pm+dmx1c+nmx2c, тогда для перемещений внутри элемента имеем матричную запись: u1i \u2J 0f u = {5}= 1^ 0 fj 0 fm 0 j m М0 fi 0 f = {5} T T. u2 J 187 В результате квазиэнергетическое произведение примет вид: (-Cu,u) = ? {Pl}{ ? }, {Pl} = (c 11 ? 11 c 22 ?22). Согласно левой части формулы (12), сначала используется сомножитель {8}T , поэтому поменяем местами сомножители T и перенесем сомножитель {Pl}: (-Cu,u)=А{5}T}T{Pl}T. Необходимо согласовать выражение T { Pl } T с выражением {8} из формулы (12). За счет отсутствия касательных напряжений в жидкой фазе матрица-столбец {Pl } T по размерности на единицу меньше матрицы {с}. Дополним матрицу {Pl } T нулевым элементом и введем нулевые элементы в матрицу , что позволит переписать {Pl } T сначала в виде произведения двух матриц, а затем при учете {?} = {8} — в виде произведения трех матриц: c 1 ? 1 11 c 1 0 0??11 {Pl } = = = 0 c2 0 ? 22 0 0 0 0JU T = 2? За счет нулевых элементов в матрице T добьемся одинаковой размерности с T (di 0 nЛ fi 0 0^ 0 ni di0 fi 0 dj 0 nj 0 nj dj, T = 2?fj 0 0 0 fj 0 dm 0 nmfm 0 0 l0 nm dmjV0 fm 0j Квазиэнергетическое произведение получило матричную запись: (-Cu,u) = ? {?}T{?}, что позволяет переписать уравнение равенства работ внутренних и внешних сил (12) в двух вариантах: {? } T (? + ){?} = {F }, и ввести матрицу жесткости для жидкой фазы или 0dj 0 d0 0ni 0 nj 0 ndn d n № 4А (fi00Г 0fi0n00г fj000n 0 0f j000 0 f m00V 0 Vfm0 188 С помощью суммарной матрицы жесткости записана система линейных алгебраических уравнений МКЭ для нахождения искомых узловых перемещений {?}. Описанная методика построения матрицы жесткости для двухфазного элемента может быть перенесена на прямоугольный и другие конечные элементы. Литература Мальцев, модель грунта и биоматериалов / , , . — СПб.: Стройиздат, 2002. — 320 с. Мальцева, двухфазного тела с учетом несущей способности жидкой фазы / , // Математическое моделирование. — 2004. -Т. 16, №11. — С. 47-60. Михлин, методы в математической физике / . — М.: Наука, 1957. — 476 с. Мальцева, функционала для решения обобщенной системы уравнений Ламе / // Вестник ТГУ. — 2003. — №5. — С. 196-201. Победря, методы в теории упругости и пластичности / . — М.: Изд-во МГУ, 1981. — 344 с. Тимошенко, упругости / , Дж. Гудьер. — М.: Наука, 1975. -575 с.