К вопросам построения эффективных алгоритмов расчета системы «Сооружение-грунт» Текст научной статьи по специальности «Физика»

К вопросам построения эффективных алгоритмов расчета системы «сооружение-грунт»

М.И. Кадомцев, А.А. Ляпин, С.И. Тимофеев

Ростовский государственный строительный университет, г. Ростов-на-Дону

При реализации методов расчета поведения поверхностных строительных объектов на основе совместного использования методов конечных (МКЭ) и граничных элементов (МГЭ) [1] важными вопросами являются построение эффективных схем расчета элементов напряженно-деформированного состояния полуограниченных структур, относящихся к многослойному основанию, и согласование методов в области их сопряжения.

1. Для расчета основания предлагается метод полуплоскостей: как при

определении полей от поверхностного источника, так и при нахождении фундаментальных решений для сосредоточенного источника с целью реализации МГЭ. Метод обладает высокой эффективностью, так как позволяет разделить волновые поля по их отражению от различных границ слоистой структуры, в том числе при анализе плотности потока энергии упругих колебаний.

Проиллюстрируем отмеченное на примере задачи плоской деформации многослойного полупространства.

Пусть область Б, занимаемая средой, представляет собой N -слойное упругое полупространство: Б = Д и Б2 и... и DN, описываемое в декартовой системе координат (х,у,2) как (рис. 1):

Рис. 1 - Область в декартовой системе координат Б = {х е (0,+ю), у, 2 е(- да,+да)} - полупространство;

Б] = {х е(- х.,-х.-1) у, 2 е(-гс,+гс)} х. = ^И; (¿1= 0) - >й слой 0=2,...^

/=1

толщины И] .

Упругие свойства сред в Б],] = 1,..., N описываются плотностью р . и коэффициентами Ламе X ., ^ . или соответственно модулем упругости Е. и коэффициентом Пуассона V. :

3Х ] + 2 и, ] X ]

Е. = ^ ^, V ,■ =---------]---.

] ] Х] + ^] 2(Х] ])

На поверхности среды в области О задана система распределенных усилий:

t(1) = ^, F = X(y) x = 0, y g W .

В случае однородной полуплоскости с применением преобразования Фурье по переменной y

f (^) = j f( y)exp (iay)dy, f( y) = — J f (a)exp (iay)da

—œ —œ

для функций перемещений точек данной области получим интегральные представления:

u(1) (x, y) = -1J P(1)(x, a) • X(a)e~iayda . (1)

Элементы матрицы P(1) имеют вид ( j = 1):

р1)(x, a) = 7^ Г Cfe1 + 2a^e2 }, P1(2)(x, a) = 4^ {“ 2a j1aj2e1 + j }

An A

\R

'R

P21)(Xa) = It j + 2C j1G j2e2 } P22} (Xa) - ^ {2a2el - С)e2 }

AR ar

Ar =-4a2G jiG j2 + d, ek - e U , к - 1,2 .

>j1Gj 2 t

f2 2 2 2 2 n.2 r. ® ^ ®

cj -a +aj2; Gjk -a -0jk; 0ji ; 0j2 •

V pj V sj

Аналогично для вектора напряжений на линиях х = const, найдем, что Т(1)(х, y) = co/|ax, т^}1 \х, у) = | W(l)(x, a) • X(a)e~“yda

wii) (Xa) -^~ fee - 4a jiG j2a2e2 }, wi2 ) (X a) -A r

( j )^.4_2zaC/2<C J

A

{e1 - e2 }

R

W2( j )( x, a) -

2/aa jiC j

A

{e1 - e2 },

w2(2j)( x, a) - {r 4a jiG j 2a2ei + r 4

XR

R

cb }.

(2)

k Ima

011 012 a R

Re a

-012 -011 ■ 'l

Рис. 2 - Контур интегрирования

Контур интегрирования Г в представлениях (1), (2) определяется применением принципа предельного поглощения [2]: при отсутствии диссипации энергии в среде обходит положительный корень уравнения Рэлея: АК (а) = 0 - снизу, отрицательный - сверху, а на остальной части совпадает с вещественной осью, как показано на рис. 2.

При наличии малой диссипации энергии в среде интегрирование можно проводить непосредственно по вещественной оси.

Решение для одного слоя при заданных на его гранях векторах напряжений:

10)(0,у) = ц] Код)(у), 1°’2)(И],у) = ц] Я0,2)(у)

строится способом суперпозиции решений для двух полуплоскостей.

Пусть в локальной системе координат для ] -го слоя: (х,у): х е (0,И]), у е(-да,да)

амплитудные функции перемещений имеют вид: и (])(х,у) = у),и2])( х у)}= {7х°) (x, уХи^х у)}.

Функции и(])(х, у), удовлетворяющие уравнениям движения, согласно предлагаемому методу будем разыскивать в виде суммы решений для полупространств

