Предисловие
1 | РАЗРАБОТАНЫ И ВНЕСЕНЫ | Акционерным обществом "Казахстанский дорожный научно-исследовательский институт" (АО "КаздорНИИ") |
2 | УТВЕРЖДЕНЫ И ВВЕДЕНЫ В ДЕЙСТВИЕ | Приказом Председателя Комитета автомобильных дорог Министерства по инвестициям и развитию Республики Казахстан № 122 от 21 декабря 2018 года |
3 | СОГЛАСОВАНЫ | Акционерным обществом "НК "ҚазАвтоЖол" № 03/14-2-2591-И от 12 ноября 2018 года |
4 | СРОК ПЕРВОЙ ПРОВЕРКИ | 2023 год |
ПЕРИОДИЧНОСТЬ ПРОВЕРКИ | 5 лет | |
5 | ВВЕДЕНЫ ВПЕРВЫЕ |
Документ доступен к просмотру в информационно-правовой системе нормативно-правовых актов Республики Казахстан "Әдiлет"
Настоящие рекомендации не могут быть полностью или частично воспроизведены, тиражированы и распространены без разрешения Комитета автомобильных дорог Министерства по инвестициям и развитию Республики Казахстан
Содержание
1 | |
2 | |
3 | |
4 | |
4.1 | |
4.2 | |
4.3 | |
4.4 | |
4.5 | Вычисление матрицы градиентов четырехугольного элемента с восемью узлами |
4.6 | |
4.7 | Краткая инструкция по использованию программы BASIC_NDS_MKE_8_uzlov |
4.8 | Примеры расчетов по программе BASIC_NDS_MKE_8_uzlov |
Приложение А (обязательное) Исходный код программы BASIC_NDS_MKE_8_uzlov Приложение Б (обязательное) Исходный код подпрограммы BtDBnds | |
Введение
Настоящие методические рекомендации разработаны согласно плану работ, принятому в рамках выполнения темы Б01.02 "Разработка теоретических решений и программ для определения напряженно-деформированного состояния многослойной конструкции дорожной одежды".
Документ содержит постановку типовой задачи по определению напряженно-деформированного состояния многослойной дорожной конструкции, статически нагруженной весом автомобиля, и основные формулы из алгоритма расчета методом конечных элементов. В методе применяется высокоточный квадратичный прямоугольный конечный элемент с восемью узлами. В приложении А приведен текст (листинг) расчетной программы BASIC_NDS_MKE_8_uzlovна языке MATLAB [1-4], снабженный необходимыми комментариями. В приложении Б приводится листинг подпрограммы BDB0, предназначенный для вычисления матрицы жесткости элемента. Приводится также краткая инструкция по использованию программы BASIC_NDS_MKE_8_uzlov. В конце документа приведен список использованной литературы.
1 Область применения
1.1 Настоящие рекомендации распространяются на сеть автомобильных дорог общего пользования Республики Казахстан и предназначены для решения вопросов, связанных с проектированием многослойных автомобильных дорог общего пользования.
1.2 Рекомендациями следует руководствоваться при проектировании конструкций дорожных одежд для автомобильных дорог общего пользования, для расчета дорожных одежд на стадиях проектирования и эксплуатации, а также при решении инженерно-экономических задач применительно к автомобильным дорогам.
2 Нормативные ссылки
Для применения настоящих рекомендаций необходимы следующие ссылочные нормативные документы:
СП РК 3.03-103-2014 "Проектирование жестких дорожных одежд"
СП РК 3.03-104-2014 "Проектирование нежестких дорожных одежд"
Примечание - При пользовании настоящими рекомендациями целесообразно проверить действие ссылочных документов по ежегодно издаваемому информационному указателю "Нормативные документы по стандартизации", составленному по состоянию на текущий год и соответствующим ежемесячно издаваемым информационным указателям, опубликованным в текущем году. Если ссылочный документ заменен (изменен), то при пользовании настоящими рекомендациями следует руководствоваться замененным (измененным) стандартом. Если ссылочный документ отменен без замены, то положение, в котором дана ссылка на него, применяется в части, не затрагивающей это ссылку.
3 Термины и определения
В настоящих рекомендациях применяются следующие термины с соответствующими определениями:
3.1 Дорожная одежда: Многослойная конструкция в пределах проезжей части автомобильной дороги, воспринимающая нагрузку от автотранспортного средства и передающая ее на грунт.
3.2 Земляное полотно: Конструктивный элемент, служащий основанием для размещения дорожной одежды, а также технических средств организации дорожного движения и обустройства автомобильной дороги.
3.3 Деформация: Относительная величина, характеризующая изменение линейных размеров тела по отношению к его первоначальным размерам.
3.4 Напряжение: Относительная величина, определяемая нормальной или касательной нагрузкой, приходящейся на единицу площади.
3.5 Дорожная одежда нежесткая: Дорожная одежда со слоями, устроенными из разного вида асфальтобетонов, из материалов и грунтов укрепленных битумом, цементом, известью, комплексными и другими вяжущими, а также из слабосвязных зернистых материалов
3.6 Конструктивный слой: Каждый слой дорожной одежды, состоящий из однородных материалов и отличающийся от соседних слоев видом материалов, его прочностью и составом. Учитывается при расчете прочности дорожной одежды.
3.7 Покрытие дорожное: Одно- или многослойная верхняя часть дорожной одежды, устраиваемая на дорожном основании, непосредственно воспринимающая нагрузки от транспортных средств и предназначенная для обеспечения заданных эксплуатационных требований и защиты дорожного основания от воздействия атмосферных факторов.
3.8 Основание дорожное: Нижний несущий слой дорожной одежды, воспринимающий нагрузки от транспортных средств совместно с покрытием и предназначенный для ее распределения на дополнительные слои или непосредственно на грунт земляного полотна.
3.9 Метод конечных элементов: Численный метод решения дифференциальных уравнений математической физики.
4 Методика использования программы BASIC_NDS_MKE_8_uzlov
4.1 Общая постановка задачи
Для построения математической модели автомобильной дороги общего пользования рассмотрим конструкцию, состоящую из двухслойной асфальтобетонной дорожной одежды и основания дорожной одежды, состоящего из щебеночно-песчаной смеси с 8 % цементной добавкой, щебеночной смеси и гравийно-песчаной смеси. Основание дорожной одежды устраивается на грунтовом основании из легкого суглинка (рисунок 1).
Возможности численного метода – метода конечных элементов, применяемого для решения задачи о равновесии многослойной конструкции, позволяют сравнительно легко назначать требуемые значения геометрических размеров конструктивных слоев и физико-механических характеристик материалов, используемых в конструкции.
Поэтому, приведенные в таблице 1 геометрические размеры конструктивных слоев в рассматриваемой математической модели (рисунок 1) и физико-механические свойства материалов конструктивных элементов (в том числе модули упругости) назначены согласно [6].
1 – покрытие, асфальтобетон мелкозернистый плотный; 2 – асфальтобетон крупнозернистый пористый; 3 – щебеночно-песчаная смесь 8 %; 4 – щебеночная смесь; 5 –гравийно-песчаная смесь; 6 – грунт суглинок легкий
Рисунок 1 – Схематический вид дорожной конструкции
Значения физико-механических свойств материалов конструктивных элементов вводятся в 108-110 строках приложения А.
Таблица 1 – Характеристики конструктивных слоев
Материал конструктивного слоя | Модуль упругости Е, МПа | Коэффициент Пуассона, | Высота слоя, м |
1 Асфальтобетон мелкозернистый плотный | 3200 | 0,18 | 0,05 |
2 Асфальтобетон крупнозернистый пористый | 2000 | 0,27 | 0,10 |
3 Щебеночно-песчаная смесь 8% | 800 | 0,30 | 0,40 |
4 Щебеночная смесь | 275 | 0,30 | 0,45 |
5 Гравийно-песчаная смесь | 180 | 0,30 | 1,00 |
6 Грунт суглинок легкий | 50 | 0,35 | - |
Для компьютерной реализации алгоритма решения задачи о напряженно-деформированном состоянии дорожной конструкции методом конечных элементов разработана расчетная программа BASIC_NDS_MKE_8_uzlovна алгоритмическом языке системы MATLAB.
Исследуемая область разбивается на квадратичные прямоугольные конечные элементы с восемью узлами (рисунок 2). Местная нумерация узлов с 1 по 8 в таком конечном элементе начинается с левого нижнего узла, и осуществляется в направлении против часовой стрелки.
Рисунок 2 – Прямоугольный элемент с восемью узлами
Для создания конечно-элементной сетки применялись переменные шаги и в направлении координатных осей и (рисунок 3).
В каждом горизонтальном ряде (рисунок 3) расположено по 60 элементов (egor=60), а количество горизонтальных рядов элементов – 23 (ever=23). Следовательно, исследуемая область разбита на 1380 конечных элементов, а число узлов равно np=(3*ever+2)*egor+2*ever+1=4307 (строки 12-15 приложения А).
Рисунок 3 – Схематический вид конечно-элементной разбивки исследуемой области
Глобальная нумерация узлов начинается с крайнего левого вертикального ряда узлов, и осуществляется сверху вниз и перемещается слева направо (рисунок 3).
При этом разбивка получается нерегулярной, т.к. в восьмиузловых элементах центральный узел отсутствует. В данном случае, в одном вертикальном ряде имеются 47 узлов, а в следующем – 24 узла.
Требуется определить компоненты вектора перемещений в узловых точках конечно-элементной сетки под действием веса автомобиля, передаваемого через одно колесо и приложенного в виде вертикальных нагрузок в соответствующих узлах (рисунок 4).
L и H – ширина и высота исследуемой области соответственно; P – нагрузка от колеса транспортного средства; ui, vi – компоненты перемещения точек на границах в направлении координатных осей x и y; I – VI конструктивные слои дорожной одежды и грунтового основания
Рисунок 4 – Расчетная схема задачи
4.2 Формирование массивов координат
В программе BASIC_NDS_MKE_8_uzlovзаданы массивы dx и dy, содержащие значения шагов по координатам. Наличие таких массивов позволяет задавать переменные шаги по координатам (строки 44-63 приложения А).
В блоке формирования координат программы (строки 68-91 приложения А) глобальные номера узлов текущего конечного элемента обозначены через n1, n2, …, n8, соответствующие локальным номерам 1,2, …, 8, приведенным на рисунке 2. Конкретные значения глобальных номеров узлов, вычисляемые по алгоритмам вычисления n1, n2, …, n8 в программе BASIC_NDS_MKE_8_uzlovопределяются номером элемента в текущем горизонтальном ряде узлов (m=1:egor) и номером элемента в текущем вертикальном ряде узлов (n=1:ever).
4.3 Граничные условия
Из рисунка 4 видно, что граничные условия в перемещениях по бокам и на нижней границе исследуемой области задаются в виде нулевых значений соответствующих компонентов вектора перемещений (u=0 и u=v=0). Чтобы задавать граничные условия в виде известных перемещений некоторых узлов также нужно создать специальный массив Mz, содержащий номера узлов задаваемых известных перемещений (строки 23-43 приложения А).
Граничные условия в нагрузках задаются в виде сосредоточенных вертикальных усилий, рассчитанных от веса автомобиля, приходящегося на одно колесо. Методика расчета граничных условий в нагрузках приведена в приложении А (строки 112-166). Конечный результат расчета представляет собой вектор правой части уравнения равновесия системы (1).
4.4 Определяющие соотношения уравнения равновесия
Уравнение равновесия системы в случае применения квадратичных элементов, как и в случае линейных треугольных конечных элементов, в матричной форме имеет вид [7]:
(1)
где [K] – матрица жесткости системы размерности (2*np,2*np);
– вектор узловых перемещений размерности 2*np;
– вектор узловых нагрузок размерности 2*np;
np – общее число узлов в конечно-элементной сетке.
Порядок вычисления значений компонентов вектора , соответствующие граничным условиям в нагрузках приведены в строках 150-165приложения А. Остальные компоненты этого вектора, в силу условий равновесия узлов, равны нулю.
В уравнении (1) матрица жесткости системы определяется в виде суммы
где – общее число элементов в конечно-элементной сетке; – матрица жесткости текущего элемента "".
Матрица жесткости элемента вычисляется с помощью объемного интеграла (2):
(2)
где алгоритм вычисления матрицы упругости приведены в строках 204-207 приложения А.
В выражении (2) так называемая матрица градиентов формируется путем дифференцирования функций формы по координатам .
В случае нелинейной зависимости функций формы от координат , что имеет место в случае квадратичных элементов с восемью узлами, производные от функций формы по координатам сами являются функциями координат , и вычисление интеграла представляет собой самостоятельную задачу с громоздкими математическими выкладками.
4.5 Вычисление матрицы градиентов четырехугольного элемента с восемью узлами
Для численного интегрирования объемного интеграла в соотношении (2) предлагается [7] переходить к локальной системе координат (рисунок 5). Здесь замену переменных интегрирования можно сделать с помощью соотношения (3):
(3)
где – единичная толщина элемента, а – определитель матрицы преобразования координат Якоби. Тогда объемный интеграл в (2) приводится к виду (4):
(4)
Переход к переменным интегрирования и упрощает пределы интегрирования, что позволяет разработать единый алгоритм вычисления интегралов через квадратуру Гаусса-Лежандра.
Рисунок 5 – Система координат для линейного прямоугольного элемента
4.6 Алгоритм определения деформации и напряжений
Выражение перепишем в следующем виде (5):
(5)
где
,
, ,
,
, (6)
,
, .
Теперь из выражения (5) компоненту тензора деформации определим следующим образом:
(7)
или учитывая равенства (5)
(8)
В выражении (8) координаты и относятся к линейному прямоугольнику (рисунок 5), и им при рассмотрении восьмиугольных прямоугольников будут соответствовать и .
Аналогично определяются и остальные компоненты тензора деформации. Из формулы (5) для компоненты имеем:
или, учитывая (4) и , получим:
(9)
Деформация сдвига определяется аналогичным образом.
Для определения компонентов напряжений воспользуемся формулами вычисления напряжений [7]:
, и , (10)
где, для плоской деформации ;
- модуль упругости;
- коэффициент Пуассона.
Программа численной реализации приведенного алгоритма вычисления компонентов вектора деформации и вектора напряжений приведена в строках 412-527 приложения А настоящего документа. Кроме этого, в приложении А предусмотрены фрагменты программ построения графиков горизонтальных перемещений точек характерных вертикальных сечений (строки 297-312 приложения А), вертикальных перемещений точек горизонтальных сечений, расположенных на различных глубинах (строки 316-408 приложения А) и значений компонент векторов деформаций и напряжений (строки 441-475, 529-602 приложения А).
4.7 Краткая инструкция по использованию программы
BASIC_NDS_MKE_8_uzlov
Целью создания настоящего нормативного документа является разработка расчетной программы на языке MATLAB.
Расчетная программа предназначена для вычисления значений компонентов перемещений, деформаций и напряжений в точках исследуемой дорожной конструкции (рисунок 4).
В мировой практике существуют и другие пакеты программ, способных решать задачу о напряженно-деформированном состоянии в твердом деформируемом теле. Наиболее известным среди них является программный комплекс ANSYS. Однако, все они являются коммерческими, алгоритм их реализации скрытый, и простому пользователю использовать их для расчетов своих задач, тем более их видоизменить просто невозможно. Настоящий программный продукт BASIC_NDS_MKE_8_uzlovявляется открытым, алгоритм его реализации расписан открытым текстом, программа легко адаптируется и всегда есть возможность ее усовершенствовать. Например, ее можно дополнить отдельными блоками и решить задачи об условиях образования низкотемпературных трещин в дорожном покрытии, усталостных трещин в конструктивных слоях дорожной одежды.
Для запуска программы BASIC_NDS_MKE_8_uzlovдля решения конкретной задачи,в первую очередь необходимо, чтобы в компьютере была установлена одна из версии программного комплекса MATLAB,и в директорию этого комплекса был внесен текст (листинг) программы BASIC_NDS_MKE_8_uzlov. На рабочий стол компьютера выносится ярлык программного комплекса MATLAB (рисунок 6).
Рисунок 6 – Ярлык программного комплекса MATLAB
В приложении А настоящего документа приведен исходный код программы BASIC_NDS_MKE_8_uzlov. Для облегчения обращения к тексту листинг снабжен построчной нумерацией, которую нужно удалить перед установкой программы в директорию программного комплекса MATLAB.
Вызов программного комплекса MATLAB для расчетов осуществляется двойным щелчком на ярлыке (рисунок 6).
На рабочем столе появится главная страница MATLAB (рисунок 7).
Рисунок 7 – Вид главной страницы программы наMATLAB
На рисунке 7 на главной странице MATLAB приведен фрагмент программы BASIC_NDS_MKE_8_uzlov. Если же этот фрагмент будет отсутствовать, а на главной странице присутствует фрагмент другой программы, то поиск необходимой программы BASIC_NDS_MKE_8_uzlov будет произведен простым нажатием на кнопку OPEN, расположенной в верхней части главной страницы. После нажатия кнопки OPEN получим список программ, находящихся в директории MATLAB (рисунок 8).
Рисунок 8 – Вид списка программ в директории
Простым нажатием на название нужной программы BASIC_NDS_MKE_8_uzlovиз приведенного списка вызывается эта программа на главную страницу (рисунок 7).
Рисунок 9 – Строка 109-110 для введения физико-механических характеристик материалов слоев
По умолчанию в программе установлены физико-механические характеристики автомобильной дороги.Для изменения их значений данные необходимо ввести в строки 109-110 (рисунок 9). После введения физико-механических характеристик вводим значение прилагаемой вертикальной силы на поверхность в строку 150,а так же вводим номера поверхностных узлов куда необходимо приложить нагрузку (рисунок10).
Рисунок 10 – Строка для введения значения прилагаемой силы на поверхность и номеров узлов приложения сил
Следующим этапом при решении задачи является получение графических результатов по необходимым для пользователя требованиям. Ниже приведены примеры запроса графиков вертикальных перемещений узлов на разных глубинах (рисунки 11-17).
Рисунок 11 – Получение графиков вертикальных перемещений узлов на глубине h= 5 см
Рисунок 12 – Получение графиков вертикальных перемещений узлов на глубине h= 15 см
Рисунок 13 – Получение графиков вертикальных перемещений узлов на глубине h= 25 см
Рисунок 14 – Получение графиков вертикальных перемещений узлов на глубине h= 45 см
Рисунок 15 – Получение графиков компонент напряжений SigmaX
Рисунок 16 – Получение графиков компонент напряжений SigmaY
Рисунок 17 – Получение графиков компонент напряжений TauXY
Просмотрев вес текст программы на главной странице, и установив все необходимые исходные данные, а при необходимости, вносив изменения в тело программы можно будет запустить программу на выполнение. Для этого достаточно нажать на кнопкуRun (рисунок 18).
Рисунок 18 – Кнопка Run на панели инструментов для запуска программы
В тексте программы присутствуют команды, согласно которым на экран будет выводиться необходимая информация о результатах решения задачи. Такая информация может иметь графическую форму, может быть в виде таблиц, или в виде текстового файла. Все это зависит от воли заказчика.
Ниже, на рисунках 13-18, в качестве примера приведены результаты решения тестовой задачи, где приведены картины формирования поля горизонтальных и вертикальных перемещений точек характерных сечений исследуемой области и компоненты напряжений в характерных горизонтальных сечениях. В тестовой задаче дорожная конструкция нагружена распределенной вертикальной нагрузкой от одиночного колеса автомобиля, действующей в центре поверхности дорожного покрытия.
4.8 Примеры расчетов по программе BASIC_NDS_MKE_8_uzlov
Ниже приводятся результаты компьютерной реализации программы BASIC_NDS_MKE_8_uzlovв виде графиков перемещений и напряжений.
1–горизонтальные смещения точек вертикального сечения, отстоящего от вертикальной оси симметрии слева на 30 см; 2–горизонтальные смещения точек вертикальной оси симметрии; 3–горизонтальные смещения точек вертикального сечения, отстоящего от вертикальной оси симметрии справа на 30 см
Рисунок 19 – Горизонтальные смещения точек вертикального сечения
Рисунок 20 – Вертикальные смещения точек горизонтальных сечений, расположенных на глубинах: 1– h=5 см; 2–h=15 см; 3–h=25 см; 4–h=45 см; 5–h=80 см.
Рисунок 21 – Распределение напряжений на глубинах h
Рисунок 22 – Распределение напряжений на глубинах h
Рисунок 23 – Распределение напряжений на глубинах h
Приложение А
(обязательное)
Исходный код программы BASIC_NDS_MKE_8_uzlov
1 | % БАЗОВАЯ ВЕРСИЯ ПРОГРАММЫ BASIC_NDS_MKE_8_uzlov |
2 | % РЕШЕНИЕ ТЕСТОВОЙ ЗАДАЧИ О НДС С ПОМОЩЬЮ |
3 | %ЧЕТЫРЕХУГОЛЬНОГО |
4 | % КВАДРАТИЧНОГО ЭЛЕМЕНТА С 8-Ю УЗЛАМИ МЕТОДОМ КОНЕЧНЫХ |
5 | % ЭЛЕМЕНТОВ |
6 | |
7 | % Начало разработки 8.07.2017 г |
8 | % Окончание 1.09.2018 г |
9 | |
10 | Tic |
11 | |
12 | ever=23; egor=60; % количество конечных элементов в вертикальном (ever) и |
13 | % горизонтальном (egor) рядах элементов в конечно-элементной |
14 | % разбивке |
15 | np=(3*ever+2)*egor+2*ever+1; % общее число узлов |
16 | x=zeros(np); y=zeros(np); |
17 | |
18 |
% ПОРЯДОК НУМЕРАЦИИ ЭЛЕМЕНТОВ - СЛЕВА НАПРАВО (РЯДЫ |
19 |
% ПОРЯДОК НУМЕРАЦИИ УЗЛОВ - СЛЕВА НАПРАВО. РЯДЫ УЗЛОВ ТАКЖЕ |
20 |
% СНИЗУ. ПРОМЕЖУТОЧНЫЕ РЯДЫ УЗЛОВ, ПРОХОДЯЩИЕ ЧЕРЕЗ СЕРЕДИНУ |
21 | % СОДЕРЖАТ egor+1 УЗЛОВ |
22 | |
23 | % ФОРМИРОВАНИЕ МАССИВА НОМЕРОВ УЗЛОВ С ЗАДАННЫМИ ЗНАЧЕНИЯМИ |
24 | % ПЕРЕМЕЩЕНИЙ ДЛЯ ЗАДАНИЯ ГРАНИЧНЫХ УСЛОВИЙ В ПЕРЕМЕЩЕНИЯХ |
25 | |
26 | Mz=zeros(400); |
27 | np=ever*(3*egor+2)+2*egor+1; |
28 | for i=1:(2*egor+1) |
29 | Mz(i)=i; |
30 | End |
31 | k=2*egor+2; |
32 | for i=1:ever |
33 | Mz(k)=(2*egor+1)+1+(i-1)*(3*egor+2); |
34 | Mz(k+1)=Mz(k)+egor; |
35 | Mz(k+2)=Mz(k+1)+1; |
36 | Mz(k+3)=Mz(k+2)+2*egor; |
37 | k=k+4; |
38 | End |
39 | for i=1:(2*egor+1) |
40 | Mz(k)=Mz(i)+np; |
41 | k=k+1; |
42 | End |
43 | |
44 | % МАССИВЫ ШАГОВ ПО КООРДИНАТНЫМ ОСЯМ В МЕТРАХ |
45 | format long |
46 | dx=zeros(2*egor); |
47 | dy=zeros(2*ever); |
48 | |
49 | dx=[ |
50 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 ... |
51 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 ... |
52 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 ... |
53 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 ... |
54 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 ... |
55 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 ... |
56 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 ... |
57 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 0.250e-01 ... |
58 | 0.250e-01 0.250e-01 0.250e-01 0.250e-01]; |
59 | |
60 | dy=[0.1000 0.1000 0.1000 0.1000 0.1000 0.1000 0.1000 ... |
61 | 0.1000 0.1000 0.1000 0.1000 0.1000 0.0750 0.0750 ... |
62 | 0.0625 0.0625 0.0500 0.0500 0.0500 0.0500 0.0250 ... |
63 | 0.0250 0.0250]; |
64 | |
65 | yy=[0.000 0.200 0.400 0.600 0.800 1.000 1.200 1.400 1.600 1.800 2.000 2.200 ... |
66 | 2.400 2.550 2.700 2.825 2.950 3.050 3.150 3.250 3.350 3.400 3.450 3.500]; |
67 | |
68 | % БЛОК ФОРМИРОВАНИЯ КООРДИНАТ |
69 | % Положительное направление оси х - слева направо |
70 | % Положительное направление оси у - снизу вверх |
71 | k=1; |
72 | NN=zeros(ever); |
73 | for n=1:ever |
74 | for m=1:egor |
75 | n1=2*m-1+(3*egor+2)*(n-1); n2=n1+1; n3=n1+2; |
76 | n4=n3+2*egor-m+1; n5=n3+(3*egor+2); |
77 | n6=n5-1; n7=n5-2; n8=n4-1; |
78 | |
79 | x(n1)=(m-1)*2*dx(m); x(n2)=x(n1)+dx(m); x(n3)=x(n1)+2*dx(m); |
80 | x(n4)=x(n3); x(n5)=x(n3); x(n6)=x(n2); x(n7)=x(n1); x(n8)=x(n1); |
81 | |
82 | y(n1)=yy(n); y(n2)=y(n1); y(n3)=y(n1); y(n4)=y(n3)+dy(n); |
83 | y(n5)=y(n3)+2*dy(n); y(n6)=y(n5); y(n7)=y(n5); y(n8)=y(n4); |
84 | if(m==36) |
85 | NN(2*k-1)=n3; |
86 | NN(2*k)=n4; |
87 | k=k+1; |
88 | end |
89 | End |
90 | End |
91 | NN(1:2*ever+1); |
92 | |
93 | % ИСХОДНЫЕ ДАННЫЕ И ОБНУЛЕНИЕ МАТРИЦ |
94 | np2=2*np; |
95 | K=zeros(np2,np2); |
96 | U=zeros(np2+1); |
97 | U0=zeros(np2+1); |
98 | F=zeros(np2); |
99 | |
100 | % НАИМЕНОВАНИЕ МАТЕРИАЛОВ КОНСТРУКТИВНЫХ СЛОЕВ (СВЕРХУ ВНИЗ) |
101 | % 1. АБ МЗ ПЛОТНЫЙ |
102 | % 2. АБ КЗ ПОРИСТЫЙ |
103 | % 3. ЩЕБЕНОЧНО-ПЕСЧАНАЯ СМЕСЬ 8% |
104 | % 4. ЩЕБЕНОЧНАЯ СМЕСЬ |
105 | % 5. ГРАВИЙНО-ПЕСЧАНАЯ СМЕСЬ |
106 | % 6. ГРУНТ СУГЛИНОК ЛЕГКИЙ |
107 | |
108 |
% ФИЗИКО-МЕХАНИЧЕСКИЕ ХАРАКТЕРИСТИКИ МАТЕРИАЛОВ СЛОЕВ(СВЕРХУ |
109 | E=[6933e+05 4142e+05 1000e+05 250e+05 180e+05 58e+05]; % Модуль упр.кг/м^2 |
110 | nu=[0.18 0.27 0.30 0.30 0.30 0.35]; % КОЭФФИЦИЕНТЫ ПУАССОНА |
111 | |
112 | % ГРАНИЧНЫЕ УСЛОВИЯ В НАГРУЗКАХ |
113 | |
114 | % Нагрузка от одиночного колеса |
115 | |
116 | % ДАННЫЕ ИЗ ИНТЕРНЕТА |
117 | % 1 Па=0.102 кгс/м^2 |
118 | % 0.6 МПа=6.12 кгс/см^2=60000.0 кгс/м^2 |
119 | |
120 |
% РАСПРЕДЕЛЕННАЯ ПО УГЛОВЫМ УЗЛАМ ЭЛЕМЕНТА НАГРУЗКА ОТ |
121 |
% ПРИКЛАДЫВАЕТСЯ ПО ЦЕНТРУ РАСЧЕТНОЙ СХЕМЫ. При этом промежуточные |
122 | % игнорируются |
123 | |
124 |
% НАГРУЗКА НА КОЛЕСО РАССЧИТЫВАЕТСЯ ИЗ РАСЧЕТА НОРМАТИВНОЙ |
125 | % НАГРУЗКИ q=0.6 MPa |
126 | % ЕСЛИ 1 Па=0.102 кгс/м^2 (из Интернета) ТО |
127 | % 0.6 МПа=6.12 кгс/см^2 = 60000.0 кгс/м^2. |
128 | % Таким образом, если на 1 погонный метр полотна приходится 60000 кг веса, то |
129 | % на одно расстояние 0.05 м – 3000 кг веса. |
130 | % Если ширину колеса принять равной 40 см, а ширину одного конечного |
131 | % элемента в горизонтальном направлении – 5 см, то транспортная нагрузка |
132 | % будет распределена по 9 узлам в пределах 8 конечных элементов. |
133 | % Расстояние между угловыми узлами - 0.05 м (5 cм - ширина конечного элемента). |
134 | % При этом надо учесть, что |
135 | % два крайних узла (1-ый и 9-ый узлы) будут нагружены наполовину. |
136 | % В таком случае, нагрузки в узловых точках будут распределены |
137 | % следующим образом: |
138 | % P1=P9=1500 кг, P2=P3= … =P8=3000 кг |
139 | |
140 | % РАСЧЕТ НАГРУЗКИ НА УЗЕЛ |
141 | % РАССТОЯНИЕ МЕЖДУ УЗЛАМИ 0.025 М (учитываются и промежуточные узлы) |
142 | % НА 1 М ПРИХОДИТСЯ 60000 КГ ВЕСА, А НА ОДНО РАССТОЯНИЕ 0.025 М |
143 | % ПРИХОДИТСЯ 1500 КГ ВЕСА |
144 | |
145 | % РЕЗУЛЬТАТЫ БУДУТ ПОЛУЧЕНЫ В КГ И МЕТРАХ |
146 | |
147 | % ВЕЛИЧИНЫ ВЕРТИКАЛЬНЫХ СИЛ, ПРИЛОЖЕННЫХ ПО УГЛОВЫМ УЗЛАМ |
148 | % В ЦЕНТРАЛЬНОЙ ЧАСТИ ЭПЮРЫ |
149 | |
150 | P=0.3e+04; % кг |
151 | |
152 | % НОМЕРА УГЛОВЫХ УЗЛОВ ПРИЛОЖЕНИЯ СИЛ |
153 | |
154 | nF(1:9)=[4239 4241 4243 4245 4247 4249 4251 4253 4255]; |
155 | for i=1:9 |
156 | j=nF(i); |
157 | xP(i)=x(j); |
158 | end |
159 | |
160 | % ПРИЛОЖЕНИЕ СИЛ К УГЛОВЫМ УЗЛАМ |
161 | F(nF(1)+np)=-P/2; |
162 | F(nF(9)+np)=-P/2; |
163 | for i=1:7 |
164 | F(nF(1+i)+np)=-P; |
165 | end |
166 | |
167 | % ФОРМИРОВАНИЕ МАТРИЦЫ ЖЕСТКОСТИ ЭЛЕМЕНТА [ke] И |
168 | % ГЛОБАЛЬНОЙ МАТРИЦЫ ЖЕСТКОСТИ СИСТЕМЫ [K] |
169 | |
170 | mpe=[-0.577350 0.577350]; % КООРДИНАТЫ ТОЧЕК ИНТЕГРИРОВАНИЯ |
171 | mH=[1 1]; % ВЕСОВЫЕ КОЭФФИЦИЕНТЫ |
172 | K=zeros(np2,np2); |
173 | for n=1:ever |
174 | for m=1:egor |
175 | n1=2*m-1+(3*egor+2)*(n-1); n2=n1+1; n3=n1+2; |
176 | n4=n3+2*egor-m+1; n5=n3+(3*egor+2); |
177 | n6=n5-1; n7=n5-2; n8=n4-1; |
178 | |
179 | X1=x(n1); Y1=y(n1); |
180 | X2=x(n3); Y2=y(n3); |
181 | X3=x(n5); Y3=y(n5); |
182 | X4=x(n7); Y4=y(n7); |
183 | |
184 |
% ФОРМИРОВАНИЕ МАТРИЦЫ УПРУГОСТИ ДЛЯ СЛУЧАЯ ПЛОСКОЙ |
185 | if((n>=1)&(n<=7)) |
186 | e1=E(6); nu1=nu(6); % ГРУНТ СУГЛИНОК ЛЕГКИЙ |
187 | end |
188 | if((n>=8)&(n<=12)) |
189 | e1=E(5); nu1=nu(5); % ГРАВИЙНО-ПЕСЧАНАЯ СМЕСЬ |
190 | end |
191 | if((n>=13)&(n<=16)) |
192 | e1=E(4); nu1=nu(4); % ЩЕБЕНОЧНАЯ СМЕСЬ |
193 | end |
194 | if((n>=17)&(n<=20)) |
195 | e1=E(3); nu1=nu(3); % ЩЕБЕНОЧНО-ПЕСЧАНАЯ СМЕСЬ 8% |
196 | end |
197 | if((n>=21)&(n<=22)) |
198 | e1=E(2); nu1=nu(2); % АБ КЗ ПОРИСТЫЙ |
199 | end |
200 | if(n==23) |
201 | e1=E(1); nu1=nu(1); % АБ МЗ ПЛОТНЫЙ |
202 | end |
203 | |
204 | D=zeros(3,3); |
205 | D(3,3)=e1/(2*(1+nu1)); D(2,2)=2*D(3,3)*(1-nu1)/(1-2*nu1); |
206 | D(1,1)=D(2,2); D(1,2)=2*D(3,3)*nu1/(1-2*nu1); |
207 | D(2,1)=D(1,2); |
208 | |
209 | dJ=abs((X1-X2)*(Y2-Y4))/4; % якобиан |
210 | |
211 |
% ИСПОЛЬЗОВАНИЕ КВАДРАТУРЫ ГАУССА ДЛЯ ЧИСЛЕННОГО |
212 | % ДВОЙНОГО ИНТЕГРАЛА В МАТРИЦЕ ЖЕСТКОСТИ ЭЛЕМЕНТА [Ke] |
213 | |
214 | % BtDBnds1(psi,eta,D,X1,X2,Y1,Y3) - подпрограмма вычисления |
215 | % матрицы жесткости элемента |
216 | Ke=zeros(16,16); |
217 | psi=mpe(1); eta=mpe(1); |
218 | k11=BtDBnds1(psi,eta,D,X1,X2,Y1,Y3); |
219 | psi=mpe(1); eta=mpe(2); |
220 | k12=BtDBnds1(psi,eta,D,X1,X2,Y1,Y3); |
221 | psi=mpe(2); eta=mpe(1); |
222 | k21=BtDBnds1(psi,eta,D,X1,X2,Y1,Y3); |
223 | psi=mpe(2); eta=mpe(2); |
224 | k22=BtDBnds1(psi,eta,D,X1,X2,Y1,Y3); |
225 | Ke=dJ*(mH(1)*(mH(1)*k11+mH(2)*k12)+mH(2)*(mH(1)*k21+mH(2)*k22)); |
226 | |
227 | ne8=[n1 n2 n3 n4 n5 n6 n7 n8]; |
228 | for i=1:8 |
229 | for j=1:8 |
230 | K(ne8(i),ne8(j))=K(ne8(i),ne8(j))+Ke(i,j); |
231 | K(ne8(i),ne8(j)+np)=K(ne8(i),ne8(j)+np)+Ke(i,j+8); |
232 | K(ne8(i)+np,ne8(j))=K(ne8(i)+np,ne8(j))+Ke(i+8,j); |
233 | K(ne8(i)+np,ne8(j)+np)=K(ne8(i)+np,ne8(j)+np)+Ke(i+8,j+8); |
234 | end |
235 | end |
236 | end |
237 | end |
238 | |
239 | disp('Матрица К сформирована') |
240 | |
241 | % ПРЕОБРАЗОВАНИЕ СЛАУ |
242 | |
243 | np2=2*np; |
244 | i=1; |
245 | for n=1:np2 |
246 | if(n==Mz(i)) |
247 | F(n)=K(n,n)*U(n); |
248 | for m=1:np2 |
249 | if(m~=n) |
250 | K(n,m)=0.0; |
251 | end |
252 | end |
253 | i=i+1; |
254 | end |
255 | end |
256 | |
257 | i=1; |
258 | for n=1:np2 |
259 | if(n==Mz(i)) |
260 | for m=1:np2 |
261 | if(m~=n) |
262 | F(m)=F(m)-K(m,n)*U(n); |
263 | K(m,n)=0.0; |
264 | end |
265 | end |
266 | i=i+1; |
267 | end |
268 | end |
269 | |
270 | % РЕШЕНИЕ СЛАУ |
271 | U=K\F; |
272 | disp('Система [K]{U}={F} решена') |
273 | |
274 | % Номера узлов характерных вертикальных сечений |
275 | nom12=[ |
276 | 49 146 231 328 413 510 595 692 ... |
277 | 777 874 959 1056 1141 1238 1323 1420 ... |
278 | 1505 1602 1687 1784 1869 1966 2051 2148 ... |
279 | 2233 2330 2415 2512 2597 2694 2779 2876 ... |
280 | 2961 3058 3143 3240 3325 3422 3507 3604 ... |
281 | 3689 3786 3871 3968 4053 4150 4235]; |
282 | nom15=[ |
283 | 61 152 243 334 425 516 607 698 ... |
284 | 789 880 971 1062 1153 1244 1335 1426 ... |
285 | 1517 1608 1699 1790 1881 1972 2063 2154 ... |
286 | 2245 2336 2427 2518 2609 2700 2791 2882 ... |
287 | 2973 3064 3155 3246 3337 3428 3519 3610 ... |
288 | 3701 3792 3883 3974 4065 4156 4247]; |
289 | nom18=[ |
290 | 73 158 255 340 437 522 619 704 ... |
291 | 801 886 983 1068 1165 1250 1347 1432 ... |
292 | 1529 1614 1711 1796 1893 1978 2075 2160 ... |
293 | 2257 2342 2439 2524 2621 2706 2803 2888 ... |
294 | 2985 3070 3167 3252 3349 3434 3531 3616 ... |
295 | 3713 3798 3895 3980 4077 4162 4259]; |
296 | |
297 | % ГРАФИКИ ГОРИЗОНТАЛЬНЫХ ПЕРЕМЕЩЕНИЙ УЗЛОВ В СРЕДНЕМ СЕЧЕНИИ |
298 | Y=zeros(2*ever+1); |
299 | X150=zeros(2*ever+1); |
300 | X120=zeros(2*ever+1); |
301 | X180=zeros(2*ever+1); |
302 | for i=1:2*ever+1 |
303 | j120=nom12(i); |
304 | j150=nom15(i); |
305 | j180=nom18(i); |
306 | Y(i)=y(j150); |
307 | X120(i)=U(j120)-0.0001; % СЛЕВА ОТ ВЕРТИКАЛЬНОЙ ОСИ СИММЕТРИИ НА 30 СМ |
308 | X150(i)=U(j150)*10; % НА ВЕРТИКАЛЬНОЙ ОСИ СИММЕТРИИ |
309 |
X180(i)=U(j180)+0.0001; % СПРАВА ОТ ВЕРТИКАЛЬНОЙ ОСИ СИММЕТРИИ НА |
310 | end |
311 | hPlot=plot(X150,Y,'b-',X120,Y,'r-',X180,Y,'k-.');grid on;xlabel('Ux, m');ylabel('y, m') |
312 | set(hPlot,'LineWidth',2); |
313 | |
314 | figure |
315 | |
316 |
% ФОРМИРОВАНИЕ МАССИВА НОМЕРОВ ГОРИЗОНТАЛЬНЫХ РЯДОВ УЗЛОВ |
317 | |
318 |
% ОБЩЕЕ ЧИСЛО УЗЛОВ В КАЖДОМ ПОЛНОМ ГОРИЗОНТАЛЬНОМ РЯДУ РАВНО |
319 | % А ЧИСЛО ТАКИХ ГОРИЗОНТАЛЬНЫХ РЯДОВ РАВНО ever+1 |
320 |
% СЕЧЕНИЯ, ПРОХОДЯЩИЕ ЧЕРЕЗ СЕРЕДИНУ ЭЛЕМЕНТОВ, НЕ |
321 |
% ОНИ ИСПОЛЬЗУЕТСЯ ДЛЯ ПОСТРОЕНИЯ ГРАФИКОВ ИСКОМЫХ ВЕЛИЧИН НА |
322 | % ГОРИЗОНТАЛЬНОМ УРОВНЕ (ГЛУБИНЕ) |
323 | |
324 | Nv=zeros(ever+1,2*egor+1); |
325 | for n=1:(ever+1) |
326 | for m=1:(2*egor+1) |
327 | I=2*n-1; |
328 | Nv(I,m)=m+(3*egor+2)*(n-1); |
329 | end |
330 | end |
331 | for n=1:ever |
332 | for m=1:(egor+1) |
333 | I=2*n; |
334 | Nv(I,2*m-1)=m+(2*egor+1)*n+(egor+1)*(n-1); |
335 | end |
336 | end |
337 | |
338 | format short |
339 | Set=zeros(2*ever+1,4); |
340 | for n=1:2*ever+1 |
341 | i=Nv(n,1); |
342 | Set(n,1)=n; |
343 | Set(n,2)=i; |
344 | Set(n,3)=y(i); |
345 | Set(n,4)=3.50-y(i); |
346 | end |
347 | for m=1:47 |
348 | % fprintf('%6u %6u %10.4f %10.4f\n',Set(m,1),Set(m,2),Set(m,3),Set(m,4)); |
349 | end |
350 | |
351 | % ГРАФИКИ ВЕРТИКАЛЬНЫХ ПЕРЕМЕЩЕНИЙ УЗЛОВ |
352 | % ГООРИЗОНТАЛЬНЫХ СЕЧЕНИЙ НА ГЛУБИНАХ h |
353 | X=zeros(2*egor+1); |
354 | Y=zeros(2*egor+1); |
355 | for i=1:(2*egor+1) |
356 | j=Nv(45,i); % ГЛУБИНА h=5 СМ |
357 | X(i)=x(j); |
358 | Y(i)=U(j+np); |
359 | end |
360 | hPlot=plot(X,Y,'-');grid on;xlabel('X, m');ylabel('Uy, m') |
361 | set(hPlot,'LineWidth',2); |
362 | hold on |
363 | |
364 | X=zeros(2*egor+1); |
365 | Y=zeros(2*egor+1); |
366 | for i=1:(2*egor+1) |
367 | j=Nv(41,i); % ГЛУБИНА h=15 СМ |
368 | X(i)=x(j); |
369 | Y(i)=U(j+np); |
370 | end |
371 | hPlot=plot(X,Y,'-');grid on;xlabel('X, m');ylabel('Uy, m') |
372 | set(hPlot,'LineWidth',2); |
373 | hold on |
374 | |
375 | X=zeros(2*egor+1); |
376 | Y=zeros(2*egor+1); |
377 | for i=1:(2*egor+1) |
378 | j=Nv(39,i); % ГЛУБИНА h=25 СМ |
379 | X(i)=x(j); |
380 | Y(i)=U(j+np); |
381 | Ux325(i)=U(j); |
382 | Uy325(i)=Y(i); |
383 | end |
384 | hPlot=plot(X,Y,'--');grid on;xlabel('X, m');ylabel('Uy, m') |
385 | set(hPlot,'LineWidth',2); |
386 | hold on |
387 | |
388 | for i=1:(2*egor+1) |
389 | j=Nv(35,i); % ГЛУБИНА h=45 СМ |
390 | X(i)=x(j); |
391 | Y(i)=U(j+np); |
392 | X345(i)=X(i); |
393 | Ux305(i)=U(j); |
394 | Uy305(i)=Y(i); |
395 | end |
396 | hPlot=plot(X,Y,'--');grid on;xlabel('X, m');ylabel('Uy, m') |
397 | set(hPlot,'LineWidth',2); |
398 | hold on |
399 | |
400 | for i=1:(2*egor+1) |
401 | j=Nv(29,i); % ГЛУБИНА h=80 СМ |
402 | X(i)=x(j); |
403 | Y(i)=U(j+np); |
404 | Ux270(i)=U(j); |
405 | Uy270(i)=Y(i); |
406 | end |
407 | hPlot=plot(X,Y,'-.');grid on;xlabel('X, m');ylabel('Uy, m') |
408 | set(hPlot,'LineWidth',2); |
409 | |
410 | figure |
411 | |
412 | % ВЫЧИСЛЕНИЕ КОМПОНЕНТОВ ДЕФОРМАЦИИ В СРЕДИННЫХ ТОЧКАХ |
413 | % ЭЛЕМЕНТОВ (eps=0.0, eta=0.0) |
414 | |
415 | epsX0=zeros(ever*egor); |
416 | epsY0=zeros(ever*egor); |
417 | gamXY0=zeros(ever*egor); |
418 | |
419 | for n=1:ever |
420 | for m=1:egor |
421 | n1=2*m-1+(3*egor+2)*(n-1); n2=n1+1; n3=n1+2; |
422 | n4=n3+2*egor-m+1; n5=n3+(3*egor+2); |
423 | n6=n5-1; n7=n5-2; n8=n4-1; |
424 | |
425 | dJ=abs((x(n1)-x(n3))*(y(n3)-y(n7)))/4; % якобиан |
426 | i=m+egor*(n-1); |
427 | epsX0(i)=(U(n4)-U(n8))/(x(n4)-x(n8)); |
428 | epsY0(i)=(-U(n2+np)+U(n6+np))/(y(n6)-y(n2)); |
429 | gamXY0(i)=-(-U(n2)+U(n6))/(x(n4)-x(n8))+(U(n4+np)-U(n8+np))/(y(n6)-y(n2)); |
430 | end |
431 | end |
432 | |
433 | % НОМЕРА ЭЛЕМЕНТОВ ПОСТРОЧНО |
434 | NEv=zeros(ever,egor); |
435 | for n=1:ever |
436 | for m=1:egor |
437 | i=m+egor*(n-1); |
438 | NEv(n,m)=i; |
439 | end |
440 | end |
441 | |
442 | % НОМЕР СРЕДИННОГО УЗЛА, ДЛЯ КОТОРОГО ВЫЧИСЛЯЮТСЯ ДЕФОРМАЦИИ % И НАПРЯЖЕНИЯ, |
443 | % СОВПАДАЕТ С НОМЕРОМ ЭЛЕМЕНТА, РАСПОЛОЖЕННОГО В МАССИВЕ NEv(n,m) |
444 | |
445 | X=zeros(egor); Y25=zeros(egor); Y200=zeros(egor); Y400=zeros(egor); |
446 | for m=1:egor |
447 | n=Nv(1,2*m); |
448 | X(m)=x(n); |
449 | i25=NEv(23,m); |
450 | i200=NEv(20,m); |
451 | i400=NEv(18,m); |
452 | Y25(m)=epsX0(i25); % 2,5 cм |
453 | Y200(m)=epsX0(i200); % 20,0 cм |
454 | Y400(m)=epsX0(i400); % 40,0 см |
455 | end |
456 | hPlot=plot(X,Y25,'-',X,Y200,'--',X,Y400,'-.');... |
457 | grid on;xlabel('X, m');ylabel('epsX') |
458 | set(hPlot,'LineWidth',2); |
459 | figure |
460 | |
461 | X=zeros(egor); Y25=zeros(egor); Y200=zeros(egor); Y400=zeros(egor); |
462 | for m=1:egor |
463 | n=Nv(1,2*m); |
464 | X(m)=x(n); |
465 | i25=NEv(23,m); |
466 | i200=NEv(20,m); |
467 | i400=NEv(18,m); |
468 | Y25(m)=epsY0(i25); % 2,5 cм |
469 | Y200(m)=epsY0(i200); % 20,0 cм |
470 | Y400(m)=epsY0(i400); % 40,0 см |
471 | end |
472 | hPlot=plot(X,Y25,'-',X,Y200,'--',X,Y400,'-.');... |
473 | grid on;xlabel('X, m');ylabel('epsY') |
474 | set(hPlot,'LineWidth',2); |
475 | figure |
476 | |
477 | % ВЫЧИСЛЕНИЕ КОМПОНЕНТОВ НАПРЯЖЕНИЙ Sigma В СРЕДИННЫХ ТОЧКАХ |
478 | % ЭЛЕМЕНТОВ (eps=0.0, eta=0.0) |
479 | SigX=zeros(ever); SigY=zeros(ever); |
480 | k=0; |
481 | for n=1:ever |
482 | for m=1:egor |
483 | n1=2*m-1+(3*egor+2)*(n-1); n2=n1+1; n3=n1+2; |
484 | n4=n3+2*egor-m+1; n5=n3+(3*egor+2); |
485 | n6=n5-1; n7=n5-2; n8=n4-1; |
486 | |
487 |
% ФОРМИРОВАНИЕ МАТРИЦЫ УПРУГОСТИ ДЛЯ СЛУЧАЯ ПЛОСКОЙ |
488 | if((n>=1)&(n<=7)) |
489 | e1=E(6); nu1=nu(6); % ГРУНТ СУГЛИНОК ЛЕГКИЙ |
490 | end |
491 | if((n>=8)&(n<=12)) |
492 | e1=E(5); nu1=nu(5); % ГРАВИЙНО-ПЕСЧАНАЯ СМЕСЬ |
493 | end |
494 | if((n>=13)&(n<=16)) |
495 | e1=E(4); nu1=nu(4); % ЩЕБЕНОЧНАЯ СМЕСЬ |
496 | end |
497 | if((n>=17)&(n<=20)) |
498 | e1=E(3); nu1=nu(3); % ЩЕБЕНОЧНО-ПЕСЧАНАЯ СМЕСЬ 8% |
499 | end |
500 | if((n>=21)&(n<=22)) |
501 | e1=E(2); nu1=nu(2); % АБ КЗ ПОРИСТЫЙ |
502 | end |
503 | if(n==23) |
504 | e1=E(1); nu1=nu(1); % АБ МЗ ПЛОТНЫЙ |
505 | end |
506 | |
507 | D=zeros(3,3); |
508 | D(3,3)=e1/(2*(1+nu1)); D(2,2)=2*D(3,3)*(1-nu1)/(1-2*nu1); |
509 | D(1,1)=D(2,2); D(1,2)=2*D(3,3)*nu1/(1-2*nu1); |
510 | D(2,1)=D(1,2); |
511 | |
512 | X1=x(n1); X2=x(n3); Y2=y(n3); Y4=y(n7); |
513 | dJ=abs((X1-X2)*(Y2-Y4))/4; % Якобиан |
514 | G=e1/(2*(1+nu1)); |
515 | |
516 | i=m+egor*(n-1); |
517 | % В МПа |
518 | SigmaX0(i)=(D(1,1)*epsX0(i)+D(1,2)*epsY0(i))/100000; |
519 | SigmaY0(i)=(D(2,1)*epsX0(i)+D(2,2)*epsY0(i))/100000; |
520 | TauXY0(i)=D(3,3)*gamXY0(i)/100000; |
521 | if(m==30) |
522 | k=k+1; |
523 | SigX(n)=SigmaX0(i); |
524 | SigY(n)=SigmaY0(i); |
525 | end |
526 | end |
527 | end |
528 | |
529 | % ПОСТРОЕНИЕ ГРАФИКОВ КОМПОНЕНТ НАПРЯЖЕНИЙ SigmaX0, SigmaY0 И |
530 | % TauXY0 В СРЕДИННЫХ УЗЛАХ ЭЛЕМЕНТОВ |
531 | |
532 | X=zeros(egor); Y23=zeros(egor); Y22=zeros(egor); Y20=zeros(egor); |
533 | Y18=zeros(egor); Y14=zeros(egor); Y12=zeros(egor); |
534 | for m=1:egor |
535 | n=Nv(1,2*m); |
536 | X(m)=x(n); |
537 | i23=NEv(23,m); % h=2,5 см |
538 | i22=NEv(22,m); % h=7,5 см |
539 | i20=NEv(20,m); % h=20,0 см |
540 | i18=NEv(18,m); % h=40,0 см |
541 | i14=NEv(14,m); % h=87,5 см |
542 | i12=NEv(12,m); % h=120,0 см |
543 | Y23(m)=SigmaX0(i23); % 2,5 cм |
544 | Y22(m)=SigmaX0(i22); % 7,5 cм |
545 | Y20(m)=SigmaX0(i20); % 20,0 cм |
546 | Y18(m)=SigmaX0(i18); % 40,0 cм |
547 | Y14(m)=SigmaX0(i14); % 87,5 cм |
548 | Y12(m)=SigmaX0(i12); % 120,0 cм |
549 | end |
550 | hPlot=plot(X,Y22,'-',X,Y20,'--',X,Y18,'-.',X,Y14,'-r');... |
551 | grid on;xlabel('X, m');ylabel('SigmaX, МПа') |
552 | set(hPlot,'LineWidth',2); |
553 | figure |
554 | |
555 | X=zeros(egor); Y23=zeros(egor); Y22=zeros(egor); Y20=zeros(egor); |
556 | Y18=zeros(egor); Y14=zeros(egor); Y12=zeros(egor); |
557 | for m=1:egor |
558 | n=Nv(1,2*m); |
559 | X(m)=x(n); |
560 | i23=NEv(23,m); % h=2,5 см |
561 | i22=NEv(22,m); % h=7,5 см |
562 | i20=NEv(20,m); % h=20,0 см |
563 | i18=NEv(18,m); % h=40,0 см |
564 | i14=NEv(14,m); % h=87,5 см |
565 | i12=NEv(12,m); % h=120,0 см |
566 | Y23(m)=SigmaY0(i23); % 2,5 cм |
567 | Y22(m)=SigmaY0(i22); % 7,5 cм |
568 | Y20(m)=SigmaY0(i20); % 20,0 cм |
569 | Y18(m)=SigmaY0(i18); % 40,0 cм |
570 | Y14(m)=SigmaY0(i14); % 87,5 cм |
571 | Y12(m)=SigmaY0(i12); % 120,0 cм |
572 | end |
573 | hPlot=plot(X,Y22,'-',X,Y20,'--',X,Y18,'-.',X,Y14,'-r');... |
574 | grid on;xlabel('X, m');ylabel('SigmaY, МПа') |
575 | set(hPlot,'LineWidth',2); |
576 | figure |
577 | |
578 | X=zeros(egor); Y23=zeros(egor); Y22=zeros(egor); Y20=zeros(egor); |
579 | Y18=zeros(egor); Y14=zeros(egor); Y12=zeros(egor); |
580 | for m=1:egor |
581 | n=Nv(1,2*m); |
582 | X(m)=x(n); |
583 | i23=NEv(23,m); % h=2,5 см |
584 | i22=NEv(22,m); % h=7,5 см |
585 | i21=NEv(21,m); % h=12,5 см |
586 | i20=NEv(20,m); % h=20,0 см |
587 | i18=NEv(18,m); % h=40,0 см |
588 | i16=NEv(16,m); % h=61,25 см |
589 | i14=NEv(14,m); % h=87,5 см |
590 | i12=NEv(12,m); % h=120,0 см |
591 | Y23(m)=TauXY0(i23); % 2,5 cм |
592 | Y22(m)=TauXY0(i22); % 7,5 cм |
593 | Y21(m)=TauXY0(i21); % 12,5 cм |
594 | Y20(m)=TauXY0(i20); % 20,0 cм |
595 | Y18(m)=TauXY0(i18); % 40,0 cм |
596 | Y16(m)=TauXY0(i16); % 61,25 cм |
597 | Y14(m)=TauXY0(i14); % 87,5 cм |
598 | Y12(m)=TauXY0(i12); % 120,0 cм |
599 | end |
600 | hPlot=plot(X,Y22,'-',X,Y20,'--',X,Y14,'-.');... |
601 | grid on;xlabel('X, m');ylabel('TauXY, МПа') |
602 | set(hPlot,'LineWidth',2); |
603 | |
604 | toc |
Приложение Б
(обязательное)
Исходный код подпрограммы BtDBnds
1 | function[BDB]=BtDBnds(psi,eta,D,X1,X2,Y1,Y3); |
2 | A=(-2/(X1-X2)); |
3 | B1=A*(1-eta)*(2*psi+eta)/4; B2= A*(-psi*(1-eta)); |
4 | B3= A*(1-eta)*(2*psi-eta)/4; B4= A*(1-eta^2)/2; |
5 | B5= A*(1+eta)*(2*psi+eta)/4; B6= A*(-psi*(1+eta)); |
6 | B7= A*(1+eta)*(2*psi-eta)/4; B8= A*(-(1-eta^2)/2); |
7 | B=(-2/(Y1-Y3)); |
8 | C1=B*(1-psi)*(psi+2*eta)/4; C2=-B*(1-psi^2)/2; |
9 | C3=B*(1+psi)*(-psi+2*eta)/4; C4=-B*eta*(1+psi); |
10 | C5=B*(1+psi)*(psi+2*eta)/4; C6=B*(1-psi^2)/2; |
11 | C7=B*(-(1-psi)*(psi-2*eta)/4); C8=-B*eta*(1-psi); |
12 | |
13 | |
14 | B=[B1 B2 B3 B4 B5 B6 B7 B8 0 0 0 0 0 0 0 0; |
15 | 0 0 0 0 0 0 0 0 C1 C2 C3 C4 C5 C6 C7 C8; |
16 | C1 C2 C3 C4 C5 C6 C7 C8 B1 B2 B3 B4 B5 B6 B7 B8]; |
17 | Bt=B'; |
18 | BtD=Bt*D; |
19 | BDB=BtD*B; |
Библиография
[1]Мартынов Н.Н., Иванов А.П. MATLAB 5.Х. Вычисления, визуализация, программирование. – М.:КУДИЦ-ОБРАЗ, 2000.-336 с.
[2] Коткин Г.Л., Черкасский В.С. Компьютерное моделирование физических процессов с использованием MATLAB:Учебное пособие/Новосиб. ун-т. Новосибирск, 2001. 173 с.
[3] Потемкин В.Г. Система инженерных и научных расчетов MATLAB 5.Х. В 2-х томах. –М.:ДИАЛОГ-МИФИ, 1999. 670 с.
[4] Тынкевич М.А. Численные методы.-Кемерово.:КузГТУ.1977.-122 с.
[5] Сегерлинд Л. Применение метода конечных элементов. М: Мир, 1979.– 392 с.
[6] СТ РК 1293-2004 Методы определения модуля упругости дорожных одежд нежесткого типа и их классификация.
