CC BY
20
1/2008
ВЕСТНИК
МГСУ
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ УПРУГО-ПЛАСТИЧЕСКОЙ РАБОТЫ СИСТЕМЫ: ШТАМП - ГРУНТОВОЕ
ОСНОВАНИЕ
Постановка задачи
Исследуется несущая способность грунтового основания под действием штампа, загруженного вертикальной силой Р и горизонтальной силой Н, в условиях плоской задачи (плоской деформации).
Грунтовое основание представляется в виде упруго-пластической среды с поверхностью текучести:
где т и а - касательные и нормальные напряжения на площадке сдвига;
тпр - сдвиговая прочность грунта. Штамп моделируется упругим материалом.
На контакте штампа и грунта основания моделируется возможность раскрытия (или закрытия) контактного шва и упруго-пластический характер его работы.
Задача решается на основе метода конечных элементов (МКЭ) с применением разработанных компьютерных программ.
Математическая модель грунта основания
Зависимость т. р = / (V) устанавливается на основе экспериментальных иссле-
дований с образцами грунта и может быть аппроксимирована нелинейной зависимостью по параболической кривой (предложение Фейерхерста, рис. 1а) либо линейной зависимостью Кулона-Мора (рис. 16).
Орехов В.Г., Толстиков В.В. (МГСУ)
(1.1)
а)
(*-<П \
*
г
Рис.1 Огибающие кривые: а) - парабола Фейерхерста; б) - зависимость Кулона - Мора.
Зависимость Фейерхерста
Нелинейную зависимость Тпр = У(<г) Фейерхерст предложил аппроксимировать параболой у2 = 2рх; где р - фокальный параметр. Зависимость Укак огибающая предельных кругов Мора, в этом случае представляет-
п р
ся в виде:
х2р = 2р{яр -4
(2)
где Яр - прочность материала на растяжение. Параметр р определяется из рассмотрения предельного круга Мора при одноосном сжатии (а! = 0; а3 = Яс). Откуда:
Р = —р
или обозначая
+!
V — У
—
■2 + ! + !
я
- —; т = Vп + ! = I— + !
Яр
п -
получим
р =
т - !)2 2
я
• —
(3)
(4)
(5)
и,следовательно,
*2 р =(т - !)2 —Р (—Р"4
Указанная зависимость в величинах главных напряжений а! и а3 записывается в
виде:
(в 1 3 )2
+а з
-2(т -1)2 Яр
1 + -
а 1 +а з
т -1 2
-1
где а1 - максимальное главное напряжение;
а3 - минимальное главное напряжение (сжатие имеет знак «-»). Зависимость (2.6), как огибающая кругов Мора, справедлива при условии
г(т - 2)ах + а3 (0;
т(
для условия
т
(т - 2)а1 + а3 > 0
прочность материала ограничивается величинои
^ = Яр.
(6)
(7)
(8) (9)
2
В зоне т < Тр принимается упругий характер работы материала. Идеализированные графики работы материала представлены на рис. 2 и рис. 3.
ту;
т
-р-о
>2
-
/
/1
ю
ъ-р
Л? о
в -о
Рис.2. Идеализированные зависимости осевой деформации е3 от напряжений <г3 при различном
боковом обжатии р:
а - идеально пластическое поведение материала; б - разуплотняющийся материал с бесконечной хрупкостью; в - разуплотняющийся материал с конечной хрупкостью
Рис.3. Графики, характеризующие свойства разуплотняющейся упругопластической среды в условиях плоской деформации
На рис. 3 показаны графики, характеризующие работу разуплотняющейся упруго-пластической среды в условиях плоской деформации: зона I - область упругой работы материала; зоны 11-У - запредельная область деформирования материала; кривая А-В-С-Е-Э - предельная поверхность текучести; кривая М-Е-Ы - остаточная поверхность текучести.
Характер работы материала в запредельной области деформирования в явном виде не выражается и отыскивается на основе алгоритма, реализующего итерационный процесс решения задачи.
Модель Кулона-Мора
Зависимость Кулона-Мора записывается в виде:
р = с<Р (10)
или через главные напряжения а1 и а3 в виде:
-03) + (о*1 + ^3 ^Бтф — 2еСовф = 0 . (11)
Указанная зависимость справедлива для случая О^ (Яр; в противном случае
происходит разрушение материала при О^ > Яр .
Решение задачи в запредельной области производится на основе итерационного процесса.
Моделирование нарушения сплошности основания
При решении задачи по методу МКЭ для моделирования нарушения сплошности материала на контакте сооружения с основанием, а также в заранее известных областях массива основания, где возможны проявления указанных нарушений, применяются специальные контактные элементы.
Как указывалось выше, контактные элементы позволяют моделировать возможность раскрытия (или закрытия) рассматриваемого сечения расчетной области и упруго-пластический характер работы материала в указанных сечениях.
В модель поведения контактного элемента входят следующие параметры
а) при действии нормальных напряжений (а):
Яр - прочность на растяжение:
Т'тах - максимальное значение закрытия элемента:
кп -удельный модуль нормальной деформации (нормальная жесткость);
б) при действии касательных напряжений (т):
Тр и Т— - предельное и остаточное значения сдвиговой прочности элемента;
к3 - удельный модуль сдвиговых деформаций (касательная жесткость).
В модели поведения контактного элемента учитывается также а0 и т0 - начальные значения нормальных и касательных напряжений в элементе, обусловленных начальным напряженным состоянием массива основания, или их значения на предыдущем этапе при поэтапном нагружении сооружения.
Указанный характер поведения контактного элемента воспроизводится в расчетах с помощью итерационной процедуры.
Итерационная процедура
На первом шаге итерационного процесса перемещения и напряжения вычисляются по характеристикам упругой части зависимостей (Гп — /(V) и Т3 — /(тл) .
Рассмотрим случай действия растягивающих напряжений (Сп )0). Если растягивающие напряжения (Рис. 4, правая часть графика зависимости
превышает предел прочности на растяжение ^, контакт нарушается и образуется
раскрытая трещина, которая уже не несет никакой нагрузки. Это приводит к перераспределению напряжений в окружающих элементах. Такое перераспределение достигается путем вычисления узловых сил, ликвидирующих напряжения в трещине, которые прикладываются на следующей итерации и являются, по существу, новым случаем нагружения.
Очередное распределение напряжений дает новое значение нормального напря-
жения СТ
^ , которое равно
а
{2) = кп^(2) +^0
(12)
аналогично касательное напряжение Т3^ для данного случая загружения [ап)Яр )
равно
Г5(2) _ к8и{2) + Т0 + А ^(1) .
(13)
Невязки напряжений вычисляются как разница между упругими (точки 1', 2', 3', 4' и т.д. рис. 4) и истинными (в данном случае нулевыми) напряжениями при достигнутых перемещениях:
п(2 )=~кпУ(2 0;
5(2) = _к5и(2) ~Х0.
(14)
Рис.4.Схема итерационного процесса а =А(у)
Для достижения истинного напряжения требуется несколько итераций. Зависимости (12^14) можно записать в общем виде для «1 - ой» итерации
°п{1 )= knv (г) +^0 А°п(, )=- knv (г)
(г) = к5Ы{г) + Т0 + (г-1) ;
(г) = ~к8и(г) _ Г0.
(15)
(16)
Рис.5.Схема итерационного процесса т =1(и)
Сходным образом учитываются нелинейные свойства среды при сдвиге для случая, когда Сп (Яр и Т$ )Гр (рис. 5). Невязка касательных напряжений вычисляется в
этом случае по следующим зависимостям: -при «правостороннем» сдвиге (т3 >0)
(г) _ ТЯ к8и(г) Г0;
-при «левостороннем» сдвиге (тх (0 )
(г) _ ТЯ к8Л(г)
(17)
(18)
где Т— — С— ~ О- остаточная прочность при сдвиге.
1/2008_иГВЕСТНИК
Если в процессе деформирования перемещения превышают заданную величину ^тах (Рис. 4, левая часть графика зависимости ип = f (у)), то на каждой итерации при достигнутом уровне напряжений С^ корректируются относительные перемещения
. Невязка напряжений при этом равна:
Аап(,п(,)- кпУтах 0. (19)
Напряжения (ДСТц^ + СХ0) и (Аг^^^ + Г0)являются начальными напряжениями для следующей «1+1» итерации. Вектор узловых усилий ^ | элемента, подсчитанный
по значениям этих напряжений, добавляется к вектору сил системы и производится следующее упругое решение с прежней матрицей жесткости, но с новым набором узловых сил. Добавление сил, обусловленных начальными напряжениями, увеличит упругие напряжения в элементе на следующей итерации, однако на величину меньшую, чем начальные напряжения, по которым были рассчитаны узловые силы, поскольку в ансамбле элементов добавленные силы распределяются также и на другие элементы расчетной области. В физическом смысле это означает итерационный поиск таких дополнительных нагрузок, которые сообщают линейно-деформируемому телу перемещения, равные перемещениям нелинейно-деформируемого тела при заданной нагрузке. Итерации продолжаются до стабилизации решения.
Алгоритм расчета для упруго-пластической модели грунта по условию Кулона-Мора и условию Фейерхерста
Условие Кулона-Мора
В расчетах воспроизводится поэтапность нагружения системы «сооружение-основание». На первом этапе нагружения определяется вектор узловых сил системы
{-Р} , обусловленный напряженным состоянием грунта основания от начальных напряжений (вызванных действием собственного веса грунта, обозначаемое вектором напряжений {<7*01). На следующих этапах нагружения в вектор узловых сил системы включается
также вектор } , соответствующий действию внешних поверхностных и сосредоточенных сил и вектор , учитывающий невязку напряжений |До"| при развитии пластических деформаций в грунте.
В общем случае нагружения вектор узловых сил системы {-Р} записывается в следующем виде:
{Р}, = {Р» }, +{Р} + {РР}1.1. (20)
Рассмотрим в общем случае решение нелинейной задачи на «1» шаге итерационного процесса.
1. Решается упругая задача
{Я }=[Ктп ]{и}, (21)
где \Ктп j -матрица жесткости системы; |и| - вектор узловых перемещений системы.
2. Далее для каждого конечного элемента системы выполняются следующие операции (п.п. 3^10)
(У У У | ^ y,a xy ) и расчетных значений, с учетом пластических деформаций, (р'Р xy )
}=[ßS ]-{u} + {a 0};
a y(-'
(22)
'} = {а у } + {Ла г ^ где j - матрица связи напряжений в элементе с перемещениями его узлов. 4. Вычисляются значения главных упругих напряжений |(7 у | и расчетных :
- !а Р
пряжение Г-' гл
а = 0,5 а = arctg
-+ ау )±v(a - -аy) +(2х)2
2т
(23)
а - а
- У
5. Назначаются расчетные значения главных напряжений |С) j , Q ^ j :
если aj > Rp, то af = 0;
ay <Rp, то a? =aУ;
a 3 > Rp, то a p = 0; (24)
a^ <Rp, то aP =a3.
6. Проверяется условие текучести и определяются теоретические значения глав-
„ / Т т\ ных напряжении IQ j , Q 3 I
R -a i (1 + Sino) -а 3 (1 - Sina) - 2 S0, где S0 - c ■ Cosa.
(25)
Если R < 0, то
Если Я >0, то
ст„ = 0,5
ст[ = *?; а I - а >.
(ар + аР ) + (аР - а£ )• Япа
_ Т _ СТ и + $ °1 _
_ г _ СТ п - $0
°3 _
1 + £ша 1 - £ша
7. Определяются теоретические значения осевых напряжений
(26)
(27)
(28)
СТ 1 = 0,5 [К + ст [ ) + (ст[
- У = 0,5 [К + а 1 )-(а[
X г = 0,5(а[- -а з V &'п2а
(29)
8. Вычисляются невязки напряжений (Асх , Дс , Дх)
{Да} = {а ' }-{<т•}.
(30)
9. Формируется вектор дополнительных узловых сил элементов от невязки напряжений Ас^:
лт
(31)
где - матрица дифференцирования перемещений.
10. Определяется вектор узловых сил системы для следующего «1+1» шага итерации:
},+1 = № },+1 +{Р} + {РР}Г (32)
11. Далее идти кп.1.
Итерационный процесс ограничивается стабилизацией величины перемещения характерной точки системы
[и - и<Д!, (33)
либо стабилизацией величины напряжений
к),-(«),1 ](А 2, <«>
где Д1 и Д 2 - величины допускаемой погрешности расчета.
Условие Фейерхерста
При нелинейной зависимости X пр = /(ст), которая по предложению Фейерхерста аппроксимируется параболой
У
■ 2рх или х2р = 2р[Ир - ап)
(35)
в приведенном выше алгоритме п.6 записывается в следующем виде:
6. Определяются теоретические значения главных напряжений (о ^, С ^ ^ .
6.1. Вычисляется значение Cos2fi (см. рис. 1,а) из решения уравнения
(а? + )-(2—р + р_
а ?
(36,а)
Cos 2р + ^-^^-^Со^2р + = 0
а1 !, < -аР
и значение нормального напряжения С п на площадке текучести
а п = + СО^2Р. (36,б)
6.2. Проверяется условие текучести гР ^Р
Я = Яп2р2р(яр -ап).
(37)
Если Я < 0, то Если Я)0, то
а [ =а [; а тъ =а 3. (38)
а 1 =(ап -р) + 42р{Яр -ап) + р2;
° Зг=(а п - р)-ЩЯЯ^ф?.
(39)
Рассмотренная методика реализована в программе для расчета напряженно-деформированного состояния систем «сооружение-основание» методом конечных элементов являющейся дальнейшим развитием программного комплекса"Сгаск" /1,2/.
Выдача расчетной информации
Разработанная математическая модель и программный комплекс позволяют получить следующую информацию:
-вертикальные V и горизонтальные и перемещения для всех узловых точек расчетной схемы конечных элементов;
-компоненты напряжений ах, ау, т (или а1, а3 и направления их действия а) для каждого конечного элемента;
-информацию о характере работы материала сооружения и грунта в пределах каждого элемента ( упругое или пластическое поведение);
-информацию о характере работы контактных элементов (их возможное раскрытие, сдвиг при т > тпред или упругое поведение);
- определить коэффициент устойчивости по любой поверхности, представленной в расчетной схеме явно с помощью контактных элементов с учетом реального характера распределения напряжений;
1/2008_М ВЕСТНИК
-имеется возможность вывода и другой интересующей информации.
Для повышения эффективности работы пользователя с программой дополнительные сервисные возможности комплекса реализуют:
a) диалоговый режим работы при создании и коррекции исходных данных и при работе всех программ комплекса;
b) хранение исходных данных и результатов расчетов в виде общей базы данных;
c) графическое представление расчетных схем и результатов расчетов;
d) возможность конвертации исходных данных и результатов расчетов в графический пакет AutoCAD для последующего редактирования и распечатки на принтере или плоттере.
Графические программы, предназначенные для визуализации расчетных схем и результатов расчетов позволяют просматривать в нескольких режимах визуализации исходную и результирующую информацию, сохранять ее в виде графических файлов и распечатывать на принтере.
Пример расчета
По представленному программному комплексу выполнен методический пример расчета упруго-пластической работы системы в виде упругого штампа (Еш =2,1 106 кг/см2) и песчаного грунта основания (Е = 850 кг/см2).
Размер штампа 194x194x70 см, размер моделируемой области грунта основания ВхН = 1180x600 см.
Расчетная область системы состоит из 1826 конечных элементов - из них 1806 изопараметрических элементов с квадратичной аппроксимацией перемещений моделировали грунт и штамп, и 20 контактных элементов моделировали контакт между штампом и основанием.
Расчет выполнены при поэтапной загрузке штампа вертикальной силой (Ршах= 64 т) и последующей его загрузке горизонтальной силой Н до величины, при которой происходит неустойчивое развитие пластических деформаций, приводящих к потере несущей способности системы вследствие сдвига.
Результаты расчетов представлены в виде: общей картины деформированного состояния системы при Р = 64 т и Н =20 т (рис. 6); эпюр горизонтальных перемещений по глубине основания (рис. 7, 9); схем площадок сдвиговых деформаций и возможной схемы образования глубинного сдвига в основании штампа (рис. 8, 10).
На рис. 6 показано, что при Р = 64 т, Н =20 т в основании под верховой гранью штампа произошли раскрытие контактных элементов на глубину 41 см и сдвиговые деформации по контакту (тупр>тпр) лишь на ограниченном участке подошвы штампа, 18 (см. рис. 6), что свидетельствует об отсутствии условий для плоского сдвига штампа по основанию. Об этом свидетельствует также анализ эпюры горизонтальных перемещений в системе (рис. 7). Из рисунка видно, что при увеличении горизонтальной силы с Н =16 тдо Н =20 т по мере увеличения числа итераций в расчетах происходит практически полная стабилизация величины перемещения штампа.
Для суждения о возможности смешенного сдвига штампа вместе с участком основания обратимся к рис. 8, где показана "возможная" кривая сдвига, которая начинается под штампом на участке упругой работы контакта и далее следующей по направлению наклона площадок сдвига в основании. Как видно из рисунка на заключительном участке "возможная" кривая сдвига проходит в зоне упругой работы материала основания, что, следовательно, исключает возможность проявления глубинного сдвига.
Рис.6. Картина деформированного состояния Р — 64 т, Н — 20 т (масштаб перемещений 1:20) 1 ■ зона раскрытия контактных элементов;--зона сдвига ( %упр > хкр)
Рис.7. Эпюры горизонтальных перемещений Р— 64 т, Н—20 т: Значения перемещений (см) указаны в следующем порядке--перемещения накопленные на предыдущем этапе нагружения
Н—16 т;
--перемещения при Н—20 т и числе итераций п —300;--то же при п —600
11::
I ^ % 111 #1,
и
Рис.8. Схема площадок сдвиговых деформаций в основании штампа Р— 64 т, Н—20 т:--грани-
^ пр
ца между областью пластических деформаций -—1,0 и областью упругих деформаций
т
^ пр
>1,0 /// - наклон площадок сдвига;
Рис.9. Эпюры горизонтальных перемещений Р — 64 т, Н — 24 т: Значения перемещений (см) указаны в следующем порядке: - перемещения на предыдущем этапе нагружения Н—20 т; - перемещения при Н—24 т и числе итераций п —500; - то же при п —1000
Дальнейшее увеличение сдвигающей нагрузки с Н —20 т до Н —24 т, как показано на рис. 10, приводит к потере несущей способности системы, так как кривые возможного сдвига
при данной нагрузке полностью располагаются в зоне пластических деформаций грунта основания.
I р=е1г
Рис.10. Схемы глубинного сдвига в основании штампа Р= 64 т, Н=24 т: кривые возможных сдвигов:_- расчетная Р= 64 т, Н=24 т.....- по СНиП 2.02.02-85 (Р= 64 т,
Н=24,87 т)
^ нр
--граница между областью пластических деформаций -=1,0 и областью упругих
т
деформаций —— >1,0 /// - наклон площадок сдвига
Указанный вывод о потере несущей способности системы может быть сделан также о их данных, приведенных на рис. 9, где показан незатухающий характер развития пластических деформаций при увеличении числа итераций в расчетах.
Выводы
1. Представленная математическая модель по сравнению с традиционными методами расчета дает качественно новую информацию (оценки предельного состояния различных швов и трещин при их различных прочностных и деформативных характеристиках; количественные оценки глубины и величины раскрытия различных швов, в том числе и контактного шва; размеры зон разуплотнения основания и областей пластических деформаций и их влияние на состояние контактного шва и т.д.). Это позволяет более объективно и детально оценивать напряженно-деформированное состояние сооружения и основания при выработке критериев безопасной работы сооружения.
Список литературы
1. Зерцалов М.Г., Толстиков В.В., Иванов В.А. «Реферат программного комплекса "Трещина" (Crack)» - Журнал «Механика грунтов, основания и фундаменты» №5, 1988.
2. Толстиков В.В. «Моделирование работы швов и трещин в расчетах напряженно-деформированного состояния бетонных плотин» - Научно-технический журнал «Вестник МГСУ» №2, 2006.