х > 0 - и(],1); и(],2) - х < И ]:

гг0') гг0’2)

и^ = и^ + и.

(3)

я'

( 7,1)

К(Л2)

Рис. 3

Здесь слагаемые и(],п)(х, у), п = 1,2 в (3) являются решениями уравнений Ламе для однородной полуплоскости с удовлетворением граничных условий:

1(]Д)(0, у) = ц ] Х(],1)(у), 1и,2)(И],у) = ц ] Х(],2)(у) .

Вектор перемещений и(],1'\х, у), представим через трансформанты вектора напряжений Х(],1)( у) в виде:

и а,1) (х, у) = — | Р(],1) (х, а) • Х(],1) (а)в~1ауёа . (4)

Здесь функции р(],1)(х,а) получены из Р(1)(х,а) заменой упругих параметров полуплоскости на параметры ] -го слоя.

Аналогично формуле (4) определяются перемещения для полуплоскости х < И] через функции Х(],2)(а), где для элементов Р<(Ц^",(х, а) справедливы соотношения:

Ршг’2\х,а) = (-^Р^И -х,а), п,т = 1,2, Ъпт - символ Кронекера.

Определяя напряженное состояние слоя в виде суммы соответствующих решений для двух полуплоскостей, получим:

2

(5)

Т(])( х, у) = ^ Т(],к)(х, у) к=1

где ^ х, а) имеют вид (2).

Для второй группы слагаемых найдем: ,,пт (х,а) = ( 1) ,,пт с

При рассмотрении далее общей краевой задачи для N -слойной полуплоскости

используются граничные условия и условия сцепления слоев между собой и

подстилающей полуплоскостью, что в пространстве преобразований Фурье по

Гп<т-2)(х,а) = (-1)5пт+1 шпт'ф, - х,а).

0

п

переменной у приводит к системе ЛАУ относительно неизвестных функций напряжений

По найденным компонентам напряжений на гранях слоя можно восстановить осредненный за период колебаний поток энергии, проходящий через границы раздела сред х =const :

Здесь, при введении малой диссипации в слоях конструкции в качестве контура Г+ выбиралась вещественная ось.

Подставляя далее выражения (3), (5) в преобразованном по Фурье виде в соотношение (6), получим:

Пki имеют смысл плотности потока энергии излучаемых и отраженных волн от плоской

поверхности х = const.

1. Важным моментом при численной реализации совместного использования методов конечных и граничных элементов для системы «сооружение-грунт» является выбор системы аппроксимирующих функций, а также применение формул численного интегрирования на элементах.

Стыковка МКЭ и МГЭ предполагает равенство векторов узловых перемещений и усилий в области контакта S фундамента здания или сооружения и окружающего грунтового массива и не требует согласования данных характеристик вне узлов. Отсюда выбор закона полиномиального распределения перемещений в области контакта может быть независимым для каждого из методов (безусловно, предпочтительнее выбирать аппроксимирующие функции одной и той же степени).

Не ограничивая общности, можно считать, что область контакта S является плоской с введением локальной системы координат (х, у). Соответствующая система

узлов £г- =(xi,yi) i = 1,...,N определяется разбиением конечной части системы

«сооружение-грунт» в методе конечного элемента и выбором типа конечного элемента.

Так для восьмиузловых твердотельных элементов (рис. 4) область £ разбивается на четырехугольные граничные элементы (1ЖЬ). При использовании опции элементов: призма или тетраэдр, граничные элементы имеют форму треугольников (ЦК). Таким образом, наиболее общей формой граничных элементов для применения МГЭ является

X(J,k)(а), j = 1,2; к = 1,2,...,N.

(6)

Syitcri

КеУО:

X

Рис. 4

треугольная. К этому же приводит разбиение области контакта неплоской формы путем триангуляции соответствующей поверхности в пространстве.

Таким образом, считаем, что область £ разбита сетью граничных треугольных элементов с узлами (1,2,3). В каждом узле вектор перемещений и(х, у) имеет значение и«, / = 1,2,3.

2

Рис. 5

Применим на элементе линейную аппроксимацию:

3

и( х у)=2 ф« (х у)и-, ф/(х у)=а«х+ь«у+с (7)

/ =1

Константы а, Ь, в1 определяются из условий ф« (х^, у^ ) = 8^, 8 у - символ Кронекера. В результате получим:

а1 =

у2 - у3

А

х2 - х3

а2 =

У ~ У3 А

с1 =

х2 У3 - х3 У2

х1 ~ х3

А

с1 =

А

х1У3- х3 У А

(8)

а3 =

У1 - У2 А

А

с3 =

х1 У2 - х2 У1 А

А = х1(у2 -у3) + х2(у3 -у1) + х3(у1 -у2) .

