CC BY
26
Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки. 2017. Т. 21, № 1. С. 137-159
ISSN: 2310-7081 (online), 1991-8615 (print) d http://doi.org/10
УДК 517.958:539.3(1)
Изучение переходных процессов в нелинейно деформируемых средах на основе интегральных представлений и метода дискретных областей
В. А. Петушков
Институт машиноведения им. А. А. Благонравова РАН,
Россия, 101990, Москва, М. Харитоньевский пер., 4.
Аннотация
Для изучения переходных процессов в объемных нелинейно деформируемых средах разработаны методы моделирования, основанные на интегральных представлениях трехмерной краевой задачи упругой динамики, численных схемах высокого порядка аппроксимации границ и коллокационного приближения решения. Представлены обобщенные формулировки метода граничных интегральных уравнений, использующие фундаментальные решения статической упругости, уравнения состояния упругопластических сред с анизотропным упрочнением и разностные методы интегрирования по времени. Учитываются сложные истории комбинированного медленно меняющегося во времени и ударного нагружения составных кусочно-однородных сред при наличии зон локального возмущения решения. С использованием разработанного метода дискретных областей получены решения прикладных задач о распространении нелинейных волн напряжений в неоднородных средах. Приведены сравнения с решениями, полученными методом конечных элементов. Они подтверждают вычислительную эффективность разработанных алгоритмов, а также общность и полезность для практических целей предлагаемого подхода.
Ключевые слова: неоднородные среды, распространение волн, нелинейное деформирование и разрушение, метод граничных интегральных уравнений, метод конечных разностей, коллокационное приближение, метод дискретных областей, математическое моделирование.
Получение: 26 июня 2016 г. / Исправление: 15 октября 2016 г. / Принятие: 9 декабря 2016 г. / Публикация онлайн:
.14498/vsgtu1498
Статья
3 ©® Контент публикуется на условиях лицензии Creative Commons Attribution 4.0 International (https://creativecommons.org/licenses/by/4.0/deed.ru) Образец для цитирования
Петушков А. В. Изучение переходных процессов в нелинейно деформируемых средах на основе интегральных представлений и метода дискретных областей // Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки, 2017. Т. 21, № 1. С. 137-159. doi: 10.14498/vsgtu1498. Сведения об авторе
Владимир Алексеевич Петушков А http://orcid.org/0000-0002-9552-1399 доктор физико-математических наук, профессор; профессор; лаб. математического моделирования; e-mail: pva_imash@bk.ru
Введение. Изучение переходных процессов в трехмерных нелинейно деформируемых средах с локальными особенностями в геометрии, структуре или режимах нагружения является актуальной и достаточно сложной проблемой. Под действием ударных волн происходят изменения геометрии и объема, образование и локализация повреждений в зонах концентрации с последующим катастрофическим разрушением, страгивание и распространение существующих дефектов или образование множественных очагов разрушения для хрупких материалов и др. [1].
В большинстве случаев реальных нагрузок изменения в объеме (плотности) сред оказываются при этом незначительными и возникновение указанных предельных состояний происходит вследствие нелинейных деформаций сдвига. Очевидно, что анализ таких состояний возможен только с использованием численных методов, иногда в сочетании с экспериментальными методами, если ударным воздействиям подвергаются особо ответственные объекты техники.
Наиболее широко для решения подобных задач применяется метод конечных элементов (МКЭ), его достоинства неоспоримы и хорошо известны. Однако связанные с ним вычислительные затраты оказываются чрезмерными, если требуется детальное описание геометрии и представление решения с большими градиентами изменения полевых функций, например для составных конструкций или их систем при наличии зон концентрации в виде локальных особенностей в геометрии или дефектов, трещин и т. п. В силу малости этих зон соответствующие краевые задачи часто переходят в разряд внешних для неограниченных областей с граничными условиями «на бесконечности», требующими для их учета в МКЭ неординарных подходов.
Более эффективным в этих случаях оказывается применение метода граничных интегральных уравнений (МГИУ), чаще именуемого в зарубежной литературе методом граничных элементов (МГЭ). В основе его лежат интегральные представления соответствующих краевых задач, что позволяет свести их решения к многообразиям меньшей размерности и определять с более высокой точностью, чем в МКЭ, искомые полевые функций и их градиенты внутри расчетных областей [2-5]. Граничные условия на «бесконечности», или условия Зоммерфельда, выполняются при этом автоматически.
Основные трудности в применении МГИУ для рассматриваемых задач связаны, прежде всего, с наличием фундаментальных решений, а также необходимостью многократного по времени вычисления граничных интегралов и обращения несимметричных плотно упакованных матриц для формирования и решения систем разрешающих алгебраических уравнений большого порядка.
Поэтому очень важна разработка эффективных схем вычисления таких интегралов, особенно почти сингулярных, которые возникают при анализе переменных физического поля вблизи границы или изучении тонких пограничных структур. С этой целью могут быть использованы численные схемы, разработанные для решения не зависящих от времени краевых задач эллиптического типа [2, 5] и др. На гладких явно параметризованных границах обычно применяются численные схемы с квадратурами высокого порядка.
Все эти схемы становятся более эффективными при их использовании совместно с так называемыми быстрыми методами формирования и реше-
ния больших систем уравнений [3,6]. Когда эти методы применимы, число операций 2), обычно необходимых на формирование, и 3) на решение систем уравнений удается снизить до 0(Ж), где N — порядок системы. Тем самым расширяются возможности МГИУ для моделирования переходных процессов в реальных объектах техники.
К настоящему времени сформировались два основных направления в применении интегральных представлений для решения задач волновой динамики нелинейно деформируемых сред. В первом используются формулировки с зависящими от времени фундаментальными решениями Стокса и их аппроксимация по пространству и времени, во втором — с фундаментальными решения Кельвина для задач статики совместно с шаговыми методами по времени. Детальный обзор состояния и перспектив развития этих направлений представлен в недавних работах [3, 7-9].
Строгое математическое обоснование соответствующих формулировок метода, берущее начало с работ Купрадзе [10], изложено в ряде монографий [5,7]. Однако примеры успешного применения метода ограничиваются в основном двухмерными приближениями.
Для ядер, зависящих от времени, формулировка метода оказывается достаточно сложной и имеет проблемы с устойчивостью. Использование не зависящих от времени статических фундаментальных решений совместно с шаговыми методами аппроксимации по времени делает формулировку менее сложной, но трудоемкой с вычислительной точки зрения, поскольку требует многократного вычисления объемных интегралов с инерционными силами и нелинейными функциями физического поля внутри расчетной области. Однако это компенсируется простотой ее реализации на ЭВМ [11].
В предлагаемой статье представлены теория и основные элементы компьютерной реализации МГИУ, определяющие эффективность его применения для моделирования нелинейной волновой динамики неоднородных деформируемых сред, подвергаемых ударным воздействиям. Формулировка метода основана на успешно использованных в [12], не зависящих от времени фундаментальных решениях в интегральном представлении соответствующих краевых задач и шаговых методах интегрирования во времени.
В решении используются методы дискретных областей со схемами высокого порядка аппроксимации границ и коллокационного приближения. Такой подход органически вписывается в современные вычислительные технологии для многопроцессорных систем [3]. В подтверждение его эффективности приведены решения модельных задач с результатами, полученными МКЭ, и задач, имеющих практическое значение.
1. Постановка задачи и основные соотношения МГИУ. Рассматривается изотропная вязко-упругопластическая среда объема V, занимающая в М3 произвольную, в общем случае многосвязную область Б, ограниченную поверхностью Липшица
п+1
дБ = У Бг,
г
где Б"+1 — поверхность, внешняя по отношению ко всем остальным.
Область Б (тело, конструкция) может быть неоднородной, составленной из различных материалов, и включать в себя локальные геометрические осо-
бенности, например типа углов, отверстий или поверхностей разрыва типа трещин 0е(ха), где ха, а = 1, 2,—локальная система координат, связанная с поверхностью разрыва, рис. 1. Здесь и далее используются обозначения и соглашения, принятые в тензорном исчислении.
Пусть среда деформируется под действием быстро меняющихся во времени массовых X](хг,Ь), хг € Б, поверхностных р](хг,Ь), хг € дБа сил и перемещений Uj (хг, Ь), заданных на части поверхности с хг € дБи, а также известного температурного поля Т(хг, Ь), хг € Б. Следовательно, дБ = дБиидБа, при этом дБи П дБа = 0, а Ь € Dt = (0,т), где т — продолжительность истории динамического нагружения.
Состояние среды в момент времени Ь = 0 может быть как недеформи-рованным, так и включать в себя начальные или остаточные напряжения и деформации, в том числе технологического происхождения.
Поведение среды при указанных воздействиях характеризуется полями перемещений и](хг,Ь), деформаций ¿]к(хг,Ь) и напряжений ](хг,Ь) или их скоростями, соответственно, и], ], . Для его описания ограничимся так называемым квазиакустическим приближением.
В этом случае актуальное давление р на фронте ударной волны полагается функцией только объема, которую можно принять в следующем виде:
Р = К(р - ро)/р,
где К = (2р+3А)/3 — так называемый объемный модуль, а р и Л — параметры изотропной среды (параметры Ламе). Ударная адиабата совпадает с изоэн-тропой расширения, а необратимость процессов деформирования обусловлена в основном вязкостью и пластическими деформациями сдвига.
Полагая деформации малыми, тензор скорости деформаций Коши ]
« -е „ -о
представим в виде суммы упругой еек и начальной составляющих:
1 2
причем ¿0°к, в свою очередь, может включать в себя как пластические деформации ], так и любые другие, например технологические ё*к. Отнесем к числу перечисленных и температурные деформации
¿¿к = аЬ]к 6,
— о (uj,k + uk,j) — j + j,
(1)
Xi
Sn+1
Рис. 1. Схема к постановке задачи [Figure 1. The problem statement scheme]
где а — коэффициент температурного расширения, 5гз — символ Кронекера, 9 — разность температур.
Соответствующее поле скоростей напряжений определяется обобщенным законом Гука
(?]к = Л]к1т (¿1т - ¿ш) , (2)
где тензор упругости для изотропной среды
Л]к1т = Л5гз 5к1 + У(5гк 5Ц + 5И5]к),
а параметры Л, у в общем случае могут зависеть от температуры. Представим их в форме, удобной для последующего использования:
Еи Е Л = Тл-ТГл-, У
(1 + v)(1 - 2ру ^ 2(1 + V)'
где Е и V — обычные (технические) модуль Юнга и коэффициент Пуассона.
Для вычисления скоростей пластических деформаций, входящих в (1), используем модифицированные соотношения теории течения упругопласти-ческих сред с анизотропным упрочнением [13,14]. В этом случае полная деформация ¿ук с учетом классических условия течения Мизеса и условий на-гружения и разгрузки представляется в виде
¿]к = ¿% + ¿¿зк + ¿зк + к,
где приращения пластической деформации —¿рк определяются из соотношения
—¿р
д¿3k = Е3к. (3)
Здесь Е3к — девиатор тензора деформаций:
Езк = езк - 2>у1Р^к = 2^ (Бзк - рзк) + —Ррк, Ег — интенсивность так называемой смещенной деформации:
(2 Л1/2
Ег = Е3к Е3к] , —¿3 — интенсивность приращения пластических деформаций: —¿3 = (3/лЕг - аП-1) (Зу + Ип-1)-1,
причем
11
е-Зк = ¿зк - з¿ц5зк, Бзк = азк - р5зк, а р = -згц — давление,
ТТ -р ТТ ¿аг
рзк = ^¿зк, Н = ¿¿¿р
3 \ 1/2
О с \ 2
аг, — I . БзкБз^ .
£4 =£
Интенсивность приращения пластических деформаций Ае1? определяется на каждом n-ном шаге нагружения деформируемой среды диаграммой деформирования Oi ~ ¿p для соответствующей температуры в и/или скорости деформаций ¿i для сред, чувствительных к скорости деформирования. При этом полагаются справедливыми следующие соотношения:
g(xi,t) = dg/dt = dg/dt и Ад = gAt, (4)
где g(xi,t) —любая переменная физического поля.
Следуя законам сохранения, уравнение волновой динамики нелинейно деформируемой среды запишем в обобщенном виде с учетом сил вязкого демпфирования:
(cf - cf) ui,ij + c2uj,ii + Xj = üj + Yüj, (xi,t) e D x Dt. (5)
Здесь запятая перед индексом обозначает производную по соответствующей координате, üj — вектор ускорения, y = c/p — коэффициент вязкого демпфирования. Свойства материала деформируемой среды представлены волновыми скоростями cf = (Л + 2ß) /р и cf = ß/p — соответственно для волн расширения/сжатия и сдвига, а начальные и нелинейные составляющие деформирования отнесены к массовым силам Xj(xi,t):
Xj(x%,t) = Xj — 22ekk,j — (l — 2f) ¿¿kj,i.
Для однозначного решения уравнения (4) необходимо располагать граничными и начальными условиями. Граничные условия, записанные в форме Неймана и Дирихле для усилий и перемещений, соответственно имеют вид
Pj = Ojknk = pj (xi, t) на dDa x Dt, üj = uj (xi, t) на dDu x Dt,
(6)
где п] — компонента вектора внешней нормали к поверхности дБ в точке хг.
Выражения для тензора напряжений а^к и вектора усилий р] можно записать в виде, удобном для последующего использования:
3 = р(о\ - 2с%)5]кщ,г + рс2(и],к + пк,з), (7)
Рз = О^к Пк = р(о\ - 2с2 )5^кПк Щ,г + рС22 (из,к Пк + Пк,] Пк). (8) Неоднородные, отличные от нуля, начальные условия имеют вид
и](хг, 0) = и0(хг), и](хг, 0) = У0(хг), хг е Б. (9)
Для описания медленно меняющихся во времени или квазистатических процессов деформирования, когда и](хк,Ь) ^ са, а = 1, 2, соответствующая нелинейная краевая задача следует из (4)—(9), если положить рй] = 0. В этом случае время Ь играет роль только индекса, характеризующего изменения в истории нагружения.
При наличии в области Б внутренних известной формы границ Г, обусловленных соединением разнородных материалов или составных тел, условия (5) дополняются следующим соотношением:
üj = üf, pj = pf на Г x Dt, (10)
где цифрами 1, 2 обозначены тела (материалы), находящиеся по обе стороны от границы Г, которая в данном случае не является обычной физической. Если соединение допускает относительное смещение тел, то внутренняя граница Г между ними рассматривается в виде первоначально совпадающих поверхностей Г1, Г2 тел в зоне контакта, рис. 1. В отличие от статики, в задачах динамического деформирования вместо усилий рз слева и справа от границы Г используются их импульсы.
Отрыву (расслоению) общих точек поверхностей Г1, Г2 соответствуют условия равенства в них нулю усилий:
р1 = р2 = 0,
а скольжению с трением — следующие условия на (Г1'2 х Б^):
иП = й2п, рП + 'р2п = 0, р\ + р2т = 0, -крП2 < Р1, р2 < крП2, (11)
где п и т отмечены соответственно нормальное и касательные направления к поверхностям Г1, Г2 в точке контакта, к — коэффициент трения. В этом случае краевая задача (4)—(10) оказывается нелинейной и при упругом (линейном) поведении деформируемой среды.
2. Интегральные представления краевой задачи. Используя метод взвешенных невязок или обобщенную теорему взаимности Бетти (иначе — тождество Грина), краевой задаче (4)—(10) можно поставить в соответствие пространственно-временное интегральное уравнение с зависящими от времени фундаментальными решениями для смещений и сил. В основе его лежит представление решения щ(хг,Ь) задачи в виде зависящего от времени поверхностного потенциала с заданной плотностью и свойства непрерывности подобных потенциалов [3,5,10]:
гь
Ог3 (х)щ (х,г) = - / РгЗ (х, у,т) * Щ (у,т)(,8у(И + .70 Л дБ
+ и 13(х,Ь; у,т) * рз(у,т)(8у(М+
%з\
о
+ / / Игз (х,г; у,т) * {Хз (у,г) - из (у,т) - чиз (у,т)} (Ьу(Ъ+
1,
+ / / ^зкг(х,Ъ у, г) * £%(у, т)(х,г; у,т)(Уу(П. (12)
/о
Напряжения (х,Ь) для точек внутри области Б определяются соответствующим интегральным представлением, однако они могут быть получены с использованием (10) и уравнения состояния для деформируемой среды (2), что оказывается более эффективным с вычислительной точки зрения. Здесь у € дБ — точка поверхности; х € М3\дБ — точка (источник) внутри области Б; Ур, V — соответственно области пластического деформирования и/или ненулевой разности температур и отличными от нуля объемными силами; €13(х) = 5%з, если х € Б, и €13(х) = 0.5813, если х € дБ, при условии ее гладкости; знак (*) означает риманову свертку по времени, которая в этом случае
обладает всеми групповыми свойствами и определяется соотношениями
Динамические фундаментальные решения Р^, и ^^кг, входящие в (12), являются сингулярными при |х — у| ^ 0, обладают свойствами симметрии, взаимности, причинности и не зависят от начала отсчета по времени. Кроме того, они автоматически удовлетворяют условию излучения Зоммерфельда для бесконечных областей [5]. В явном виде выражения для них достаточно громоздкие и здесь в целях экономии места не приводятся, их можно найти, например, в [2].
Решения интегральных уравнений для перемещений (хг,г) и напряжений а^ (х,г) ищутся в пространстве Ь2(В х В^) интегрируемых с квадратом функций и их обобщенных производных и строятся численно МГИУ с использованием коллокационного приближения и аппроксимаций по пространству и времени в В х В^.
Однако этот обычный в МГИУ подход, основанный на зависящих от времени фундаментальных решениях, имеет известные проблемы с устойчивостью и оказывается весьма трудоемким и малоэффективным с вычислительной точки зрения [3,9,11]. Причина заключается в достаточно сложных выражениях для ядер, необходимости многократного вычисления интегралов свертки Римана (12) с обеспечением при этом причинно-следственной связи и устойчивости в течение всего времени интегрирования, а также объемных интегралов для изменяющихся во времени нелинейных деформаций, инерционных и вязких сил.
Более простым и эффективным оказывается подход, основанный на применении фундаментальных решений трехмерной задачи статики и шаговых методов интегрирования по времени [11]. Примеры его успешного использования для решения трехмерной задачи упругой динамики во временной и частотной областях были показаны в [12]. Здесь мы обобщим этот подход, включив учет нелинейного поведения деформируемых сред, как в [14,15].
Интегральное представление для краевой задачи (5)-(10) в этом случае формально имеет вид, аналогичный (12):
Однако входящие в него ядра — тензоры , Р^, — теперь соответствуют фундаментальному решению Кельвина—Сомильяны и определяются достаточно простыми выражениями:
f (хг, г) * д(хг, г) = 0 У(хг, г) е в х в1.
(13)
Щ, = [Ш»(1 — и)т]-1 [(3 — 4и)6г, + Т,гТ,з] ,
1 (дт
Рг, = — [8п(1 — V)т2 - { тп [(! — )5г3 + 3т,гт,з]
— (1 — 2и)(т,гп, — т,, пг)}, (14)
Ъ]кг = — [8п(1 — и)т2] 1 [(1 — 2и)(§г,т,к + 6гкт,, — 5,кт,г) + 3т,гт,,т,к],
дт
т = Ух — у\\ и т,г = — = т-1(Уг — хг).
Существование и единственность решений уравнения (13) с условиями (6) для задач упругости, в нашем случае на каждом шаге по времени, установлено в классе непрерывных по Гельдеру функций С 1,а(В) С Ш^В), где Ш2^(В) — пространство Соболева с параметром 0 < а ^ 1 [7,16]. Для нерегулярных областей, в том числе с поверхностями разрыва типа трещины, такие исследования выполнены, например, в [16] на основе теории псевдодифференциальных операторов и метода Винера—Хопфа.
Следует отметить, что уравнение (13) с выражениями для ядер (14) применимо только для ограниченных сред, поскольку условия излучения Зо-ммерфельда, как известно, в этом случае не выполняются. Для рассматриваемых неоднородных сред уравнение (13) записывается отдельно для каждой из подобластей, входящих в область В и имеющих общие границы. В целом для области соответствующие уравнения объединяются с помощью условий (10), (11) на этих границах.
Выражение для напряжений а,к(хг,г) может быть получено дифференцированием (13) с учетом (2)-(4) и для точек внутри области В имеет вид
аг, (х,г) = щ,к (х,у)рк (у,г)йв — Ркг, (х,у)йк (у,г)(в+ ■¡до .¡до
+ щ,к(х,у)р{йк (у,г) + чй к (у,г) — Хк (у,г)} (V+ .Уу
7 — 5у р _ 4/л(1 + V) 15(1 — V) ^ 3(1 — V)
+ Ъцш (х,у)е°ы(у,1)дУ — ^— 1Л 2^. — 5г, аб, (15)
где
йкг, = —^кг] = г, йк1,1, Ркг, = г,Рк1,1. Тензор определяется соотношением
= 2р [3(1 — 2v)(5г, т,к т,\ + 6ыт,гт,,) +
+ 3v(6lгT,j т,к + 6,к т,ьт,г + 6гк т ,1т,, + т ,гт ,к) — 15т,гт,, т,к т,\+
+ (1 — 2V)(6гк6ц + 6,кб1г) + (1 — ^)6г,5Ы] (8п(1 — v)т3)-1.
Для вычисления напряжений на границе дВ области В выражение (15) неприменимо из-за сингулярного характера входящих в него интегралов. Однако они могут быть получены в точках хг е дВ через граничные усилия и производные от граничных перемещений по касательным направлениям [14]. Для
этого необходима параметризация геометрии границы дБ и аппроксимация решения на ее поверхности. Аппроксимация решения также необходима и внутри области Б для определения последних двух интегралов по внутренним объемам V, Ур в выражениях (13), (15) и параметров механики разрушения в окрестности трещины (ха).
Полагая рассматриваемые поверхности области Б кусочно-гладкими, численное решение уравнения (13) выполним с использованием М граничных элементов Ак так, чтобы
м
дБ и дБг = у Ак и Ак П Аг = 0, к = I. к=1
Параметрическое представление геометрии и искомых функций щ(хз,г) и р1(хз ,г) на каждом из них выберем подобно МКЭ в следующем виде:
хг (Сч)хг '; ,
1,1, о (16)
иг :=КГ(^К, рг := N(^Ж, <-у € Ак С М2.
Здесь N(^), г = 1, 2,..., 8, — квадратичная функция формы локальных координат € (-1,1), принадлежащая к лагранжеву семейству конечных элементов Б1к'т(дБн) в смысле [17].
Аналогичное параметрическое представление используется и для вычисления объемных интегралов внутри области Б. В этом случае функция формы N в (16) включает в себя третью локальную координату п € (-1,1), нормальную к двум первым , и принадлежит тому же семейству:
ь
N ((7, п) € Б^Б), Б и Бг = у Ак, г = 1, 2,..., 8,... 20,
к=1
для объемных элементов с 8 и/или 20 узлами. Следует отметить, что размер области с пластическими деформациями заранее не известен и устанавливается в процессе решения.
Заменяя далее интегралы по границе и объему соответствующими суммами по граничным и объемным элементам, вместо (13), (15) получаем их численные аналоги:
ни (г) + си(г) + ми(г) = вх (г) + ср (г) + Бео(г),
а(г) = с'р (г) - Ни (г) + в'х (г) + Б'ео(г) + Qe0(t) - с 'и (г) - м 'и (г), ( )
где М и С — соответственно матрицы масс и вязкого демпфирования, остальные матрицы и векторы являются стандартными в МГИУ.
Вычисление входящих в эти выражения интегралов основано на численных квадратурах Гаусса. Объем и точность вычислений при этом зависят от порядка формулы и минимального расстояния между точкой наблюдения и узлом рассматриваемого элемента Ак. Если точка наблюдения совпадает с одним из узлов граничного или объемного элементов, т. е. г ^ 0, интегралы в (13), (15) становятся сингулярными. Для их вычисления, как и в задачах
статики, используются переход соответственно к полярным и сферическим координатам и квадратуры Гаусса по радиальным и угловым координатам.
Сильно сингулярные интегралы с особенностью О (г-3) рассматриваются в смысле главного значения Коши. Их вычисление по объемным элементам строится путем перехода к сферическим координатам и выделения локальной зоны вблизи сингулярной точки, по которой проводится расчет аналитически по радиусу и численно с использованием квадратуры Гаусса по угловым координатам [14].
Прямое численное интегрирование уравнений (17) во времени может быть осуществлено с помощью стандартных разностных схем; среди них выберем неявную двухслойную разностную в-схему Ньюмарка, обычно применяемую в МКЭ для решения задач нелинейной динамики. Предпочтение такому выбору обычно отдают из-за наличия у схемы контролируемых характеристик устойчивости и возможности простого подключения различных итерационных схем для учета нелинейностей.
Согласно этой схеме на Б вводится равномерная разностная сетка с шагом Аг таким образом, что
А = {гп : гп = пАг; п € 0 < пАг < т}.
На каждом интервале времени [гп,гп+1 ] ускорение, входящее в (17), обычно полагается постоянным и равным среднему значению (осредненному ускорению):
и(г) = (1 - 7 Ук + 7&п+1,
где 7 — свободный параметр, принимаемый здесь равным 1/2.
Ускорения и скорости на временном слое (п + 1) могут быть вычислены через перемещения и(г) следующим образом:
1 \ . 2 •• 2и
ип+1 = ШРп+1 - Ш2> [ип + ипАг +{2> - в) Аг2и
вАг2 п+1 вАгК ^ , ,
1 1.. 1 Г ( 1 \ 2.. 1 ()
ип+1=2вАип+1 + ип + ^ипАг - ^^ ип + ипАг + - в) Аг2иг .
где в — второй свободный параметр, управляющий точностью и устойчивостью схемы Ньюмарка. Для безусловной устойчивости схемы необходимо, вообще говоря, выполнение условия 2в ^ 7 ^ 1/2.
Выражения (18) включают в себя текущие (актуальные) значения вектора перемещений. Тем самым учитывается нелинейный, зависящий от решения характер изменения инерционных и вязких сил в процессе динамического деформирования среды. Подставим эти выражения вместе с граничными условиями в уравнения (17) и сведем в левую часть все входящие в них неизвестные. Полученную разрешающую систему нелинейных уравнений МГИУ можно переписать в следующем виде:
[А]Х,п+1 = К+1 + [Б]{ёо}п+1, (19)
& = -[А']Хп+1 + К+1 + [Б']{ёо}п+1 + [С]гп+1. (20)
Таким образом, на каждом шаге по времени Аг требуется решение обычной нелинейной задачи квазистатики, детали которой подробно обсуждаются в [14,15].
Входящие в (19), (20) [А] и [А] —матрицы коэффициентов при векторе граничных неизвестных, Хп+1 — вектор, содержащий неизвестные усилия и перемещения на границе, векторы Рп+1 и Р'п+1 содержат члены, определяемые отличными от нуля заданными значениями граничных перемещений и усилий, а также действием объемных сил, включая инерционные и вязкие. Матрицы [Б] и [Б'] включают в себя объемные интегралы с начальными и нелинейными деформациями в (13), (15), квазидиагональная матрица [С] определяется наличием свободных от интеграла членов в (15). Выражения для перечисленных матриц и векторов стандартные в МГИУ и здесь не приводятся. Их можно найти, например, в [2].
Для неоднородных (составных) сред уравнения (19), (20) записываются отдельно для каждой из сред и на общей поверхности между ними Г = Г1 иГ2 вводится два слоя узлов и элементов. В случае идеального контакта на первом из них в качестве граничных условий задаются неизвестные перемещения и нулевые усилия, на втором — неизвестные усилия и нулевые перемещения. В случае проскальзывания с трением на первом слое задаются неизвестные компоненты перемещений — нормальные и касательные, а на втором — неизвестные нормальные усилия и касательные перемещения. При этом условия (10), (11) учитываются при формировании матриц [А] и [А] .
Моделирование нелинейного деформирования сред осуществляется здесь, как и в [13,14], с использованием на каждом шаге по времени уравнений состояния (2), (3) и классического метода последовательной линеаризации в форме начальных деформаций с применением известной схемы «упругий предиктор - пластический корректор». Приращения пластической деформации на каждой итерации определяются из решения (19), (20) с использованием соотношений (3) и экспериментально полученных диаграмм деформирования для материала каждой из рассматриваемых сред.
Большая размерность, плотная упакованность и несимметричная структура матриц и векторов в уравнениях (19), (20) делают вычислительные операции над ними чрезвычайно трудоемкими и представляют серьезное препятствие для применения МГИУ в математическом моделировании переходных процессов в реальных крупномасштабных объектах. В связи с этим для ускорения формирования и решения систем уравнений МГИУ в последние годы разработаны так называемые быстрые методы, такие как РЫМ, АСА и другие [6], а также используемый ниже метод дискретных областей [14,18].
В основе этого метода лежит представление исходной расчетной области Б в виде совокупности подобластей с искусственными границами между ними и граничными условиями типа (10). В результате матрицы уравнений становятся блочно-симметричными с ленточной структурой, что позволяет эффективно использовать прямые алгоритмы гауссовского типа для их обращения [19].
3. Вычисление напряжений на границе тела, метод дискретных областей. Интегральное выражение (14) неприменимо для вычисления напряжений на границе нелинейно деформируемого тела из-за сингулярного характера входящих в него интегралов с особенностью О(г-2, г-3). Вместо него
воспользуемся вычисленными на границе усилиями и касательными производными от граничных перемещений и соотношениями (7), (8).
Перемещения й,(хг,г) и их производные по локальным координатам с учетом аппроксимации (16), а также граничные усилия р,(хг,г) определяются из решения уравнений (19), (20) для каждого шага по времени. Через них с использованием (7) могут быть вычислены напряжения на границе дВ рассматриваемой области В.
Такой прием оказывается весьма эффективным с вычислительной точки зрения, поскольку обладает высокой точностью и не требует дополнительного вычисления производных от интегралов по области пластического деформирования, нет также необходимости вычисления матрицы В' в выражении (20) в процессе решения нелинейной задачи. Его легко обобщить на нелинейно деформируемые среды, расчетные области которых состоят из элементов с сильно различающимися свойствами и/или размерами или включают в себя локальные особенности, например, типа произвольно ориентированных трещин.
В этом случае область В, занимаемая деформируемой средой, рассматривается в виде совокупности N дискретных подобластей 0,г
N
В = У пг и П1 П пт = 0, I = т,
г=1
образованных произвольно ориентированными фиктивными границами (см. рис. 1). На каждой из них задаются условия типа (10) или (11) и при необходимости учитываются, например, взаимодействия берегов трещины [20]. Кусочно-неоднородные среды с различными свойствами материалов и реальными разделяющими их границами учитываются при этом естественным образом.
При таком подходе появляется возможность объединения предлагаемого метода с другими, например, МКЭ для областей с сильной нелинейностью, или МГИУ в его прямой формулировке (12) с областью, включающей «бесконечность», для учета граничных условий Зоммерфельда [11], или даже с экспериментальным методом определения переменных физического поля на поверхности.
Изложенные в статье методы и разработанные на их основе алгоритмы реализованы в компьютерной программе, как развитие программы [11] разработанной ранее применительно к изучению квазистатических процессов нелинейного деформирования трехмерных сред — конструкций и их элементов. В качестве пре- и постпроцессора в пакет включена программа ОГО [21] с широкими возможностями автоматического построения расчетных схем МГИУ и визуализации результатов.
4. Численные результаты. Для обоснования вычислительной эффективности предложенного подхода и компьютерной программы [10], развитой применительно к моделированию нелинейной волновой динамики, вначале приведем результаты решения МГИУ классической задачи о распространении упругой и упругопластической волны вдоль неоднородного стержня при его мгновенном нагружении в момент времени г = 0.
Задача рассматривается в объемной постановке (рис. 2), где показаны геометрия, условия закрепления и функция нагружения стержня. Длина стерж-
ня Ь = 10 см, поперечное сечение — квадрат с длиной стороны I = 0.2Ь. Стержень состоит из двух подобластей и с общей границей В, каждая со своим материалом, диаграммы деформирования которых приведены на рис. 2, Ь.
Механические свойства материалов: модули Юнга Е = 2.15 ■ 105 МПа, Е2 = 2.05 ■ 105 МПа; коэффициенты Пуассона и1 = 0.3, = 0.33, начальные плотности материалов р1 = р2 = 7.85 ■ 103 кг/м3. Функция нагружения — равномерно распределенное на торце давление = р0Н(¿), где Н(¿) — функция Хевисайда, значение р0 = 220 МПа, рис. 2, с.
В случае упругого поведения материала задача имеет аналитическое решение. Для нелинейного деформирования воспользуемся сравнительными результатами, полученными МКЭ с использованием программы [22]. Для аппроксимации решения на поверхности стержня выбраны квадратичные вось-миузловые граничные элементы, при этом объемные интегралы вычислялись с использованием аналогичных двадцатиузловых элементов. Изучались вопросы точности и вычислительной устойчивости МГИУ при последовательном увеличении числа граничных элементов для каждой подобласти (материала) от 88 до 1426 и выборе соответствующих значений шагов по времени с использованием известного критерия КФЛ — ^ДЬ/Н < 1/3, где Л — размер равномерной сетки, и в, выбираемого из диапазона в = 0.25 ^ 0.10, в соотношениях (18). Результаты, представленные ниже, получены с использованием в аппроксимации стержня 352 квадратичных элементов, в = 1/6 и шагом по времени 1.5 ■ 10-4 с.
Strain, £i (percentages) Time
Ь с
Рис. 2. Ударное нагружение неоднородного стержня: a — геометрия, b — диаграммы нелинейного деформирования материалов, c — функция нагружения [Figure 2. Impact loading of the inhomogeneous rod: a—the geometry, b—the diagrams of nonlinear deformation of the materials, c—the loading function]
Распространение упругопластической волны напряжения вдоль неоднородного стержня, полученное с использованием уравнения состояния (3) и диаграмм деформирования составляющих его материалов (рис. 2, b), показано на рис. 3. Здесь же приведены результаты упругого решения в сечениях A и C стержня и для сравнения результаты, полученные МКЭ. Для демонстрации объемности возникающих в стержне напряженных состояний в точке C приведены также другие значимые компоненты тензора напряжений. Из сравнения полученных результатов следует хорошее качество получаемых МГИУ решений.
Деформированная геометрия стержня и распределение в нем упругопла-стического напряжения 033 для момента времени t = 6.0 ■ 10-3 c от начала ударного нагружения показаны на рис. 4, a. Видно, что основная концентрация напряжения имеет место в зоне соединения разнородных материалов (сечение B на рис. 2, a). Здесь, по-видимому, и происходит накопление повреждений и разрушение стержня. Практический интерес представляет и характер взаимодействия распространяющих волн напряжения в стержне с границей этого соединения (рис. 4, b), где также приведено аналогичное решение, полученное МКЭ. Сравнение результатов и в этом случае свидетельствует о хорошей точности и разрешающей способности МГИУ в описании больших градиентов изменения полевых функций.
С использованием результатов вычислительного эксперимента, выполненного на объемном стержне, МГИУ решена более сложная задача о кинетике напряженно-деформированных и возникновении предельных состояний в зонах резкого локального изменения геометрии в трехмерных конструкциях при ударных воздействиях.
В объемной постановке выполнено моделирование переходных нелинейных процессов деформирования в биметаллическом образце с острым надрезом под действием внезапно приложенной к одному из торцов ударной нагрузки (рис. 5, a). Надрез расположен вблизи зоны соединения разнородных материалов, которые сильно отличаются друг от друга своими деформационными свойствами (рис. 5, b). Условия нагружения представлены на рис. 5, с. Уровни нагружения и механические свойства материалов аналогичны приведенным в предыдущем примере. На этом же рисунке показаны сечения образца, представляющие интерес для анализа указанных состояний, и распределение зон пластических деформаций в центральном сечении ABCD образца для момента времени t = 8.83 ■ 10-4 c от начала ударного нагружения. В аппроксимации расчетной области полосы использовано 784 объемных и 398 граничных элементов.
На рис. 6 представлена информация о максимальном раскрытии берегов надреза (по существу трещины) в точке N образца (см. рис. 5, d) для упругого и упругопластического поведения материалов.
На рис. 7, a приведено распределение упругопластического напряжения в самом образце с сечением вблизи вершины надреза для момента времени t = 8.83 ■ 10-4 c и изменение этого напряжения вдоль траектории S, построенное для различных моментов времени от начала нагружения (рис. 7, b). Распределение по форме характеризует собой известный режим «с обострением», возникающий в зоне разнородного соединения материалов с концентратором напряжений, и результаты его очень важны для оценки прочности
Рис. 3. Распространение волн напряжений в точках A и C стержня: 1 — упругое решение МГИУ для <7зз, 1' —пластическое решение МГИУ для 733; 2 — упругое решение МКЭ для 733, 22 —пластическое решение МКЭ для 733; 3 — решение МГИУ для 7ц, 33 —решение
МГИУ для 722 (онлайн в цвете)
[Figure 3. Stress waves distribution in the points A and C of the rod: the line 1 is an elastic
BIEM solution for the stress 733, the line 1' is a plastic BIEM solution for the stress 733; the line 2 is an elastic FEM solution for the stress 733, the line 2' is a plastic FEM solution for the stress 733; the line 3 is a BIEM solutions for the stress 7ц, the line 33 is a BIEM solutions for the stress 722 (color online)]
С В_ _ _A
Ж ' 1 1 .198
-63? Ж ■578 1=0 006c J
a
Time, t (secounds) b
Рис. 4. Распределение упругопластических напряжений 733 (МПа) в стержне: a — геометрия и напряжения в момент времени t = 6 • 10-3 c от начала нагружения, b — изменение напряжения во времени в точке B по МГИУ (линия 1) и МКЭ (линия 1) (онлайн в цвете)
[Figure 4. Distribution of the elastoplastic stresses 733 (MPa) in the rod: a—the geometry and the stresses at time t = 6.0 • 10-3 s from the beginning of loading; b—stress evolution in time at point B by BIEM (line 1), and by FEM (line 2) (color online)]
to-
rn
1.9 cm
0.1 cm
2 cm
JL
n2
1 cm
0 4 8 12 Strain, £i (percentages)
Ei = 2.15 • 105 MPa лВ VI = 0.3
pi = 7.85 ■ 103 kg/m;
B2 = 2.05 • 105 MPa vï = 0.33
P2 = 7.85 ■ 103 kg/m3
Рис. 5. Образец с надрезом при мгновенном нагружении: a — геометрия и условия закрепления, b — функция нагружения, c — диаграммы деформирования материалов, d — зоны пластических деформаций в центральном сечении для момента времени t = 8.83 • 10-4 c
(онлайн в цвете)
[Figure 5. A bi-material plate with edge notch/crack along interface at instantaneous loading: a—the geometry and the fixing conditions, b—the loading function, c—the diagrams of deformation of the materials, d—the plastic deformation zones in the central section at time t = 8.83 • 10-4 s (color online)]
-12
0.5 • 10
2.5 • 10"
1 • 10"3 1.5 • 10"3 Time, t (secounds)
Рис. 6. Максимальное раскрытие берегов надреза u^ во времени: пластическое (линия 1) и упругое (линия 2) решения по МГИУ [Figure 6. Maximum opening of the bi-material sample edge notch uN in time: the plastic (line Î) and elastic (line 2) solutions by BIEM]
400
(2 200
сп
-200
x...................
\|7\ \
5 17 \ \
0 0.4 0.8 1.2 1.6 2.0
Cross-section, Si (cm) с
Рис. 7. Максимальное напряжение 733 в вершине надреза образца: a — распределение для момента времени t = 8.83 • 10-4 c, b — кинетика распределения напряжений для различных моментов времени t вдоль траектории S, c — кинетика распределения напряжений для различных моментов времени t вдоль сечения Si; 1 — t = 3.65 • 10-4 c, 2 — t = 7.28 • 10-4 c, 3 — t = 8.83 • 10-4 c, 4 — t = 1.46 • 10-3 c, 5 — t = 2.18 • 10-3 c (онлайн в цвете)
[Figure 7. The maximum stress 733 in the top of the edge notch of the bi-material specimen: a—distribution at time t = 8.83 • 10-4 s, b—stress distribution kinetics for various times t along the trajectory S, c—stress distribution kinetics for various times t along the cross-dection S1; line 1 for t = 3.65 • 10-4 s, line 2 for t = 7.28 • 10-4 s, line 3 for t = 8.83 • 10-4 s, line 4 for t = 1.46 • 10-3 s, line 5 for t = 2.18 • 10-3 s (color online)]
самого соединения.
Не меньший интерес представляет распределение напряжений и в другом сечении образца на границе Г соединения разнородных материалов. Изменения напряжений во времени вдоль выбранной траектории Si (рис. 7, с) также свидетельствуют о возможности разрушения образца в этом сечении.
Заключение. Для изучения переходных процессов в объемных нелинейно деформируемых средах при динамических, в том числе ударных, воздействиях разработаны методы моделирования, основанные на интегральных представлениях краевых задач, численных схемах высокого порядка аппроксимации границ и коллокационного приближения. Представлены соответствующие формулировки МГИУ с фундаментальными решениями статической упругости и разностными методами по времени, позволяющие учитывать неоднородности и локальные особенности сред.
Также обсуждаются основные аспекты численной реализации выполненных разработок. Изучены вопросы точности, вычислительной устойчивости и эффективности метода и компьютерной программы. Приведены примеры моделирования линейных и нелинейных процессов распространения волн деформаций и напряжений в трехмерных средах и дано их сравнение с численными решениями, полученными МКЭ. Они подтверждают вычислительную эффективность разработанных алгоритмов, а также общность и полезность для практических целей предлагаемого подхода. Конкурирующие интересы. У меня нет конкурирующих интересов. Авторская ответственность. Я несу полную ответственность за предоставление окончательной версии рукописи в печать. Окончательная версия рукописи мною одобрена.
Финансирование. Исследование не имело финансирования.
Библиографический список
1. Майборода В. П., Кравчук А. С., Холин Н. Н. Скоростное деформирование конструкционных материалов. М.: Машиностроение, 1986. 264 с.
2. Wrobel L. C., Aliabadi M. H. The Boundary Element Method. vol. 1: Applications in Thermo-Fluids and Acoustics. New York: John Wiley & Sons, Ltd., 2007. 1066 pp.
3. Liu Y. J., Mukherjee S., Nishimura N., Schanz M. at all Recent Advances and Emerging Applications of the Boundary Element Method// Appl. Mech. Rev., 2012. vol.64, no. 3, 030802. 38 pp. doi: 10.1115/1.4005491.
4. Hayami K. Variable transformations for nearly singular integrals in the boundary element method // Publ. Res. Inst. Math. Sci., 2005. vol.41, no. 4. pp. 821-842. doi: 10.2977/prims/ 1145474596.
5. Hayami K., Costabel M. Time-dependent problems with the boundary integral equation method/ Encyclopedia of Computational Mechanics. New York: John Wiley & Sons, Ltd., 2004. pp. 703-721. doi: 10.1002/0470091355.ecm022.
6. Rjasanow S., Steinbach O. The Fast Solution of Boundary Integral Equations / Mathematical and Analytical Techniques with Applications to Engineering. Heidelberg: Springer, 2007. xi+279 pp. doi: 10.1007/0-387-34042-4.
7. Hsiao G. C., Wendland W. L. Boundary Integral Equations / Applied Mathematical Sciences. vol.164. Berlin: Springer, 2008. xix+618 pp. doi: 10.1007/978-3-540-68545-6.
8. Hatzigeorgiou G. D. Dynamic Inelastic Analysis with BEM: Results and Needs / Recent Advances in Boundary Element Methods; eds. G. D. Manolis, D. Polyzos. Berlin: Springer, 2009. 193-208 pp.. doi: 10.1007/978-1-4020-9710-2_13.
9. Hatzigeorgiou G. D., Beskos D. E. Dynamic inelastic structural analysis by the BEM: A review// Engineering Analysis with Boundary Elements, 2011. vol.35, no. 2. pp. 159-169. doi: 10.1016/j.enganabound.2010.08.002.
10. Купрадзе В. Д., Бурчуладзе Т. В. Динамические задачи теории упругости и термоупругости/ Итоги науки и техн. Сер. Соврем. пробл. мат., Т. 7. М.: ВИНИТИ, 1975. С. 163-294.
11. Soares D. Dynamic analysis of elastoplastic models considering combined formulations of the time-domain boundary element method // Engineering Analysis with Boundary Elements, 2015. vol.55. pp. 28-39. doi: 10.1016/j.enganabound.2014.11.014.
12. Петушков В. А., Потапов А. И. Численные решения трехмерных динамических задач теории упругости / Сб. докладов Седьмого Всесоюзного съезда по теоретической и прикладной механике. М.: МГУ, 1991. С. 286-287.
13. Petushkov V. A., Shneiderovich R. M. Thermoelasticplastic deformation of corrugated shells of revolution at finite displacements // Strength of Materials, 1979. vol. 11, no. 6. pp. 578585. doi: 10.1007/bf00770100.
14. Петушков В. А. Численная реализация метода граничных интегральных уравнений применительно к нелинейным задачам механики деформирования и разрушения объемных тел / Моделирование в механике: Сб. научных трудов ИТПМ СО АН СССР. Т. 3(20). Новосибирск, 1989. С. 133-156.
15. Петушков В. А. Моделирование нелинейного деформирования и разрушения неоднородных сред на основе обобщенного метода интегральных представлений // Матем. моделирование, 2015. Т. 27, №1. С. 113-130.
16. Costabel M. Boundary Integral Operators on Lipschitz Domains: Elementary Results // SIAM J. Math. Anal., 1988. vol.19, no. 3. pp. 613-626. doi: 10.1137/0519043.
17. Strang G., Fix G. An Analysis of the Finite Element Method, 2nd edition, 2008. 400 pp.
18. Hsiao G. C., Steinbach O., Wendland W. L. Domain decomposition methods via boundary integral equations // Journal of Computational and Applied Mathematics, 2000. vol. 125, no. 1-2. pp. 521-537. doi: 10.1016/s0377-0427(00)00488-x.
19. Петушков В. А., Зысин В. И. Пакет прикладных программ МЕГРЭ-3Д для численного моделирования нелинейных процессов деформирования и разрушения объемных тел. Алгоритмы и реализация в ОС ЕС / Сб. Пакеты прикладных программ: Программное обеспечение математического моделирования. М.: Наука, 1992. С. 111-126.
20. Петушков В. А. Метод граничных интегральных уравнений в моделировании нелинейного деформирования и разрушения трехмерных неоднородных сред // Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки, 2014. №2(35). С. 96-114. doi: 10.14498/ vsgtu1292.
21. GiD—The Personal Pre and Post Processor (ver. 11). Barcelona: CIMNE, 1998.
22. Петушков В. А., Фролов К. В. Динамика гидроупругих систем при импульсном возбуждении / Динамика конструкций гидроаэроупругих систем. М.: Наука, 2002. 162-202 с.
Vestn. Samar. Gos. Tekhn. Univ., Ser. Fiz.-Mat. Nauki
[J. Samara State Tech. Univ., Ser. Phys. Math. Sci.], 2017, vol. 21, no. 1, pp. 137-159 ISSN: 2310-7081 (online), 1991-8615 (print) d http://doi.org/10.14498/vsgtu1498
MSC: 74G30, 74H25, 74J20
Transient dynamics of 3D inelastic heterogeneous media analysis by the boundary integral equation and the discrete domains methods
V. A. Petushkov
A. A. Blagonravov Mechanical Engineering Institute RAS, 4, M. Khariton'evskii per., Moscow, 101990, Russian Federation.
Abstract
For the study of transients in 3D nonlinear deformable media we devel-ope modeling methods which based on integral representations of 3D boundary value problem of elastic dynamics, numerical high-order approximation schemes of boundaries and collocation approximation of solutions. The generalized boundary integral equation method formulations using fundamental solutions of static elasticity, equation of state of elastoplastic media with anisotropic hardening and difference methods for time integration are represented. We take into account the complex history of combined slowly changing over time and impact loading of composite piecewise-homogeneous media in the presence of local perturbation solutions areas. With the use of this method and discrete domains method the solutions of applied problems of the propagation of non-linear stress waves in inhomogeneous media are received. Comparisons with the solutions obtained by the finite element method are represented also. They confirm the computational efficiency of the developed algorithms, as well as common and useful for practical purposes of the proposed approach.
Keywords: inhomogeneous media, wave propagation , nonlinear deformation and failure, boundary integral equation method, finite difference method, collocation approximation, subdomains method, mathematical simulation.
Received: 26th June, 2016 / Revised: 15th October, 2016 / Accepted: 9th December, 2016 / First online:
Competing interests. I have no competing interests.
Author's Responsibilities. I take full responsibility for submitting the final manuscript in print. I approved the final version of the manuscript.
Funding. The research has not had any sponsorship. Article
Q ©® The content is published under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/) Please cite this article in press as:
Petushkov V. A. Transient dynamics of 3D inelastic heterogeneous media analysis by the boundary integral equation and the discrete domains methods, Vestn. Samar. Gos. Tekhn. Univ., Ser. Fiz.-Mat. Nauki [J. Samara State Tech. Univ., Ser. Phys. Math. Sci.], 2017, vol. 21, no. 1, pp. 137-159. doi: 10.14498/vsgtu1498 (In Russian). Author's Details:
Vladimir A. Petushkov http://orcid.org/0000-0002-9552-1399
Dr. Phys. & Math. Sci.; Professor; Lab. of Mathematical Simulation; e-mail: pva_imash@bk.ru
References
1. Maiboroda V. P., Kravchuk A. S., Kholin N. N. Skorostnoe deformirovanie konstruktsion-nykh materialov [High-Velocity Deformation of Structural Materials]. Moscow, Mashinos-troenie, 1986, 264 pp. (In Russian)
2. Wrobel L. C., Aliabadi M. H. The Boundary Element Method, vol. 1, Applications in Thermo-Fluids and Acoustics. New York, John Wiley & Sons, Ltd., 2007, 1066 pp.
3. Liu Y. J., Mukherjee S., Nishimura N., Schanz M. at all Recent Advances and Emerging Applications of the Boundary Element Method, Appl. Mech. Rev., 2012, vol. 64, no. 3, 030802, 38 pp. doi: 10.1115/1.4005491.
4. Hayami K. Variable transformations for nearly singular integrals in the boundary element method, Publ. Res. Inst. Math. Sci., 2005, vol.41, no. 4, pp. 821-842. doi: 10.2977/prims/ 1145474596.
5. Hayami K., Costabel M. Time-dependent problems with the boundary integral equation method, In: Encyclopedia of Computational Mechanics. New York, John Wiley & Sons, Ltd., 2004, pp. 703-721. doi: 10.1002/0470091355.ecm022.
6. Rjasanow S., Steinbach O. The Fast Solution of Boundary Integral Equations, Mathematical and Analytical Techniques with Applications to Engineering. Heidelberg, Springer, 2007, xi+279 pp. doi: 10.1007/0-387-34042-4.
7. Hsiao G. C., Wendland W. L. Boundary Integral Equations, Applied Mathematical Sciences, vol.164. Berlin, Springer, 2008, xix+618 pp. doi: 10.1007/978-3-540-68545-6.
8. Hatzigeorgiou G. D. Dynamic Inelastic Analysis with BEM: Results and Needs, In: Recent Advances in Boundary Element Methods; eds. G. D. Manolis, D. Polyzos. Berlin, Springer, 2009, 193-208 pp.. doi: 10.1007/978-1-4020-9710-2_13.
9. Hatzigeorgiou G. D., Beskos D. E. Dynamic inelastic structural analysis by the BEM: A review, Engineering Analysis with Boundary Elements, 2011, vol.35, no. 2, pp. 159-169. doi: 10.1016/j.enganabound.2010.08.002.
10. Kupradze V. D., Burchuladze T. V. The dynamical problems of the theory of elasticity and thermoelasticity, J. Soviet Math., 1977, vol. 7, no. 3, pp. 415-500. doi: 10.1007/BF01091837.
11. Soares D. Dynamic analysis of elastoplastic models considering combined formulations of the time-domain boundary element method, Engineering Analysis with Boundary Elements, 2015, vol.55, pp. 28-39. doi: 10.1016/j.enganabound.2014.11.014.
12. Petushkov V. A., Potapov A. I. Numerical solutions of three-dimensional dynamic problems in elasticity theory, In: Book of Abstracts of the Seventh All-Union Congress on Theoretical and Applied Mechanics. Moscow, Moscow State Univ., 1991, pp. 286-287 (In Russian).
13. Petushkov V. A., Shneiderovich R. M. Thermoelasticplastic deformation of corrugated shells of revolution at finite displacements, Strength of Materials, 1979, vol. 11, no. 6, pp. 578-585. doi:10.1007/bf00770100.
14. Petushkov V. A. Numerical realization of the method of boundary integral equations as applied to nonlinear problems of mechanics of deformation and fracture of volume bodies, Model. Mekh., 1989, vol. 3(20), no. 1, pp. 133-156 (In Russian).
15. Petushkov V. A. Simulation of nonlinear deformation and fracture of heterogeneous media based on the generalized method of integral representations, Matem. Mod., 2015, vol. 27, no. 1, pp. 113-130 (In Russian).
16. Costabel M. Boundary Integral Operators on Lipschitz Domains: Elementary Results, SIAM J. Math. Anal., 1988, vol.19, no. 3, pp. 613-626. doi: 10.1137/0519043.
17. Strang G., Fix G. An Analysis of the Finite Element Method, 2nd edition, 2008, 400 pp.
18. Hsiao G. C., Steinbach O., Wendland W. L. Domain decomposition methods via boundary integral equations, Journal of Computational and Applied Mathematics, 2000, vol. 125, no. 12, pp. 521-537. doi: 10.1016/s0377-0427(00)00488-x.
19. Petushkov V. A., Zysin V. I. "MEGRE-3D" software package for numerical simulation of nonlinear processes of deformation and fracture 3D bodies, In: Algorithms and implementation in OS ES", Programmnoye obespecheniye matematicheskogo modelirovaniya [Software for Mathematical Simulation]. Moscow, Nauka, 1992, pp. 111-126 (In Russian).
20. V. A. Petushkov Boundary Integral Equation Method in the Modeling of Nonlinear Deformation and Failure of the 3D Inhomogeneous Media, Vestn. Samar. Gos. Tekhn. Univ. Ser. Fiz.-Mat. Nauki [J. Samara State Tech. Univ., Ser. Phys. Math. Sci.], 2014, no. 2(35), pp. 96-114 (In Russian). doi: 10.14498/vsgtu1292.
21. GiD—The Personal Pre and Post Processor (ver. 11). Barcelona, CIMNE, 1998.
22. Petushkov V. A., Frolov K. V. Dynamics of hydroelastic systems with pulsed excitation, In: Dinamika konstruktsii gidroaerouprugikh sistem [Dynamics of constructions of hydroaeroe-lastic systems]. Moscow, Nauka, 2002, 162-202 pp. (In Russian)
