В настоящей работе рассматриваются вопросы расчета напряженно-деформированного состояния роговицы глаза под действием внутреннего давления. Определение напряжений и деформаций в роговице основано на решении нелинейной краевой задачи теории упругости с использованием численного метода конечных элементов. Предложенный метод решения позволяет учитывать реальную геометрию роговицы как в предоперационный, так и в послеоперационный период. Сравнение экспериментальных и расчетных величин перемещений для множества точек свидетельствует об эффективности предложенного метода расчета.
Modeling of stress-deformed status of the eye cornea under the influence of intraocular pressure
The paper concerns questions of calculation of stress-deformed state of the eye cornea under the influence of inner pressure. Definition of stresses and deformations in the cornea is based on the solution of non-linear boundary problem of elasticity theory using numerical method of finite elements. The suggested method of decision gives a possibility to take into account real corneal geometry both in preoperative and postoperative period. Comparison of experimental and calculated values of displacements for a multitude of the node points speaks for the efficacy of the suggested calculation method.
В последнее время вопросам математического моделирования различного рода процессов, протекающих в биологических системах, уделяется повышенное внимание со стороны специалистов различного направления: врачей, биологов, математиков, инженеров. Моделирование биологического процесса на основе надежной расчетной математической модели позволяет более полно понять причины влияния различных факторов на поведение исследуемого объекта, оценить степень их влияния на выбранные методы лечения. В работе рассматриваются вопросы математического моделирования роговицы глаза под действием внутреннего давления на основе нелинейной краевой задачи теории упругости с использованием численного метода решения. При этом анализируются два состояния роговицы в дооперационный и послеоперационный периоды.
Описание расчетной модели
Два фрагмента на рис. 1 характеризуют геометрию роговицы глаза. На верхнем фрагменте показано поперечное сечение роговицы, отчетливо характеризующее изменение его толщины. Роговица имеет в плане форму эллипса с большой осью 13,5 мми малой осью 12,6 мм. Толщина роговицы изменяется в интервале от 0,5 ммв центральной части до 1,2 ммна периферии. Радиус кривизны R в плоскости YOX — принимался равным7,8 мм, а в плоскости ZOY —6,5 мм.
Рисунок 1. Геометрия роговицы глаза (описание в тексте)
В процессе проведения операции на круге диаметром6 ммв плане изначально с использованием лазера вырезается лоскут определенной толщины. На рис. 2, характеризующем аксонометрический вид роговицы, эта область отмечена красным цветом. При проведении Wavefront-Guided LASIK, операции с использованием микрокератома толщина лоскута может изменяться в довольно широком диапазоне — от 146 мкм в центре роговицы до 200 мкм на краю.
Рисунок 2. Аксонометрия роговицы
При использовании Wavefront-Guided FEMTOLASIK отмеченный выше срез лоскута производится лучом лазера, поэтому толщина лоскута в этом случае является постоянной величиной, равной 100 мкм, отклонения от этой толщины не превышают 10 мкм. Затем вырезанный лоскут откладывается в сторону от области обозначенного хирургического вмешательства и с помощью эксимер-лазера производится абляция нижерасположенных слоев роговицы. После выполнения операции лоскут возвращается обратно на свое место и выполняет, таким образом, функцию «пластыря».
Рисунок 3. Формирование расчетной модели: а — изменяемая в ходе операции центральная часть роговицы; б — неизменная часть модели роговицы
Формирование расчетной конечноэлементной модели проводилось в предположении, что изменяется только центральная часть роговицы с диаметром 6 мм(рис. 3а). Неизменная часть модели роговицы изображена на правом фрагменте рис. 3б. Такой подход к созданию расчетной модели позволяет достаточно оперативно изменять ту часть роговицы, которая изменяется при выполнении операции и не менять остальную часть. В период подготовки пациента к операции, по четырем меридианам роговицы, расположенных под углом 45° по отношению друг к другу, на отрезке 6 ммопределяется ее толщина с шагом равным 1 ммвдоль диаметра. На рис. 3 эти толщины для одного из отмеченных выше сечений обозначены через h0 .. h6. Поскольку центральная толщина h0 роговицы в каждом сечении должна быть одной и той же, то в случае ее отличия при выполнении измерений в различных сечениях роговицы в модели принималась наименьшая величина. Таким образом, при моделировании роговицы в дооперационный период имеется набор из двадцати пяти точно измеренных значений толщин роговицы. Этот набор толщин является каркасом, на котором строится расчетная модель изменяемой части роговицы. Использование линейной поверхности, которая может быть натянута на обозначенный каркас, будет сопровождаться аномальными величинами, характеризующими напряженно-деформированное состояние. Этот факт обусловлен тем, что в точках измерения толщин нарушается непрерывность первой производной от линейной функции, описывающей изменение толщины слева и справа от точки. Для устранения отмеченной аномалии для моделирования использовались кубические сплайны, позволяющие обеспечить непрерывность самой функции, ее первой и второй производных. На рис. 4а изображены ломаные линии, характеризующие изменение толщины роговицы по диаметру для четырех сечений оперируемого глаза, на рис. 4б приведен набор измеренных толщин для одного из сечений и показана линейная и гладкая интерполяция с использованием кубических сплайнов.
Рисунок 4. Аппроксимация толщины роговицы с использованием кубических сплайнов (пояснение в тексте)
Поскольку в анализируемом случае значение первой (или второй) производной на границе неизвестно, то можно задать естественные ограниченные условия F»(A1)=0, F»(A2)=0, и получить естественный сплайн. Погрешность интерполяции естественным сплайном составляет O(h 2). При этом максимум погрешности наблюдается в окрестностях ограниченных узлов, во внутренних узлах точность интерполяции значительно выше. Между точками, расположенными на одинаковом удалении от центра роговицы, также использовалась аппроксимация кубическими сплайнами.
Рисунок 5. Конечноэлементная модель роговицы
Расчетная конечноэлементная модель роговицы глаза изображена на рис. 5. Для моделирования использовался объемный конечный элемент с десятью узлами и тридцатью степенями свободы. Модель содержит 14887 узлов и 8670 конечных элементов и позволяет учитывать реальную геометрию роговицы и моментное напряженное состояние, поскольку вычислительный процесс основан на решении нелинейной краевой задачи теории упругости. Видно, что в центральной части роговицы сетка расчетной конечноэлементной модели имеет меньшие размеры по сравнению с периферийной частью.
Рисунок 6. Расчетные модели роговицы в предоперационный и послеоперационный периоды (пояснение в тексте)
Два фрагмента на рис. 6 характеризуют размеры конечных элементов в поперечном сечении роговицы. На рис.6а изображена расчетная модель, построенная по результатам измерений в предоперационном периоде, до проведения операции. На рис. 6б идентичная модель после снятия лоскута и абляции. Сравнение представленных расчетных моделей позволяет сформировать наглядное представление об их заметном отличии в центральной части и незначительном отличии в периферийной части роговицы. Полагалось, что изнутри роговица была нагружена равномерным давлением величиной p =18 мм рт. ст.
Рисунок 7. Схема нагружения роговицы
Схема нагружения роговицы изображена на рис. 7. Для решения нелинейной краевой задачи необходимо сформулировать граничные условия, в связи с этим на периферийной поверхности расчетной модели роговицы накладывались ограничения на все перемещения в узлах. Материал роговицы полагался изотропным с модулем упругости E = 0,83 МПа и коэффициентом поперечной деформации μ = 0,48. На основе данных, приведенных в [1], следует вывод о том, что модуль упругости роговицы может изменяться в довольно широком диапазоне: от 0,1 до 1,0 МПа. Расчетная модель позволяет довольно просто изменять механические характеристики при выполнении расчета напряженно-деформированного состояния.
Рисунок 8. Мозаика перемещений в роговице глаза, мкм
Распределение перемещений в роговице глаза изображено на рис. 8. При этом каждому цветовому оттенку на приведенной ниже шкале соответствует строго определенный интервал перемещений. Видно, что окрашенное красным цветом пятно овальной формы располагается в центральной части роговицы и соответствует области наибольших перемещений. С увеличением расстояния от центра до анализируемой области, величина перемещений уменьшается. Этот факт объясняется тем, что изгибная жесткость роговицы с увеличением расстояния от центра заметно возрастает и в периферийных областях превышает аналогичную величину жесткости в центральной части более чем на порядок. Это обстоятельство свидетельствует о том, что использование мембранной теории тонких оболочек применительно к расчету напряженного состояния роговицы приводит к недостоверным результатам.
Рисунок 9. Изменение изгибных напряжений в диаметральном сечении роговицы
Два фрагмента на рис. 9 иллюстрируют изменение изгибных напряжений в диаметральном сечении 180°-0° роговицы. Рис. 9а характеризует изменение изгибных напряжений в роговице до проведения операции, рис. 9б — после выполнения операции. Видно, что в обоих случаях наибольшие напряжения наблюдаются на периферийной поверхности роговицы. Их максимальная величина равна 13,42 кПа в исходной роговице и 17,46 кПа в роговице после операции. В центральной части роговицы в исходном состоянии максимальная величина анализируемых напряжений составляет 8,72 кПа, а в послеоперационной роговице — 11,67 кПа. Кроме того, в центральной части роговицы напряжения по толщине не меняют знак и в любой точке сечения являются растягивающими. В краевых зонах, напротив, изгибные напряжения по толщине роговицы изменяют знак с минуса на наружной поверхности на плюс — на внутренней.
Уменьшение толщины роговицы сопровождается увеличением изгибных напряжений в центральной части роговицы на 25%, а также более чем двукратным возрастанием применительно к сжимающим изгибным напряжениям на периферийной поверхности. В работе [2] приведены экспериментальные кривые в координатах «Напряжение»-«Деформация», которые характеризуют величину напряжений в роговице глаза человека в зависимости от степени деформации. При деформации 7,38%, которая имеет место в периферийной части роговицы по данным расчета, величина напряжений разрыва, полученная по отмеченной кривой, незначительно превышает 50 кПа. На основе такого сравнения можно сделать заключение о том, что суммарное уменьшение толщины роговицы на 120..150 мкм в результате абляции и снятия лоскута не приводит к напряженному состоянию роговицы, близкому к потере прочности.
С целью оценки достоверности полученных результатов расчетов был выполнен ряд вычислений для пациентов, прооперированных посредством Wavefront-Guided FEMTOLASIK операции. В таблице 1 приведены результаты измерений и вычислений. Измерения толщины роговицы проводились до операции, в ее процессе и после. При выполнении операции измерялась толщина лоскута, величина абляции. Экспериментальные и расчетные данные приведены для четырех диаметральных сечений на диаметре роговицы 6 мм. Вычисления проводились для двух значений толщин роговицы: исходной толщины роговицы до проведения операции и толщины роговицы после абляции.
Таблица 1.
Результаты пахиметрических показателей роговицы
180°-0°, радиус, мм | Толщина до операции, мкм | Абляция, мкм | Толщина лоскута, мкм | Толщина после абляции, мкм | Измеренная толщина после операции, мкм | Расчетная толщина после операции, мкм | Эксперимен- тальное изменение толщины, мкм | Разность перемеще-ний, мкм |
-3 | 548 | 13,53 | 96 | 459 | 555 | 534,47 | 20,53 | 18,33 |
-2 | 517 | 43,25 | 97 | 413 | 510 | 473,75 | 36,25 | 34,21 |
-1 | 503 | 57,74 | 95 | 396 | 491 | 445,26 | 45,74 | 41,27 |
0 | 501 | 62,00 | 95 | 396 | 491 | 439 | 52 | 48,43 |
1 | 515 | 57,96 | 95 | 407 | 502 | 457,04 | 44,96 | 43,58 |
2 | 539 | 43,16 | 95 | 414 | 509 | 495,84 | 13,16 | 14,12 |
3 | 574 | 12,87 | 95 | 503 | 598 | 561,13 | 36,87 | 35,76 |
225°-45° | ||||||||
-3 | 539 | 13,94 | 96 | 473 | 569 | 525,06 | 43,94 | 38,27 |
-2 | 518 | 42,89 | 96 | 419 | 515 | 475,11 | 39,89 | 36,03 |
-1 | 506 | 57,25 | 96 | 349 | 485 | 448,75 | 36,25 | 33,92 |
0 | 506 | 62,00 | 95 | 402 | 497 | 444 | 53 | 48,17 |
1 | 523 | 58,1 | 94 | 382 | 476 | 464,9 | 11,1 | 12,12 |
2 | 557 | 44,2 | 95 | 402 | 497 | 512,8 | -15,8 | -11,38 |
3 | 606 | 14,14 | 95 | 488 | 583 | 591,86 | -8,86 | -6,14 |
270°-90° | ||||||||
-3 | 557 | 15,10 | 95 | 452 | 547 | 541,9 | 5,10 | 3,76 |
-2 | 519 | 43,25 | 96 | 403 | 499 | 475,75 | 23,25 | 20,13 |
-1 | 515 | 57,42 | 95 | 392 | 487 | 457,58 | 29,42 | 27,02 |
0 | 518 | 62,00 | 97 | 382 | 479 | 456 | 23 | 24,04 |
1 | 514 | 58,33 | 96 | 391 | 487 | 455,67 | 31,33 | 30,41 |
2 | 536 | 43,53 | 96 | 432 | 528 | 492,47 | 35,53 | 33,23 |
3 | 594 | 13,93 | 96 | 530 | 616 | 580,07 | 35,93 | 32,67 |
315°-135° | ||||||||
-3 | 557 | 13,14 | 94 | 504 | 598 | 543,86 | 54,14 | 49,87 |
-2 | 528 | 43,89 | 96 | 452 | 548 | 484,11 | 63,89 | 58,28 |
-1 | 510 | 57,11 | 96 | 405 | 501 | 452,89 | 48,11 | 46,07 |
0 | 511 | 62,00 | 96 | 389 | 485 | 449 | 36 | 34,43 |
1 | 517 | 57,85 | 97 | 393 | 490 | 459,15 | 30,85 | 27,96 |
2 | 538 | 43,87 | 94 | 408 | 502 | 494,13 | 7,87 | 6,18 |
3 | 585 | 13,19 | 95 | 485 | 580 | 571,81 | 8,19 | 6,91 |
При этом принималось во внимание то обстоятельство, что лоскут, возвращенный на роговицу после операции, выполняет защитную функцию пластыря и не оказывает влияние на жесткость роговицы. В предпоследнем столбце таблицы приведены величины, характеризующие экспериментальное изменение толщины, полученные как разность между измеренной и расчетной толщиной роговицы после операции. В последнем столбце таблицы приведены величины, характеризующие разность перемещений роговицы до и после проведения операции. С целью оценки корреляции между изменением толщины роговицы и разностью перемещений для всех четырех диаметральных сечений, на рис. 10 изображены восемь линий, обозначенных цифрами от 1 до 8.
Рисунок 10. Корреляция расчетных и экспериментальных значений перемещений и толщин
Нечетные линии характеризуют разность толщин роговицы, а четные — разность перемещений роговицы до и после операции. Надпись над линиями указывает их принадлежность к конкретному диаметральному сечению. Анализ приведенных линий свидетельствует о заметной корреляции между линиями, расположенными в одном диаметральном сечении и менее выраженной корреляцией между линиями в разных сечениях. Довольно значительный разброс анализируемых величин можно объяснить рядом причин, к числу важнейших из которых следует отнести анизотропию материала роговицы, неоднородность механических свойств по объему. В статье [3] отмечено, что при изменении внутриглазного давления в диапазоне от 15 до30 ммрт. ст. модуль упругости роговицы не изменяет своего значения, но с увеличением давления его величина резко возрастает. Этот факт можно объяснить только нелинейностью материала роговицы или изменением его структуры. Предложенная математическая модель позволяет проводить моделирование роговицы глаза для анизотропного материала при рассматриваемом нагружении. Вместе с тем для получения достоверных результатов расчетов необходимо располагать достаточно точно определенными механическими характеристиками материала роговицы. Имеющиеся в литературе [4-7] экспериментальные данные, относящиеся к определению предела прочности и модуля упругости роговицы, свидетельствуют об их значительном отличии в зависимости от метода испытаний и анализа.
Выводы
Использование математического моделирования на основе конечно-элементной модели роговицы глаза подтвердило его эффективность [8, 9]. В предложенном методе анализа напряженно-деформированного состояния роговицы использован объемный конечный элемент и нелинейная модель. Это позволило учитывать переменную толщину роговицы по всему её объему и строить расчетную математическую модель с учетом реальных толщин, определенных с помощью оптической измерительной техники. Определение напряженно-деформированного состояния осуществлялось на основе решения нелинейной краевой задачи теории упругости. Расчет напряжений и перемещений во множестве узловых точек роговицы позволяет хирургу сформировать более полное представление о влиянии параметров, изменяемых в процессе операции, на состояние пациента.
О.А. Костин, К.А. Гончаров
Екатеринбургский филиал МНТК «Микрохирургия глаза» им. акад. С.Н. Федорова
Уральский федеральный университет им. Б.Н. Ельцина-
Костин Олег Александрович — кандидат медицинских наук, заведующий отделением лазерной хирургии
Литература:
1. Хусаинов Р.Р., Цибульский В.Р., Якушев В.Л. Моделирование деформации глаза при измерении внутриглазного давления оптическим методом // Журнал вычислительной математики и математической физики. — 2011. — Т. 51, № 2. — С. 349-362.
2. Аветисов С.Э., Воронин Г.В. Экспериментальное исследование механических характеристик роговицы после эксимерлазерной фотоабляции // Клин. офтальмология. — 2001. — № 3. — С. 83-86.
3. Wollensak G., Spoerl E., Seiler T. Stress-strain measurements of human and porcine corneas after riboflavin–ultraviolet-A-induced cross-linking // Journal of cataract& refractive surgery. — 2003. — Vol. 29. — P. 1780-1785.
4. Bryant M.R., McDonnell P.J. Constitutive laws for biomechanical modeling of refractive surgery // J. Biomech.Eng.— 1996. — 118. — P. 473-481.
5. Buzard K.A. Introduction to biomechanics of the cornea // Refract. Corneal Surg. — 1992. — № 8. — P. 127-138.
6. Jun Liu, Cynthia J. Roberts. Influence of corneal biomechanical properties on intraocular pressure measurement: Quantitative analysis // Journal of cataract& refractive surgery. — Vol 31. — P. 127-138.
7. Бауэр С. М., Любимов Г. А., Товстик П. Е. Математическое моделирование метода Маклакова измерения внутриглазного давления // Механика жидкости и газа: Изв. РАН. — 2005.— № 1.— С. 24–39.
8. Кочина М.Л. Некоторые результаты моделирования напряженного состояния роговой оболочки глаза // Кибернетика и вычислительная техника. — 1991. — Вып. 90. — С. 97-99.
9. Кочина М.Л., Калиманов В.Г. Результаты моделирования напряженно-деформированного состояния роговицы глаза с использованием системы инженерного анализа ANSYS // Клиническая информатика и Телемедицина. — 2009. — Т. 5, вып. 6. — C. 26-30.