Следует отметить, что линейная интерполяция неизвестной функции на каждом элементе не нарушает условия непрерывности поля перемещений в целом на границе £. Это следует из условий равенства полного набора констант аппроксимации для области £ ( 3М, где М - число граничных элементов модели) и суммы числа условий непрерывности перемещений в узлах для смежных элементов и числа самих узлов как точек коллокации для определения неизвестных перемещений в них.

При использовании аппроксимаций более высокого порядка, например квадратичной (лагранжевы элементы), условия согласования перемещений в узлах и их корректное определение требуют введения дополнительных узлов по центру сторон треугольников. Это в свою очередь приводит к необходимости использования в сопрягаемом МКЭ 10-узловых пирамидальных элементов (рис. 6), что увеличивает порядок системы линейных уравнений метода конечных элементов и сложности стыковки с МГЭ. Получаемое же увеличение точности решения легко можно компенсировать уменьшением сетки разбиения при использовании линейной интерполяции на элементах.

Рис. 6

Решение граничного интегрального уравнения требует также аппроксимации в области £ вектора напряжений:

р = о • п, где о - тензор напряжений Коши, п - нормаль к области £ . Исходя из требования сохранения количества узлов сопрягаемых сеток граничных и конечных элементов, для вектора напряжений необходимо применять интерполирующие функции того же порядка, что и для вектора перемещений.

Форма сопряжения МГЭ и МКЭ по напряжениям в узловых точках в этом случае может иметь следующий вид:

и = X р* —

о,

веБ

К

здесь О, - узловое усилие в , -м узле МКЭ;

р, - значение вектора напряжений в , -м узле МГЭ;

Б - множество граничных элементов, сопряженных с , -м узлом;

1

8в - площадь элемента с номером в (рис. 5), 8в =

2

х1 - х3 х2 - х3 ^

ч У1 - У3 У2 - У3 )

К - количество узлов граничного элемента (в случае линейной интерполяции К = 3 , для квадратичной - К = 6 ).

При использовании метода граничных элементов необходимым элементом является интегрирование по двумерной области с разбиением на треугольные элементы:

АЛВС

где подынтегральная функция /(х, у, Е, л) может иметь интегрируемую особенность степенного или логарифмического характера в точках (Е, л) , совпадающих с узлами аппроксимации (вершинами треугольника или серединами его сторон). В этом случае одним из способов интегрирования, показавшим существенную эффективность, является использование квадратурных формул с узлами внутри треугольной области.

Отобразим треугольник ЛВС (рис. 5): А(х1,у1), В(х2,у2), С(х3,у3), на

равносторонний треугольник РОЯ в системе координат (и, V) (рис.7) с использованием линейного преобразования:

Г х(и, V) = аи + bv + с [ у(и, V) = ёи + ву + /

Рис. 7

Несложно получить, что

а =

2 — Л*2 — ^ Х2 — Х3 + Х2 + Х3

а = 2 У1— У2—У3, е =

■Л ’

У2 — У3

73

с = ■

, /=

3

у1 + у2 + .Уз 3

В результате имеем:

Я / (х. У, 5 , л)ахау = |з| Ц / (х(и, V), у (и, V), 5, л)амОу,

ААВС АРдЯ

2л/3

9

х1 — х2 х2 — х3

у1 — у2 у2 — у3

якобиан перехода к системе координат (и, V)

Для вычисления двойного интеграла по треугольнику РОЯ воспользуемся 7-узловой квадратурной формулой:

Ц /(u, v)dudv = т^ I]

рдя 3^3 і=1

^і/(иі ,vi) + Я,

Ардя

где весовые коэффициенты ^ и узлы (и,, V,) приведены в таблице 1.

Таблица 1.

3

3

І Wi иі V'

1 270/1200 0 0

2 155 — л/Г5 л/Ї5 +1 0

1200 7

3 155 — -Л5 — л/Ї5 +1 Л5 + ^

1200 14 14

4 155 — 715 — •>/15 +1 715+^

1200 14 14

5 155 + л/Ї5 л/15 — 1 0

1200 7

6 155 + л/Т5 л/15 -1 Л5-^

1200 14 14

7 155 + 715 л/15 -1 Л5 -1^

1200 14 14

Данная формула имеет 6 порядок точности: Я = J|3) [3].

Литература

1. Кадомцев, М.И. Исследование динамики заглубленных фундаментов методами граничных и конечных элементов / Кадомцев, М.И., Ляпин, А.А., Селезнев, М.Г. // Строительная механика и расчет сооружений. - 2010. - № 3. - С.61-64.

2. Бабешко, В.А. Динамика неоднородных линейно-упругих сред / Бабешко, В.А., Глушков, Е.В., Зинченко, Ж.З./ - М. : Наука; Главная редакция физикоматематической литературы. 1989. - 343 с.

3. Справочник по специальным функциям /Под ред. М.Абрамовица и И.Стиган/ -М.: Наука, 1979. -832 с