Mathematical models of nonlinear viscoelasticity with operators of fractional integro-differentiation

Abstract


Based on the method of structural modelling and the Boltzmann-Volterra hypothesis of the hereditary elasticity deformable solid, the article considers linear and nonlinear fractional analogues of classical rheological models, such as Newton (the so-called model of Scott Blair), Voigt, Maxwell, Kelvin and Zener using the tool of fractional integro-differentiation of Riemann Liouville. The classes of nonlinear mathematical models are distinguished, for which the solution of the creep problem can be obtained explicitly in terms of known special functions. A technique for identifying the parameters of the proposed mathematical models is developed on the basis of known experimental data on uniaxial stretching of samples at different and constant load levels. In the presence of explicit solutions to the creep problem, the parameters of mathematical models are determined from the solution of the problem of approximating the experimental values of deformation using the method of least squares with subsequent refinement by the coordinate-wise descent method at all time points for all stress values in a series of experiments. For nonlinear mathematical models of viscoelastic deformation, which do not allow to find the solution of the creep problem in an explicit form, a method for determining the model parameters based on the coordinate-wise descent method with inversion at each step to the numerical solution of the defining integral equation has been developed. The method of identifying model parameters with operators of fractional integro-differentiation is realized on the example of creep of polyvinylchloride plastic compound. The values of the parameters for all the models studied are given, their adequacy to the experimental data is checked, and errors in the deviation of the calculated data from the experimental values are analyzed. As an example, a comparative analysis of the relative error in approximating the experimental creep curves by theoretical deformation values is made in the framework of a linear, nonlinear integrable and nonlinear nonintegrable fractional analog of the Kelvin model. The article oulines the appropriateness of using the viscoelastic deformation models with the operators of fractional integro-differentiation, based on the comparison of calculations of the considered models with the calculations of the viscoelasticity models having integral operators of integro-differentiation.

Full Text

