CC BY
24
ISSN 2409-6601. Российский журнал биомеханики. 2017. Т. 21, № 1: 8-17
DOI: 10.15593/RZhBiomeh/2017.1.01 УДК 539.3
К ОЦЕНКЕ ДЕФОРМАТИВНОСТИ РЕШЕТЧАТОЙ ПЛАСТИНКИ ГЛАЗА
А.О. Ватульян1,2, О.А. Потетюнко1
1 Южный федеральный университет, Институт математики, механики и компьютерных наук, Россия, 344006, Ростов-на-Дону, ул. Мильчакова, 8а, e-mail: [email protected]
2 Южный математический институт Владикавказского научного центра Российской академии наук и Республики Северная Осетия-Алания, Россия, 362027, Владикавказ, ул. Маркуса, 22
Аннотация. Рассмотрено несколько задач о деформировании круглой упругой пластинки переменной жесткости с различными граничными условиями, в том числе и при наличии упругого опирания, характеризующегося двумя коэффициентами упругости. На основе вариационного принципа Лагранжа для пластинки был сформулирован и решен численно с помощью метода Ритца ряд вспомогательных задач, не содержащих искомых коэффициентов; получено представление прогиба как дробно-рациональной функции от этих коэффициентов. На основе измерения прогиба в нескольких точках решена обратная задача их отыскания из системы нелинейных алгебраических уравнений. Также была поставлена и решена задача о восстановлении трех параметров: коэффициентов упругости заделки и распределенной нагрузки, связанной с внутриглазным давлением. Проведен ряд вычислительных экспериментов, показавших хорошую точность реконструкции параметров. Выяснено, как характер изменения жесткости пластинки влияет на точки ее перегиба. В рамках исследования также была рассмотрена задача об оценке влияния этих коэффициентов на реконструкцию давления; установлено, что неучет упругости заделки ведет к существенным ошибкам в определении амплитуды нагрузки.
Ключевые слова: решетчатая пластинка, склера, неоднородность, упругое закрепление, метод Ритца, реконструкция.
Введение
Задачи биомеханики глаза в последние годы весьма интенсивно использовались в практической офтальмологии. Среди важных проблем биомеханики глаза отметим задачи, связанные с оценкой факторов, влияющих на внутриглазное давление, в частности, на деформирование решетчатой пластинки, и факторов, его определяющих. Как известно, задний отдел склеры в месте выхода зрительного нерва очень тонок и имеет множество отверстий в области, занятой решетчатой пластинкой. Изучение ее характеристик имеет существенное значение в патогенезе глаукомного повреждения зрительного нерва, поскольку ряд авторов устанавливает анатомические отличия между решетчатой пластинкой здорового и больного глаз [11, 14]. Одной из первых публикаций, посвященных решетчатой пластинке, является работа [10]. В ней решетчатая пластинка моделировалась однородной изотропной пластинкой, где была реализована попытка учета влияния растягивающих усилий, действующих со стороны склеры на прогиб пластинки. В современных работах принято моделировать объект неоднородной пластинкой. Учет неоднородности решетчатой пластинки принципиально
© Ватульян А.О., Потетюнко О. А., 2017
Ватульян Александр Ованесович, д.ф.-м.н., профессор, заведующий кафедрой теории упругости, Ростов-на-Дону
Потетюнко Ольга Андреевна, студентка кафедры теории упругости, Ростов-на-Дону
важен, поскольку она оказывает существенное влияние на ее деформирование при повышении внутриглазного давления [6]. Это и другие свойства решетчатой пластинки и их влияние на деформирование наиболее полно изучены в работах [4, 9]. Однако следует отметить, что модели однослойных пластинок не позволяют учесть некоторые экспериментальные данные; например, тот факт, что атрофия зрительного нерва при повышении внутриглазного давления возникает на «наружном слое» [1].
Немаловажной задачей является поиск точки перегиба решетчатой пластинки [5]. Это связано с тем, что она состоит из нескольких листов соединительной ткани и наиболее сильные деформации претерпевает последний слой. В работе [15] отмечено, что при деформации пластинки возникает точка перегиба, характеризующая смену участков выпуклости и вогнутости, исследование местоположения которой играет важную роль в практической диагностике глаукомы; в этой же работе предложен экспериментальный метод нахождения координаты точки перегиба.
Как отмечалось выше, решетчатая пластинка представляет собой задний отдел склеры, ослабленный множеством отверстий. По различным данным [7, 12], число отверстий около 700 и занимают они примерно 2/3 площади всей пластинки. От центра к краю решетчатая пластинка становится более разреженной, т.е. эффективный модуль упругости в плоскости изотропии убывает в этом направлении, и при моделировании переменности модуля упругости в [2] использован убывающий экспоненциальный закон.
В настоящий момент при исследованиях деформирования решетчатой пластинки принято моделировать ее как круглую или же близкую к круглой пластинку переменной жесткости с жестко заделанным краем [6]. Вместе с тем отметим, что непосредственный контакт решетчатой пластинки со склерой требует учета упругого закрепления, что позволит более точно оценивать прогиб и в рамках такой модели с большей точностью восстанавливать значение внутриглазного давления. На это указывает работа [13], в которой авторы установили, что деформация пластинки в значительной степени зависит от модулей упругости и толщины склеры. Моделирование упругого опирания пластинки может быть осуществлено различными способами. Так, в [3] была рассмотрена задача о колебаниях неоднородной упруго опертой пластинки с двумя коэффициентами жесткости, в рамках такой модели была решена задача восстановления этих коэффициентов по известному прогибу или по двум резонансным частотам. Подходы, предложенные в этой работе, могут оказаться полезными при совершенствовании моделей деформирования решетчатой пластинки.
В настоящей статье рассмотрено несколько задач о деформировании круглой упругой пластинки переменной жесткости с различными граничными условиями, в том числе упругого опирания. Считается, что пластинка по краю имеет упругую связь, которая характеризуется двумя коэффициентами в граничных условиях. В осесимметричном случае сформулирована краевая задача для пластинки переменной жесткости при наличии упругого опирания. Для этого на основе вариационного принципа Лагранжа для пластинки был сформулирован и решен численно с помощью метода Ритца ряд вспомогательных задач, не содержащих искомых коэффициентов; получено представление прогиба как дробно-рациональной функции, содержащей искомые параметры, построена система нелинейных алгебраических уравнений, из которой и находятся коэффициенты жесткости заделки. В рамках такой модели также была поставлена и решена задача о восстановлении трех параметров коэффициентов жесткости и распределенной нагрузки, связанной с внутриглазным давлением. Проведен ряд вычислительных экспериментов, показавший хорошую точность реконструкции параметров. Выяснено, как распределение жесткости пластинки влияет на точки ее перегиба. В рамках исследования также была рассмотрена задача об оценке влияния коэффициентов жесткости на реконструкцию давления; установлено, что неучет упругости заделки ведет к существенным ошибкам реконструкции в его определении.
Постановка задачи
Рассмотрим равновесие решетчатой пластинки, которая моделируется круглой упругой изотропной пластинкой радиусом а переменной жесткости, на нее действует распределенная нагрузка ц. Будем считать, что край пластинки упруго оперт, моделируется двумя упругими связями на краю. Безразмерный лагранжиан для такой пластинки в случае установившихся колебаний построен в [3], а в рассматриваемом статическом случае соответствующий функционал Лагранжа имеет вид
F » = 2 J f (r ) {[ w "(r )]2
" w '(r) " 2 + 2v " w'(r) w"(r) "
_ r _ _ r _
rdr - J qw(r )rdr-
+g w2(1) + g2( w'(1))2
(1)
где введены безразмерные параметры и переменные по формулам
Э(аг) Э0-1 = / (г), Э0 = тах Э( %), w = аЖ,
= С1а2Э01, #2 = С2Э-1, ц(г) = Цо(аг)аг = %а\
Э(%) - цилиндрическая жесткость; w = w(%) - функция поперечного прогиба пластинки; а - радиус пластинки, 0 < % < а; /(г) - безразмерная жесткость; w(г) - прогиб; С1, С2- коэффициенты жесткости заделки; V - коэффициент Пуассона.
Из (1) путем варьирования функционала и приравнивания к нулю коэффициентов при независимых вариациях можно получить уравнение равновесия
Lw = (>" г)(/г')'+ V [(>')''- (^ ")'] = цг (2)
и соответствующие граничные условия
M{w = {-(fw'' r)'+ fr 1w'- v( fw')'+ vfw"+ gxw} = M^w + g1w
= 0,
M2w = {fw'' r + vfw'+ g2 w'}| = M20 w + g2w
= 0.
(3)
Отметим, что граничные условия содержат два параметра, характеризующие влияние склеры на деформирование решетчатой пластинки.
Для нахождения прогиба w использован метод Ритца. Функция прогиба представляется в виде линейной комбинации
w(r) = Z ck Ф*(r)'
k=1
где ск - некоторые коэффициенты; ф^ - базисные функции вида
Ф*(r) = r
_ „2(k-1)
k = 1,2...N.
(4)
(5)
После подстановки (4) в функционал (1) и нахождения его стационарного значения построена система линейных алгебраических уравнений относительно неизвестных коэффициентов разложения ск, решая которую находим функцию прогиба согласно (4).
Таблица 1
Положение точки перегиба в зависимости от закона изменения жесткости
f (r) r
1 0,5799
1 - 0,5r2 0,6399
e -r 0,6499
1 + cos r 0,5999
e-r2 0,6599
Как отмечено во введении, важной задачей в диагностике глаукомы является нахождение точек, в которых происходит изменение характера выпуклости, - точек перегиба. Для их нахождения необходимо выяснить структуру прогиба согласно (4) и далее найти точки перегиба. В табл. 1 показано, как искомая точка перегиба г* движется в зависимости от задания закона
изменения жесткости /(г). Серия расчетов показала, что достаточно выбрать N = 8, поскольку дальнейшее увеличение вносит несущественные поправки и усложняет расчеты. Для простоты принято q = 1. Коэффициент Пуассона считается равным 0,4.
В настоящее время считается, что экспоненциальный закон наиболее адекватно описывает реальную жесткость решетчатой пластинки.
Для того чтобы убедиться в достоверности результатов расчетов прогиба пластинки, сравним численные значения прогибов для однородной пластинки с известными аналитическими результатами [8] в случае жесткой заделки (см. табл. 2).
В табл. 3 приведены прогибы в наборе точек для различных законов жесткости. Принималось, что g1 = g2 = 5 -104. Видно, что наибольшие прогибы реализуются при экспоненциальном распределении жесткости.
Таблица 2
Сравнение прогибов при жестком закреплении
r Точное решение [15] Метод Ритца
g = g2 = 104 g1 = g2 = 5-104 g1 = g2 = 105
0 0,015625 0,015681 0,015636 0,015631
0,25 0,013733 0,013789 0,013744 0,013738
0,5 0,008789 0,008844 0,008799 0,008795
0,75 0,002991 0,003043 0,003001 0,002996
0,95 0,000149 0,000199 0,000159 0,000153
Таблица 3
Прогибы для различных законов жесткости
r f (r) = 1 f (r) = 1 - 0,5r 2 f(r)=
0 0,015636 0,021941 0,029942
0,25 0,013744 0,019285 0,026318
0,5 0,008799 0,012346 0,016847
0,75 0,003001 0,004208 0,005739
0,95 0,000159 0,000219 0,000294
В табл. 4 приведено сравнение прогибов в центре однородной пластинки (I(г) = 1), вычисленных с помощью метод Ритца, с данными, имеющимися в литературе [4], причем в [4] заделка считалась жесткой.
Если оценить влияние параметров й 1 и g2 на деформативность пластины, то оказывается, что значительно более выражено влияние gх, чем g2.
Таблица 4
Сравнение прогибов однородной пластины
Параметры заделки Давление д, мм рт. ст. Прогиб, 10-2 мм Прогиб [5], 10-2 мм Расхождение, %
& = g2 = Ю5 15 0,8904 0,86 3,53
& = g2 = Ю4 15 0,8954 0,86 4,12
& = g2 = 103 15 0,9461 0,86 10,01
Й = 103, g2 = 10-2 15 0,9465 0,86 10,05
й = 5-104, Й2 =7-104 30 1,7819 1,73 2,99
Й = 8-103, Й2 =2-102 30 1,7937 1,73 3,68
Й = 9-102, Й2 =2-104 30 1,9046 1,73 10,09
Й = 9-105, Й2 =3-104 40 2,3730 2,30 3,17
Й = Й2 = 103 40 2,5228 2,30 9,69
й = Й2 = 105 40 2,3743 2,30 3,23
Задача о реконструкции параметров жесткости
Поставим обратную задачу об определении параметров и #2 по известному
прогибу w(гk) = ^к, к = 1, 2...т. Стандартный способ определения параметров состоит в
использовании метода наименьших квадратов. В то же время отметим, что возможен иной подход к решению этой проблемы, для чего сформулируем три вспомогательные задачи:
10 Lw 0 = Щ , Wо |г=1 = ^ W0 ' |г=1 =0;
' =1; (6)
20 Lw1 = 0, M0w1 \ = 0 , M20w1
r=1
30 Lw2 = 0, М0 wJ = 1, wJ = 0.
2 2 2 |г=1 2 1г=1
Оператор L определен при выводе уравнения равновесия, а М(°, М2 -операторы из краевых условий (3), не содержащие и #2. Прогиб будем искать в виде линейной комбинации описанных выше вспомогательных задач
w = ^ + С1W1 + С2W2 . (7) В динамическом случае в [3] установлено, что прогиб в точках г есть дробно-рациональная функция от параметров и #2 вида
^ ) = а0#1 #2 + а1#1 + а2#2 + а3
г а #1 #2 + а, #1 + а2 #2 +1 '
где а0, а1, а2, а0, а[, а2, а3 - известные числовые коэффициенты. В статическом случае прогиб будет иметь аналогичный вид. Таким образом, для решения задачи об определении параметров достаточно знать прогиб в двух точках пластины, тогда для определения параметров имеет место система
К #1 #2 + К2 #1 + ^2 + К3 = ^ (8)
К2 #1 # 2 + К2 #1 + К2 # 2 + К2 = 0.
Рис. 1. Область параметров, обеспечивающих погрешность реконструкции # 2 менее 10 %
(выше сплошной линии)
Таблица 5
Коэффициенты из разложения (4)
10 20 30
N = 3 N = 5 N = 7 N = 3 N = 5 N = 7 N = 3 N = 5 N = 7
С1 0,0109 0,0110 0,0110 0,3111 0,3142 0,3143 -0,2740 -0,2742 -0,2742
С2 -0,0035 -0,0037 -0,0037 -0,0447 -0,0541 -0,0549 0,0372 0,0399 0,0401
Сз 0,0010 0,0017 0,0018 0,0075 0,0195 0,0223 -0,0075 -0,0142 -0,0149
С4 -0,0008 -0,0010 -0,0066 -0,0111 0,0054 0,0074
С5 0,0002 0,0006 0,0012 0,0051 -0,0012 -0,0037
Сб -0,0002 -0,0017 0,0014
С7 0,0001 0,0003 -0,0003
Каждое из уравнений системы (8) в пространстве параметров задает ветви гипербол с осями, параллельными осям координат, и такая система, вообще говоря, имеет два решения. Проблема выбора единственного решения разрешается дополнительным требованием положительности искомых параметров. Вычислительные эксперименты показали хорошую степень точности реконструкции при значениях параметров, показанных на рис. 1, причем при восстановлении й2 погрешность на порядок выше, чем при восстановлении g1. Это объясняется тем, что g 2 вносит меньший вклад в значение прогиба (что можно увидеть также из табл.4).
В табл. 5 приведем значения коэффициентов ск из разложения (4) для
вспомогательных задач (6).
С ростом N коэффициенты убывают, что свидетельствует о корректности вспомогательных задач и устойчивости вычислительной схемы.
Задача о реконструкции трех параметров
В рамках предложенной модели также можно поставить задачу об определении не только параметров g1, g2, но и нагрузки q, моделирующей уровень внутриглазного давления. При решении такой задачи необходимо знать прогиб в трех точках ™(гк ^ = ^к, к = 1,2,3. Тогда на основе представления прогиба будем иметь систему трех нелинейных уравнений вида
q(b1 й g 2+¿11 Й + Ь1 Й 2+Ь1) = Ь1,
q(bo2 Й1 Й2 + Ь12Й1 + Ь22Й2 + ьЗ2) = Ь2, (9)
q(bo3glg 2 + ь13Й1 + ¿23Й 2+ьЗ3) = Ь3.
Видно, что выражения в левых частях пропорциональны q, и (9) приводится к системе вида (8) путем деления первого уравнения на второе и второго - на третье, что позволяет исключить q. Из полученной системы вида (8) находятся й2, аналогично предыдущему пункту, а далее полученные значения подставляются в одно из уравнений (9), откуда и находится q. Сравнение восстановленных данных с исходными показало достаточную степень точности алгоритма.
Исследование влияния упругости заделки
Как отмечалось выше, в большинстве работ решетчатая пластинка считается жестко защемленной по краю. Предлагаемая в настоящей работе модель с упругим опиранием позволяет сформулировать задачу об оценке ошибки при определении давления, которую влечет за собой неучет упругости заделки, т.е. если решетчатая пластина считается жестко защемленной (g1, g2 ^да). Для моделирования такой ситуации сначала заделка считается упругой, фиксируются значения g1, g2, q, находится прогиб ~м(г) при заданных параметрах. Далее заделка полагалась жесткой и по найденному ~м(г) восстанавливалось значение q*. Сравнение изначально заданного q и восстановленного q* позволяет оценить степень влияния коэффициентов g1, g2 на определение давления. В табл. 6 ниже приведены погрешности реконструкции давления в различных точках.
Из таблицы видно, что неучет влияния упругости при небольших значениях параметров g1, g2 приводит к значительной погрешности при реконструкции давления. По мере приближения к краю пластинки точек измерения погрешность увеличивается, так как увеличивается влияние упругости заделки, что также видно из рис. 2.
Таблица 6
Оценка погрешности при восстановлении нагрузки
Заданные параметры Погрешность восстановления q в точке a/8,% Погрешность восстановления q в точке a /2, % Погрешность восстановления q в точке 4a / 5, %
g = g2 = 105 0,04 0,06 0,26
g = g2 = 5-104 0,07 0,12 0,52
g1 = g2 = 104 0,37 0,62 2,58
g1 = g2 = 5-103 0,74 1,24 5,16
g1 = g2 = 103 3,71 6,22 25,80
g1 = g2 = 5-102 7,42 12,44 51,60
8и 50001
4000
3000-
2000 -
1000 -
4500
5000
5500 6000 £2
б
Рис. 2. Область параметров, обеспечивающих погрешность реконструкции давления
точке
2
менее 10 % (выше сплошной линии): a - в
б - в точке a / 2, f (r) = 1 + r ISSN 2409-6601. Российский журнал биомеханики. 2017. Т. 21, № 1: 8-17
4a/5, f (r) = 1 + r2;
а
Заключение
В работе были представлены результаты решения нескольких типов задач. Была решена прямая задача: получены прогибы и найдены точки перегиба пластинки. Также решена обратная задача: проведена реконструкция двух (жесткостей закрепления), трех (жесткостей закрепления и нагрузки) параметров, а также показано, как неучет упругости закрепления приводит к ошибкам реконструкции параметра нагружения.
Список литературы
1. Бауэр С.М., Воронкова Е.Б. Модели теории оболочек и пластин в задачах офтальмологии // Вестник СПбГУ. Серия 1. - 2014. - Т. 59, № 3. - С. 438-458.
2. Бауэр С.М., Зимин Б.А., Товстик П.Е. Простейшие модели теории оболочек и пластин в офтальмологии. - СПб.: Изд-во СПбГУ, 2000. - 92 с.
3. Ватульян А.О., Потетюнко О.А. О колебаниях неоднородной пластины с упруго опертым краем // Известия вузов. Северо-Кавказский регион. Естественные науки. - 2016. - № 2. -С. 28-33.
4. Зинковский П.А., Пальмов В.А. Математические модели и компьютерное моделирование в биомеханике. - СПб.: Изд-во Политехн. ун-та, 2004. - 516 с.
5. Золотухина Л.А. О деформации многослойной решетчатой пластинки диска зрительного нерва // Российский журнал биомеханики. - 2008. - Т. 12, № 4. - С. 40-46.
6. Иомдина Е.Н., Бауэр С.М., Котляр К.Е. Биомеханика глаза: теоретические аспекты и клинические приложения / под ред. В.В. Нероева. - М.: Реал Тайм, 2015. - 208 с.
7. Нестеров А.Г. Глаукома. - М.: Медицина, 1995. - 256 с.
8. Тимошенко С.П., Войновский-Кригер С. Пластинки и оболочки. - М.: Физматгиз, 1963. -635 с.
9. Bauer S.M., Voronkova E.B. On the deformation of the lamina cribrosa under intraocular pressure // Russian Journal of Biomechanics. - 2001. - Vol. 5, № 1. - P. 73-82.
10. Dongqi H., Zeqin R. A biomathematical model for pressure-dependent lamina cribrosa behavior // Journal of Biomechanics. - 1999. - № 6. - P. 579-584.
11. Jonas J.B., Holbach L. Central corneal thickness and thickness of the lamina cribrosa in human eyes // Invest. Ophthalmol. Vis. Sci. - 2005. - Vol. 46. - P. 1275-1279.
12. Quigley H.A., Addicks E.M. Regional differences in the structure of the lamina cribrosa and their relation to glaucomatous optic nerve damage // Arch. Ophthalmol. - 1981. - Vol. 99, № 1. -P. 137-143.
13. Sigal I.A., Flanagan J.G., Ethier C.R. Factors influencing optic nerve head biomechanics. // Invest. Ophthalmol. Vis. Sci. - 2005. - Vol. 46, № 11. - P. 4189-4199.
14. Spoerl E., Goehm A., Pillunat L. The influence of various substances on the biomechanical behavior of lamina cribrosa and peripapillary sclera // Invest. Ophthalmol. Vis. Sci. - 2005. - Vol. 46. -P.1286-1290.
15. Yan D.B., Flanagan J.G., Farra T., Trope G.E., Ethier C.R. Study of regional deformation of the optic nerve head using scanning laser tomography // Current Eye Research. - 1998. - Vol. 17. -P. 903-916.
ESTIMATES OF DEFORMABILITY OF THE SCLERA LAMINA CRIBROSA
A.O. Vatulyan, O.A. Potetyunko (Rostov on Don, Russia)
Several problems on the deformation of a circular elastic plate of variable rigidity with various boundary conditions, including the presence of elastic support, characterized by two coefficients of elasticity, were considered. On the basis of the Lagrange variational principle, a number of auxiliary problems that do not contain the required coefficients was formulated
and solved numerically on the basis of the Ritz method; the representation of the deflection as a fractional function of these coefficients was obtained. On the basis of measuring the deflection at several points, the inverse problem of finding them from the system of nonlinear algebraic equations was solved. The problem of restoring three parameters was also set and solved: the elasticity coefficients of the seal and the distributed load associated with intraocular pressure. A number of computational experiments was performed, which showed a good accuracy of reconstruction of the parameters. It was found out how the nature of the change of the rigidity of the plate affects the points of its inflexion. The study also examined the problem of assessing the effect of these coefficients on the reconstruction of pressure; it was established that neglecting the elasticity of the seal leads to significant errors in determining the amplitude of the load.
Key words: lamina cribrosa, sclera, inhomogeneity, elastic fixing, Ritz method, reconstruction.
Получено 30 ноября 2016