Введение В работе [1] на основе метода структурного моделирования и гипотезы Больцмана-Вольтерры о наследственно-упругом деформируемом твердом теле рассмотрены линейные дробные аналоги некоторых классических реологических моделей: - модель Скотт Блэра как аналог и обобщение модели Ньютона для вязкой жидкости , (1) а также дробные аналоги моделей: - Фойхта , (2) - Максвелла , (3) - Кельвина , (4) - Зенера , (5) где - левосторонняя дробная производная Римана-Лиувилля порядка функции ; - интеграл Римана-Лиувилля порядка [2]; E1, E2, η - постоянные величины, которые при совпадают с модулями упругих элементов и коэффициентом демпфирования соответственно в классических моделях. В указанной работе было отмечено, что если в исходном определяющем соотношении В. Вольтерры использовано ядро Абеля, то в определяющих соотношениях в дифференциальной форме возникнут дробные производные Римана- Лиувилля на отрезке временной оси. Показана корректность классической задачи Коши для дробных дифференциальных уравнений относительно некоторых линейных комбинаций функций напряжения и деформации. Найдены явные решения задачи ползучести при постоянном напряжении на стадиях нагружения и разгрузки. Отмечено, что они могут быть использованы для аппроксимации экспериментальных зависимостей деформации от времени, полученных в опытах по одноосному растяжению деформируемых тел под действием постоянных нагрузок. Решение задачи аппроксимации позволяет идентифицировать параметры реологической модели и тем самым определить вид и структуру математической модели вязкоупругого поведения изучаемой среды в форме определяющего соотношения между функциями напряжения и деформации. В работе [3] решена задача идентификации параметров линейных математических моделей для перечисленных выше механических моделей вязкоупругого деформирования в упрощенном варианте и приведены оценки погрешностей соответствующих аппроксимаций. В работах [1] и [3] были использованы данные эксперимента по нагружению образцов из поливинилхлоридного пластиката, приведенные в работе [4]. Расчеты [3] показали существенный разброс значений параметров математических моделей. Поэтому выполнялось осреднение их значений по всем реализациям. Математическая модель с усредненными параметрами затем предлагалась в качестве теоретической модели вязкоупругого поведения испытываемых образцов. При этом относительная погрешность аппроксимаций заключалась в пределах от 4 до 19 %, что хорошо видно из расчетных данных, приведенных в работе [3]. Как отмечено в работе [5], установленная в [3] зависимость идентифицируемых параметров модели от уровня напряжения в серии экспериментов и в особенности факт увеличения погрешности аппроксимаций экспериментальных точек теоретическими значениями деформации при крайних уровнях напряжения может свидетельствовать о нелинейности механической модели. Построение нелинейных моделей базируется на обобщении линейных. Основные идеи классических (линейных) наследственно-упругих моделей восходит к работам В. Вольтерры [6, 7] и X. Больцмана [8]. Они получили мощное развитие в работе Ю.Н. Работнова [9], которая, по существу, заложила начало школы механики наследственно-упругих твердых тел. Отметим здесь также несколько предшествующих работ [10-15]. К настоящему моменту опубликовано огромное количество работ, посвященных дробным реологическим моделям (см., например, библиографические списки в работах [2, 16-19]). Как отмечено в монографии [2, с. 443], большинство авторов используют в определяющих соотношениях дробные производные Капуто [20-23], производные Грюнвальда-Летникова [24] или лиувиллевскую форму дробного интегро-дифференцирования на всей числовой оси, в то время как дробные реологические моде-ли описываются в терминах дробных производных Римана-Лиувилля на отрезке временной оси [25]. Для нахождения явных решений дифференциальных уравнений дробного порядка, возникающих из определяющих соотношений, используется в основном преобразование Лапласа. В работе [1] показано, что они могут быть найдены собственно методами дробного анализа путём нахождения резольвент интегральных уравнений в терминах некоторых специальных функций. Проблема идентификации параметров дробных реологических моделей в процессе аппроксимации экспериментальных данных определёнными аналитическими зависимостями обсуждалась уже авторами ранних работ по наследственно упругим моделям (см., например, [11]). В монографии Ю.Н. Работнова [26, с. 86] приведен метод прямого определения параметров дробно-экспоненциального ядра ползучести, предложенный ранее авторами работы [27]. Однако указанный метод основан на сравнении образов преобразований Лапласа явного (аналитического) решения задачи с образом полиномиальной аппроксимации экспериментальных точек и минимизации суммы квадратов отклонений методом координатного спуска в пространстве параметров модели. В действительности проблема может быть решена непосредственно путём сопоставления экспериментальных и теоретических значений. В немногочисленных публикациях последних десяти лет проблема идентификации параметров дробных моделей в основном решается на теоретическом уровне, например, методами спектрального анализа [28]. В работе [29] параметры модели определяются исходя из нескольких характерных точек, полученных в эксперименте, путем подстановки значений деформации в аналитические решения соответствующей задачи. В связи с вышеизложенным целью настоящей работы является: 1) разработка методики идентификации параметров линейных дробных моделей (1)-(5) вязкоупругого деформирования поливинилхлоридного пластиката (ПВХП) с использованием методов наименьших квадратов и координатного спуска; 2) построение нелинейных моделей деформирования на основе обобщения линейных соотношений (1)-(5); выделение класса нелинейных математических моделей, для которых решение задачи ползучести можно получить в явном виде, и разработка численных методов решения этой задачи для неинтегрируемых моделей; 3) проверка адекватности линейных и нелинейных дробных моделей экспериментальным данным вязкоупругого деформирования ПВХП, анализ погрешностей моделей; 4) обсуждение целесообразности использования дробных моделей в теории вязкоупругости. 1. Описание эксперимента и идентификация параметров линейных математических моделей вязкоупругого тела В работе [4] исследовано поведение поливинилхлоридных трубок длиной и площадью поперечного сечения 1,2 мм2 при температуре 20 оC и постоянных напряжениях МПа ( ), действующих на образцы в течение 8 ч, после чего производилась полная разгрузка образцов ( при ч). Экспериментальные значения деформации измерялись в определенные дискретные моменты времени. При каждом уровне напряжения испытывалось от 3 до 5 образцов, а затем значения деформации усреднялись. На рисунке значками представлены усреднённые кривые экспериментальных значений деформации на стадиях нагрузки и разгрузки. Ось абсцисс отражает шкалу времени в часах, ось ординат - шкалу деформации. Рис. Экспериментальные (значки), расчетные по линейной модели (4) (штриховые линии) и расчетные по нелинейной модели (29) (сплошные линии) кривые деформации для поливинилхлоридного пластиката при пяти постоянных уровнях напряжений: 1 - s = 4,655 МПа; 2 - s = 6,288 МПа; 3 - s = 8,738 МПа; 4 - s = 10,372 МПа; 5 - s = 12,005 МПа с последующей разгрузкой Fig. Experimental (symbols), calculated by the linear model (4) (dashed line) and calculated by the nonlinear model (29) (solid line) stress curves for polyvinyl chloride plasticate at five different constant stress levels: 1 - s = 4.655 MPa; 2 - s = 6.288 MPa; 3 - s = 8.738 MPa; 4 - s = 10.372 MPa; 5 - s = 12.005 MPa with subsequent unloading Данный материал имеет мгновенную упругую деформацию при приложении и снятии нагрузки, а под действием постоянной нагрузки проявляет свойство вязкоупругости. В настоящей работе предлагается определять параметры линейных (и далее нелинейных) математических моделей путём аппроксимаций экспериментальных данных во временных точках (i = , - число измерений) по всем уровням напряжения ( - количество уровней напряжения, H(t) - функция Хевисайда) на основе метода наименьших квадратов с последующим уточнением значений этих параметров методом координатного спуска, используя аналитические или численные решения определяющих уравнений соответствующих моделей. Найденные в работе [1] аналитические зависимости деформации как функции времени и напряжения для линейных математических моделей позволяют вычислить теоретическую деформацию в тех же дискретных по времени и напряжению точках: , где - идентифицируемые параметры. Искомые параметры определяются путем минимизации функционала (6) где - экспериментальные значения деформации в дискретные моменты времени ti при указанных выше уровнях напряжений . Рассмотрим решение задачи идентификации параметров математической модели на примере использования линейного дробного аналога модели Кельвина (4). При изменении напряжения по закону (7) где ; Т - момент времени разгрузки, решение задачи ползучести для модели (4) имеет вид [1] (8) (9) - обобщенная (двухпараметрическая) дробная экспоненциальная функция [30], Eα(z,μ) - функция типа Миттаг-Леффлера [31]. Будем аппроксимировать экспериментальные данные теоретическими значениями деформации на стадии нагружения в дискретные моменты времени ti при заданных значениях напряжения по формуле (10) Дробный аналог модели Кельвина содержит четыре параметра (α, E1, E2, η), подлежащие идентификации. В начальный момент времени t1 = 0 полная деформация принимает значения, равные величинам мгновенной упругой деформации: . (11) Используя экспериментальные значения упругой деформации при каждом значении напряжения из (11) находим величину , затем усреднением получаем оценку этой величины. Остальные неизвестные параметры , E1, η предварительно находим в первом приближении. Раскладывая дробно-экспоненциальную функцию в ряд и ограничиваясь первыми двумя членами (12) где гамма-функция Эйлера [2], из (10)-(12) получим приближенное равенство (13) Обозначим , , , , . Логарифмируя (13), получим . Минимизируя далее функцию найдём значения параметров α и δ в первом приближении, используя необходимые условия экстремума Вернувшись к введенным выше обозначениям, получим приближенное значение параметра . Оставшийся параметр Е1 найдем из известной асимптотической формулы для линейного дробного аналога модели Кельвина (4), справедливой для при [1]: , (14) где . Приравнивая правую часть равенства (14) к последнему эксперимен-тальному значению деформации ( ) при различных уровнях нагрузки , получим приближенное равенство , из которого определяется E1 при каждом , а затем усреднением находится оценка величины E1. Полученные значения параметров являются начальными приближениями для дальнейшего их уточнения методом координатного спуска (методом Хука-Дживса [32]), в соответствии с которым осуществлялась вариация параметров в четырёхмерном пространстве вдоль координатных осей в направлениях, минимизирующих функционал (6). В итоге получены следующие значения параметров модели Кельвина (4): (параметр дробности), МПа, МПа, Погрешность аппроксимации экспериментальных данных по всем пяти реализациям (сплошные линии на рис. 1) в норме , (15) где - количество кривых ( ); - количество экспериментальных значений на реализации при , составила 10,34 %. Для анализа полученных результатов в табл. 1 приведены погрешности аппроксимаций в каждой из пяти реализаций модели Кельвина (4). По разработанной схеме идентификации параметров аналогичные расчеты выполнены и для остальных моделей (1)-(3) и (5). Для моделей Скотт Блэра (1) и дробного аналога модели Фойхта (2), которые не описывают мгновенную упругую деформацию, использовались экспериментальные данные без учета упругой деформации. В табл. 2 приведены погрешности всех моделей (1)-(5) для времени ч (т.е. учитывалась деформация и при нагрузке, и при разгрузке), рассчитанные в норме (15). Анализ данных табл. 2 свидетельствует, что в целом все пять линейных моделей дают соизмеримую погрешность, при этом модели Скотт Блэра и Максвелла дают несколько худшие результаты, чем остальные. Второй важный результат следует из анализа данных, представленных на рис. 1 и в табл. 1. Здесь наблюдается систематическое увеличение погрешности при относительно больших и относительно малых значениях напряжений. Это, в свою очередь, свидетельствует о наличии нелинейности между напряжениями и деформациями и ставит проблему разработки нелинейных дробных моделей вязкоупругого деформирования, что и является одной из целей дальнейших исследований. Таблица 1 Погрешности аппроксимаций для линейного дробного аналога модели Кельвина (4) для каждой реализации Table 1 Approximation errors for the linear fractional Kelvin model (4) for each implementation Stresses, MPa 4.655 6.288 8.738 10.372 12.005 Δ, % 19.103 5.406 3.212 11.897 12.101 Таблица 2 Средние погрешности аппроксимаций для различных линейных дробных математических моделей поливинилхлоридного пластиката Table 2 Average errors of approximations for different linear fractional mathematical models of polyvinyl chloride plastic Fractional models of Scott Blair (1) Voigt (2) Maxvell (3) Kelvin (4) Zener (5) 13.975 % 10.901 % 12.842 % 10.343 % 10.329 % 2. Нелинейные математические модели, допускающие аналитическое решение задачи ползучести. Идентификация параметров Существуют различные феноменологические подходы к введению нелинейных зависимостей в определяющее соотношение между напряжением и деформацией. Например, в равенстве , (16) где s(t), e(t) - напряжение и деформация в момент времени t, можно предположить, что непрерывный оператор A является нелинейным. Однако для создания «работоспособной» математической модели вид оператора A требует конкретизации. Как уже отмечалось в разд. 1 в математических моделях вязкоупругого тела c памятью возникают дифференциальные или интегральные операторы Римана-Лиувилля дробного порядка , причем параметр a - один из параметров, подлежащих идентификации по результатам экспериментальных данных. Считая оператор A в равенстве (16) линейным, можно рассматривать зависимости следующего вида: , где - монотонные, достаточно гладкие функции. Вид этих функциональных зависимостей также требует конкретизации. Обычно рассматриваются зависимости степенного характера. Поскольку в работе используется метод структурного моделирования вязкоупругого тела, будем считать, что зависимость между напряжением и деформацией в упругих элементах структурных моделей определяется равенствами вида . (17) Вязкие элементы будем моделировать зависимостями, в которых напряжение пропорционально степени дробной производной Римана-Лиувилля функции деформации (18) или пропорционально дробной производной степени деформации . (19) В приведенных равенствах и - некоторые положительные константы, совпадающие с модулем упругого элемента и коэффициентом демпфирования, если n = m = 1, α = 1. Заметим, что при α = 1 нелинейный вариант элемента Скотт Блэра (18) превращается в хорошо известную нелинейную зависимость напряжения от степени скорости деформации . Нелинейный же вариант элемента Скотт Блэра (19) при α = 1 приводит к закону, в котором напряжение пропорционально скорости деформации с коэффициентом, зависящим от величины деформации , где . Отметим, что соотношения такого типа соответствуют теории упрочнения [37], широко используемой при описании ползучести различных материалов. Выполняя последовательное соединение упругого элемента (17) и вязкого элемента (18), получим первый вариант нелинейного дробного аналога модели Максвелла в дифференциальной форме . (20) Если соединить последовательно тот же упругий элемент и вязкий элемент (19), то получим второй вариант этой же математической модели (21) Нелинейные уравнения (20) и (21) - простейшие дифференциальные уравнения с производной Римана-Лиувилля порядка . Они позволяют найти функцию деформации для любой зависимости . Действительно, с начальным условием уравнение (20) редуцируется к интегральному соотношению , (22) а уравнение (21) - к соотношению , (23) которые можно рассматривать, с одной стороны, как нелинейные дробные аналоги модели Максвелла в интегральной форме, с другой - как представления функции деформации при известном законе нагружения. Вычисляя дробные интегралы для , где , найдем решения задачи ползучести для обоих вариантов нелинейной модели Максвелла: , (24) , (25) где - гамма-функция Эйлера. Выполняя параллельное соединение упругого элемента (17) и вязкого элемента (18), получим первый вариант нелинейного дробного аналога модели Фойхта в дифференциальной форме: . (26) Это определяющее соотношение является существенно нелинейным дифференциальным уравнением с дробной производной Римана-Лиувилля. Найти аналитическое решение подобного дифференциального уравнения не представляется возможным. Вопросы, связанные с методами численного решения нелинейных дифференциальных уравнений с дробными производными и проблемой идентификации параметров такого рода моделей, будут рассмотрены в следующем разделе данной статьи. Если соединить параллельно упругий элемент (17) и вязкий элемент (19), то получим второй вариант нелинейного дробного аналога модели Фойхта: . (27) Если в этом соотношении ограничиться предположением m = n, то уравнение (27) окажется хорошо известным дифференциальным уравнением Барретта относительно функции [33]. Его решение с однородным начальным условием при в терминах дробно-экспоненциальной функции записывается в виде (28) где функция определена соотношением (9). Выполним теперь последовательное соединение второго варианта дробного аналога модели Фойхта (27) с параметрами , описываемого в дифференциальной форме равенством , с упругим элементом, определяемым зависимостью , где - напряжения, а - деформации в соответствующих структурных элементах. Учитывая, что а полная деформация получим определяющее соотношение для нелинейного варианта дробного аналога модели Кельвина в следующем виде: . (29) Обозначая получим вновь дифференциальное уравнение Барретта относительно функции u(t): . (30) С начальным условием уравнение (30) редуцируется к интегральному уравнению Вольтерры второго рода с ядром Абеля . Его решение записывается в терминах резольвентного оператора [1] . Используя свойства этого оператора [34], нетрудно найти решение интегрального уравнения для в виде . Для функции деформации получим выражение (31) где ; . Записывая закон изменения напряжения равенством , найдем зависимость для деформации на стадии нагрузки ( при ) и разгрузки ( при ) в следующем виде: (32) Продемонстрируем методику решения задачи идентификации параметров на примере нелинейного дробного аналога модели Кельвина (29) применительно к уже рассмотренным экспериментальным данным ПВХП. Сохраняя обозначения, принятые в разд. 1, аппроксимируем экспериментальные данные теоретическими значениями деформации в точках , вычисленными по формуле (31): . (33) В начальный момент времени возникает лишь мгновенная упругая деформация . (34) Логарифмируя равенство (34) и учитывая обозначения: , , , , для нахождения параметров γ и минимизируем функционал: . Остальные неизвестные параметры предварительно находим в первом приближении, подобно тому как это делалось для линейных моделей. Для этого возьмем первые два члена ряда дробно-экспоненциальной функции (12). Тогда из равенства (33) следует, что . Прологарифмируем это приближённое равенство. Получим , где , , , . Обозначая и минимизируя функционал на основании метода наименьших квадратов определяем величины и далее . Для нахождения оставшегося параметра заметим, что, как и в линейной модели, функция деформации в задаче ползучести в рамках нелинейного дробного аналога модели Кельвина имеет асимптоту. Из формулы (31) следует существование предела . (35) Приравнивая правую часть равенства (35) к последнему экспериментальному значению деформации на стадии нагрузки ( ) при различных уровнях нагрузки , получим приближенные равенства , где , из которых параметр вычисляется для каждого уровня нагрузки и далее производится усреднение этой величины. Полученные значения параметров для стадии нагружения будут являться начальными приближениями для дальнейшего их уточнения методом координатного спуска, но уже на стадиях нагрузки и разгрузки. Рассчитанные по описанной схеме значения параметров модели Кельвина (29) приведены в последней строке табл. 3. По аналогичной схеме определялись значения остальных нелинейных моделей, допускающих аналитические решения, которые также представлены в табл. 3. Отметим, что для моделей (18), (19) и (27), не описывающих мгновенно-упругую деформацию, использовались лишь экспериментальные значения вязкоупругой деформации В табл. 4 приведены оценки погрешности в соответствии с нормой (15) для всех исследованных нелинейных моделей, допускающих аналитические решения, при этом минимальную погрешность дает модель Кельвина (29). Для сравнения в первой строке приведены погрешности для линейных моделей. Модели Скотт Блэра (18) и (19) при переходе к линейному случаю (n = 1, m = 1) принимают вид уравнения (1), нелинейная модель Фойхта (27) становится моделью (2), нелинейные модели Максвелла (20) и (21) преобразуются в (3), а нелинейная модель Кельвина (29) - в (4). В качестве примера на рис. 1 сплошными линиями нанесены расчетные данные по нелинейной модели (29). Таблица 3 Значения параметров для нелинейных дробных моделей, допускающих аналитические решения задачи вязкоупругости Table 3 The values of parameter for nonlinear fractional models that admit analytical solutions of the viscoelasticity problem Fractional Models of α E1, MPa E2, MPa η, MPa∙c b γ Scott Blair (18) 0.180 - - 126.026 1.298 - Scott Blair (19) 0.173 - - 139.119 1.271 - Voigt (27) 0.274 44.413 - 84.081 1.299 - Maxvell (20) 0.180 - 331.807 127.721 1.291 1.281 Maxvell (21) 0.174 - 331.807 139.342 1.271 1.281 Kelvin (29) 0.284 47.743 331.807 81.328 1.281 1.299 Таблица 4 Погрешности аппроксимаций (Δ, %) для линейных и нелинейных дробных аналогов реологических моделей, допускающих аналитические решения задачи вязкоупругости Table 4 The approximation errors (Δ, %) for linear and nonlinear fractional analogs of the rheological models that admit analytical solutions of the viscoelasticity problem Fractional models of Scott Blair (18) Scott Blair (19) Voigt (27) Maxvell (20) Maxvell (21) Kelvin (29) Linear (n = 1, m = 1) Nonlinear 13.975 9.167 13.975 8.002 10.901 8.007 12.842 7.985 12.842 7.022 10.343 7.014 Анализ данных табл. 4 свидетельствует о том, что средние погрешности аппроксимаций значений деформации ПВХП нелинейными дробными моделями улучшены на 3,3 % - 6 % по сравнению с линейными дробными моделями. Полученный результат может служить обоснованием использования нелинейных моделей с операторами дробного интегро-дифференцирования в теории вязкоупругости. 3. Численное исследование и идентификация параметров нелинейных моделей наследственно-упругого тела дробного порядка, не допускающих аналитического решения Основная цель данного раздела заключается в разработке методики решения задачи идентификации параметров математических моделей одноосного напряженного состояния в случае, когда определяющие соотношения, связывающие напряжение и деформацию вязкоупругой среды, в дифференциальной или интегральной формах не позволяют найти решение в явном виде в терминах известных специальных функций. В этом случае неизвестные параметры математической модели могут быть найдены численно методом координатного спуска с обращением на каждом шаге к численному решению определяющего интегрального уравнения, требуя минимума среднеквадратического отклонения расчетных значений деформации от экспериментальных данных. В данном разделе для идентификации параметров использовались экспериментальные данные ПВХП без учета разгрузки. В связи с этим значения погрешностей аппроксимаций для интегрируемых линейных и нелинейных моделей пересчитаны только для случая нагрузки с целью дальнейшего сравнения с погрешностями аппроксимаций при численном исследовании нелинейных дробных моделей неинтегрируемого типа. При аппроксимации данных серии экспериментов с различными постоянными нагрузками следует обратиться к формуле (15), в которой вместо точных значений деформации должны фигурировать численные значения . Рассмотрим в качестве примера нелинейный дробный аналог модели Фойхта (26). Разрешая это уравнение относительно дробной производной, получим . (36) С однородным начальным условием дифференциальное уравнение (36) эквивалентным образом редуцируется к интегральному уравнению Вольтерры с ядром Абеля в канонической форме Гаммерштейна [35]: . (37) Обозначим в нелинейном интегральном уравнении (37) подынтегральное выражение как функцию u(t). Тогда вычислительная процедура сведется к численному решению интегро-функциональной системы уравнений . (38) Для построения итерационных процедур численного решения нелинейных интегральных уравнений с ядром Абеля воспользуемся идеями, изложенными в работе [36]. Метод основывается на замене входящего в уравнение интегрального оператора квадратурными формулами, в качестве которых используются формулы приближенного вычисления дробных интегралов, аналогичные формуле прямоугольников и формуле трапеций. Для и фиксированного выбирается разбиение отрезка на s равных частей с шагом точками , где t0 = 0, ts = . Дробный интеграл представляется в виде суммы интегралов по отрезкам , а подынтегральная функция - рядом Тейлора в окрестности точки tk. Ограничиваясь двумя первыми членами разложения для функции , в работе [36] авторы привели аналог первой (левой) формулы прямо-угольников , (39) где . Аналог второй (правой) формулы прямоугольников имеет вид , (40) где . Ограничиваясь тремя членами ряда Тейлора в разложении подынтегральной функции , в работе [36] авторы получили аналог формулы трапеций (41) Обозначая точное решение системы уравнений (38) в точке за , - за , а численные значения соответственно за и , и используя, например, аналог левой формулы прямоугольников для дробного интеграла, получим следующие итерационные формулы: , (42) , (43) (44) Отметим, что использование аналога левой формулы прямоугольников предпочтительно для нелинейных интегральных уравнений, так как на каждом шаге вычисляется через значения, найденные на предыдущих шагах. Использование двух других формул требует обращения к процедуре решения нелинейного алгебраического уравнения на каждом шаге. Для численной реализации метода построения для неинтегрируемой нелинейной модели с операторами дробного интегро-дифференцирования разработано соответствующее программное обеспечение. Его тестирование в линейном случае (n = m = 1) с использованием аналога формулы левых прямоугольников с параметрами s = 480 (количество разбиений), α = 0,4, E = 128 МПа, η = 114 МПа∙cα, s0 = 12 МПа, которые были найдены в результате идентификации параметров дробного линейного аналога модели Фойхта (2) по результатам эксперимента на растяжение ПВХП без учета мгновенной упругой деформации, показало, что среднеквадратическое отклонение численных решений от аналитического решения составляет 0,34 %, что вполне удовлетворительно для наших целей. Использование формул (40) и (41) приводит ожидаемо к более точным результатам: 0,17 % - для аналога правых прямоугольников и 0,05 % - для аналога формулы трапеций. Тестирование программного обеспечения на основе известного решения нелинейной интегрируемой модели Фойхта (27) показало практически такие же результаты, как для описанной выше линейной модели. В дальнейшем авторы использовали в расчетах аналог левой формулы прямоугольников как более алгоритмически удобную. Для того чтобы сравнить качество аппроксимаций экспериментальных данных испытываемых поливинилхлоридных трубок теоретическими значениями деформации, определяемыми из явных решений задачи ползучести в рамках линейных и нелинейных интегрируемых математических моделей, со значениями в неинтегрируемом случае, рассмотрим второй вариант нелинейного дробного аналога неинтегрируемой модели Кельвина, приведя предварительно вывод его основного определяющего соотношения. Выполним последовательное соединение первого варианта дробного аналога модели Фойхта (26) с параметрами , описываемого в дифференциальной форме равенством , с упругим элементом, определяемым зависимостью , где - напряжения, а - деформации в соответствующих структурных элементах. Учитывая, что , а полная деформация , получим определяющее соотношение для этого нелинейного варианта дробного аналога модели Кельвина в следующем виде: . (45) Обозначая , (46) из (45) получим нелинейное дифференциальное уравнение относительно функции v(t): . Разрешим это уравнение относительно дробной производной: . (47) С начальным условием уравнение (47) редуцируется к интегральному уравнению . (48) Обозначим . (49) Тогда с учетом формулы (49) интегральное уравнение (48) приводится к виду . (50) Таким образом, задача сведена к необходимости численного решения интегро-функциональной системы уравнений (46), (49), (50). Теперь продемонстрируем методику решения задачи идентификации параметров для нелинейного дробного аналога модели Кельвина (45), используя эти же экспериментальные данные по растяжению образцов из ПВХП. Параметры E2 и n2, необходимые для расчета значений мгновенной упругой деформации, найдем аналогично методике, изложенной для модели (29), учитывая, что . Остальные параметры: α, η, E1, n1 и m найдем с помощью метода координатного спуска, минимизируя функцию (15), в которой точные значения деформации заменены на численные . Численные значения деформации будем искать, используя аналог левой формулы прямоугольников для дробного интеграла (40), по следующим итерационным формулам: (51) (52) где . , а . (53) В качестве начальных параметров α, η, E1 используются значения, полученные при идентификации параметров в серии экспериментов для линейной модели (4), поэтому на первом этапе полагаем n1 = 1, m = 1. По результатам идентификации параметров модели (45) получены следующие их значения: α = 0,434, E1 = 104,472 МПа, E2 = 331,807 МПа, η = 123,045 МПа∙с, n1 = 0,882, n2 = 0,781, m = 1. В табл. 5 представлены погрешности аппроксимаций для всех трёх рассмотренных дробных аналогов модели Кельвина в линейном и нелинейном случаях в норме (15), при этом погрешность Δ для неинтегрируемого случая составила 6,379 %. Анализ данных табл. 5 свидетельствует о том, что относительная погрешность аппроксимаций экспериментальных данных аналитическим решением задачи ползучести на основе нелинейного дробного аналога модели Кельвина (29) существенно меньше, чем в случаях, когда используется линейная модель (4) или численные методы решения определяющего интегрального уравнения в нелинейном неинтегрируемом дробном аналоге модели Кельвина (45). Но обе нелинейные модели дают лучшие результаты по сравнению с линейной моделью. Таблица 5 Средние погрешности аппроксимаций (Δ, %) для различных вариантов дробного аналога модели Кельвина Table 5 The average approximation errors (Δ, %) for different variants of the fractional Kelvin model Linear model (4) Nonlinear model (29) Nonlinear model (45) 9.872 % (using analytical solution (8)) 4.256 % (using analytical solution (31)) 6.379 % (using numerical solution (51)-(53)) 4. Анализ математических моделей вязкоупругих сред с операторами дробного и целочисленного интегро-дифференцирования Вопросы целесообразности применения моделей с операторами дробного интегро-дифференцирования в теории вязкоупругости имеют долгую историю (см., например, [9]). В настоящей работе авторы выполнили сравнительный анализ решений на основе дробных линейной (4) и нелинейных (26), (27) моделей Фойхта с данными, полученными на основе классической модели теории вязкоупругости с операторами целочисленного интегрирования Ю.П. Самарина [38, 39]. В цитируемых работах автором представлена нелинейная модель вязкоупругого тела в виде , (54) где и - деформация и напряжение в момент времени t; - нормирующий коэффициент; Е - мгновенный модуль упругости; - некоторые подлежащие идентификации параметры модели. Отбрасывая мгновенно-упругую деформацию при , из (54) получим выражение для чистой вязкоупругой компоненты деформации : . (55) Посредством экспериментальных данных для серии кривых ползучести при пяти постоянных уровнях напряжения поливинилхлоридного пластиката (см. рис. 1) и методики [38, 39] получена аппроксимация для рассматриваемого материала со следующими параметрами: (1/ч), (1/ч), (1/ч), т.е. имеем в (55) 3 экспоненциальных слагаемых и 7 параметров. Приведем оценки погрешностей в норме (15) анализируемых моделей: для модели (55) (семь параметров) - Δ = 9,23 %, линейного дробного аналога модели Фойхта (2) (три параметра) - Δ = 10,47 %, нелинейного неинтегрируемого дробного аналога модели Фойхта (26) (пять параметров) - Δ = 4,41 %, нелинейного интегрируемого дробного аналога модели Фойхта (27) (четыре параметра) - Δ = 4,27 %. Отметим, что погрешности определены только по кривым ползучести при нагрузке ( ч). Таким образом, при меньшем числе параметров линейная дробная модель дает практически одинаковые результаты с моделью (55), а нелинейная дробная модель всего при четырех параметрах дает вдвое меньшую погрешность, чем модель с целочисленным интегральным оператором с семью параметрами. Отсюда можно сделать вывод о целесообразности использования математических моделей с дробными операторами интегро-дифференцирования в теории вязкоупругости. Выводы Выполненные исследования позволяют сформулировать следующие результаты. 1. Разработан ряд линейных и нелинейных математических моделей с операторами дробного интегро-дифференцирования для описания реологического деформирования вязкоупругих сред. 2. Разработана методика идентификации параметров дробных моделей вязко-упругого деформирования, реализованная применительно к экспериментальным данным по ползучести поливинилхлоридного пластиката. 3. Выполнен сравнительный анализ погрешностей аппроксимации экспериментальных данных расчётными кривыми вязкоупругого деформирования при постоянных напряжениях, полученных по линейным и нелинейным дробным моделям. 4. Выполнен сравнительный анализ данных расчетов по линейным и нелинейным дробным моделям с данными расчетов по определяющим соотношениям вязкоупругости с целочисленными операторами интегрирования применительно к поливинилхлоридному пластикату. Показано, что при меньшем числе параметров линейной дробной модели погрешность аппроксимации практически совпадает, а для нелинейной дробной модели погрешность в два раза меньше погрешности по модели с классическим интегральным представлением. Этот факт может служить обоснованием использования моделей с операторами интегро-дифференцирования дробного порядка в теории вязкоупругости.

About the authors

E N Ogorodnikov

Samara State Technical University

V P Radchenko

Samara State Technical University

L G Ungarova

Samara State Technical University

References

  1. Огородников Е.Н., Радченко В.П., Унгарова Л.Г. Математическое моделирование наследственно-упругого деформируемого тела на основе структурных моделей и аппарата дробного интегро-дифференцирования Римана-Лиувилля // Вестн. Самар. гос. техн. ун-та. Сер. Физ.-мат. науки, - 2016. - № 1(20). - С. 167-194. doi: 10.14498/vsgtu1456
  2. Kilbas A.A., Srivastava H.M., Trujillo J.J. Theory and Applications of Fractional Differential Equations / North-Holland Mathematics Studies. Vol. 204. - Amsterdam: Elsevier, 2006. - 523 p.
  3. Унгарова Л.Г. Применение линейных дробных аналогов реологических моделей в задаче аппроксимации экспериментальных данных по растяжению поливинилхлоридного пластиката // Вестн. Самар. гос. техн. ун-та. Сер. Физ.-мат. науки. - 2016. - № 4(20). - С. 691-706. doi: 10.14498/vsgtu1523
  4. Радченко В.П., Голудин Е.П. Феноменологическая стохастическая модель изотермической ползучести поливинилхлоридного пластиката // Вестн. Самар. гос. техн. ун-та. Сер. Физ.-мат. науки. - 2008. - № 1(16). - С. 45-52. doi: 10.14498/vsgtu571
  5. Огородников Е.Н., Унгарова Л.Г. Аналитические решения задачи о ползучести и идентификация параметров нелинейных математических моделей наследственно-упругого тела // Тр. Х Всерос. науч. конф. по механике деформируемого твердого тела. - Самара: Изд-во СамГТУ, 2017. - С. 120-123.
  6. Volterra V. Sulle equazioni integro-differeziali della teoria dell'elascita // Rend. Acс. Naz. Lincei. - 1909. - Vol. 5. - P. 295-301.
  7. Вольтерра В. Теория функционалов, интегральных и интегро-дифференциальных уравнений / пер. с англ. П.Н. Кузнецова. - М.: Наука, 1982. - С. 304.
  8. Boltzmann L. Theorie der elastischen Nachwirkung (Theory of elastic after effects) // Wien. Ber. - 1874. - Vol. 70. - P. 275-306; Boltzman L. Zur Theorie der elastischen Nachwirkung (On the elastic after effect) // Pogg. Ann. (2). - 1878. - Vol. 5. - P. 430-432; Boltzman L. Zur Theorie der elastischen Nachwirkung / Wissenschaftliche Abhandlungen. vol. 2 / Cambridge Library Collection, ed. F. Hasenӧhrl. - Cambridge: Cambridge University Press. - 2012. - P. 318-320. doi: 10.1017/CBO9781139381437.015
  9. Работнов Ю.Н. Равновесие упругой среды с последействием // ПММ. 1948. - Т. 12, № 1. - С. 53-62.
  10. Duffing G. Elastizität und Reibung beim Riementrieb (Elasticity and friction of the belt drive) // Forschung auf dem Gebiet des Ingenieurwesens A. - 1931. - Vol. 2. - No. 3. - P. 99-104. doi: 10.1007/BF02578795
  11. Gemant A. A Method of Analyzing Experimental Results Obtained from Elasto-Viscous Bodies // J. Appl. Phys. - 1936. - P. 311-317. doi: 10.1063/1.1745400
  12. Gemant A. On fractional differentials // Philos. Mag. VII. Ser. - 1938. - Vol. 25. - P. 540-549.
  13. Бронский А.П. Явление последействия в твёрдом теле // ПММ. - 1941. - Т. 5, № 1. - С. 31-56.
  14. Слонимский Г.Л. О законах деформации реальных материалов // ЖТФ. - 1939. - Т. 9, № 20. - С. 1791-1799.
  15. Герасимов А.Н. Обобщение лилейных законов деформирования и его применение к задачам внутреннего трения // ПММ. - 1948. - Т. 12, № 3. - С. 251-260.
  16. Булгаков И.И. Ползучесть полимерных материалов. - М.: Наука. Гл. ред. физ.-мат. лит., 1973. - 288 c.
  17. Podlubny I. Fractional Differential Equations. An introduction to fractional derivatives, fractional differential equations, to methods of their solution and some of their applications // Mathematics in Science and Engineering. Vol. 198. - San Diego: Academic Press, 1999. - 340 p.
  18. Учайкин В.В. Метод дробных производных. - Ульяновск: Артишок, 2008. - 512 с.
  19. Mainandi F. Fractional Calculus and Waves in Linear Viscoelasticity: An Introduction to Mathematical Models. - London: Imperial College Press, 2010. - 347 p. doi: 10.1142/9781848163300
  20. Caputo M., Mainardi F. A new dissipation model based on memory mechanism // Pure and Applied Geophysics. - 1971. - Vol. 91. - No. 1. - P. 134-147. doi: 10.1007/bf00879562
  21. Caputo M., Mainardi F. Linear models of dissipation in anelastic solids // La Rivista del Nuovo Cimento. - 1971. - P. 161-198. doi: 10.1007/bf02820620
  22. Bagley R.L., Torvik P.J. A theoretical basis for the Application of Fractional Calculus to Viscoelasticity // J. Rheol. - 1983. - Vol. 27. - No. 3. - P. 201-210. doi: 10.1122/1.549724
  23. Bagley R.L., Torvik P.J. Fractional Calculus - A Different Approach to the Analysis of Viscoelastically Damped Structures // AIAA. 1983. - Vol. 21. - No. 5. - P. 741-748. doi: 10.2514/3.8142
  24. Schmidt A., Gaul L. Parameter Identification and FE Implementation of a Viscoelastic Con-stitutive Equation Using Fractional Derivatives // Proc. Appl. Math. Mech. - 2002. - Vol. 1(1). - P. 153-154. doi: 10.1002/1617-7061(200203)1:1<153::AID-PAMM153>3.0.CO;2-J
  25. Lewandowski R., Chorążyczewski B. Identification of the parameters of the Kelvin Voigt and the Maxwell fractional models, used to modeling of viscoelastic dampers // Computers and Structures. - 2009. - Vol. 88. - No. 1-2. - P. 1-17. doi: 10.1016/j.compstruc.2009.09.001
  26. Работнов Ю.Н. Элементы наследственной механики твёрдых тел. - М.: Наука, 1977. - 384 c.
  27. Определение характеристик ползучести линейных упруго-наследственных материалов с использованием ЭЦВМ / Е.Н. Звонов, Н.И. Малинин, Л.X. Паперник, Б.М. Цейтлин // Изв. АН СССР. МТТ. - 1968. - № 5. - С. 76-85.
  28. Vasques С.М.A., Dias Rodrigues J., Moreira R.A.S. Experimental identification of GHM and ADF parameters for viscoelastic damping modeling // III European Conference on Computational Mechanics. - Springer, 2006. - P. 76-95. doi: 10.1007/1-4020-5370-3_173
  29. Параметрическая идентификация математической модели вязкоупругих материалов с использованием производных дробного порядка / С.В. Ерохин, Т.С. Алероев, Л.Ю. Фриштер, А.В. Колесниченко // International Journal for Computational Civil and Structural Engineering. - 2015. - Т. 11, № 3. - С. 82-86.
  30. Огородников Е.Н., Яшагин Н.С. Некоторые специальные функции в решении задачи Коши для одного дробного осцилляционного уравнения // Вестн. Самар. гос. техн. ун-та. Сер. Физ.-мат. науки. - 2009. - № 1(18). - С. 276-279. doi: 10.14498/vsgtu685
  31. Джрбашян М.М. Интегральные преобразования и представления функций в комплексной области. - М.: Наука, 1966. - 672 с.
  32. Hooke R., Jeeves Т.A. “Direct Search” Solution of Numerical and Statistical Problems // Journal of the ACM (JACM). - 1961. - P. 212-229. doi: 10.1145/321062.321069
  33. Barrett J.H. Differential equations of non-integer order // Canad. J. Math. - 1954. - Vol. 6. - P. 529-541. doi: 10.4153/cjm-1954-058-2
  34. Огородников E.H., Радченко В.П., Яшагин Н.С. Реологические модели вязкоупругого тела с памятью и дифференциальные уравнения дробных осцилляторов // Вестн. Самар. гос. техн. ун-та. Сер. Физ.-мат. науки. - 2011. - № 1(22). - С. 255-268. doi: 10.14498/vsgtu932.
  35. Манжиров А.В., Полянин А.Д. Справочник по интегральным уравнениям: методы решения. - М.: Факториал Пресс, 2000. - 384 с.
  36. Яшагин Н.С. Математическое моделирование и исследование осцилляционных явлений в системах с памятью на основе аппарата дробного интегро-дифференцирования: дис. … канд. физ.-мат. наук. - Самара, 2011. - 186 с.
  37. Работнов Ю.Н. Ползучесть элементов конструкций. - М.: Наука, 1966. - 752 с.
  38. Самарин Ю.П. Уравнения состояния материалов со сложными реологическими свойствами. - Куйбышев: Куйб. гос. ун-та, 1979. - 84 с.
  39. Самарин Ю.П. Построение экспоненциальных аппроксимаций для кривых ползучести методом последовательного выделения экспоненциальных слагаемых // Проблемы прочности. - 1974. - № 9. - С. 24-27.

Statistics

Views

Abstract - 342

PDF (Russian) - 214

Cited-By


PlumX


Copyright (c) 2018 Ogorodnikov E.N., Radchenko V.P., Ungarova L.G.

Creative Commons License
This work is licensed under a Creative Commons Attribution-NonCommercial 4.0 International License.

This website uses cookies

You consent to our cookies if you continue to use our website.

About Cookies