CC BY
16
УДК 531.51-7 Вестник СПбГУ. Прикладная математика. Информатика... 2021. Т. 17. Вып. 4
MSC 92С10, 70Е55
Математическое моделирование распределения давления при деформациях межпозвонкового диска
Д. А. Чубарое, В. П. Трегубое, Н. К. Егорова
Санкт-Петербургский государственный университет, Российская Федерация, 199034, Санкт-Петербург, Университетская наб., 7-9
Для цитирования: Чубарое Д. А., Трегубое В. П., Егорова Н. К. Математическое моделирование распределения давления при деформациях межпозвонкового диска // Вестник Санкт-Петербургского университета. Прикладная математика. Информатика. Процессы управления. 2021. Т. 17. Вып. 4. С. 381-388. https://doi.org/10.21638/11701/spbul0.2021.406
На основании анализа приведенных в литературе механических моделей межпозвонкового диска делается вывод о том, что используемые в них конечно-элементные вычислительные сетки не отражают реальную структуру межпозвонкового диска. В связи с этим была построена механическая модель межпозвонкового диска, структура которой близка к его реальной структуре. С помощью предложенной модели определялась динамика распределения давления в межпозвонковом диске при повороте одного из позвонков на заданный угол. Для оценки образующихся при этом выпуклостей каждой ячейки в структуре фиброзного кольца использовались метод Рэлея и его модификация, что позволило рационально рассчитать объемы ячеек. При вычислении давления в каждой ячейке применялась его линейная зависимость от величины прогиба ячейки. В результате работы описанного алгоритма была определена динамика давления в каждой ячейке модели межпозвонкового диска при изменении угла между позвонками.
Ключевые слова: механическая модель, межпозвонковый диск, распределение давления, математическое моделирование.
1. Введение. Межпозвонковый диск является основным звеном, связывающим соседние позвонки. Основа его структуры — студенистое (пульпозное) ядро, окруженное общим фиброзным кольцом, которое, в свою очередь, представляет собой совокупность внутренних колец, наполненных несжимаемой жидкостью (рис. 1) и заключенных между жесткими хрящевыми пластинами. Благодаря своей деформируемости межпозвонковые диски обеспечивают подвижность позвоночника.
Позвоночный сегмент состоит из двух соседних позвонков, связанных между собой межпозвонковым диском (рис. 2). Кроме того, отростки позвонков соединены межпозвонковыми связками и связками, проходящими вдоль позвоночника (рис. 3).
© Санкт-Петербургский государственный университет, 2021
Рис. 1. Общий вид структуры межпозвонкового диска
Рис. 2. Позвоночный сегмент, состоящий из двух соседних позвонков, соединенных межпозвонковым диском
Рис. 3. Связки позвоночного сегмента
Моделированию межпозвонкового диска посвящено достаточно большое количество публикаций (см., например, [1-6]). Отдельно следует упомянуть работу [7], где проведен анализ различных теорий материалов, используемых в моделях движения поясничного сегмента, а также статью [8], в которой рассмотрено применение нового метода калибровки для трехмерной конечно-элементной модели фиброзных колец поясничного отдела позвоночника. В [9] был проведен численный анализ механического поведения межпозвонкового диска с учетом структуры коллагеновых волокон. В [10] подчеркивается, что изменение давления в фиброзных кольцах при движении шейного отдела позвоночника приводит к активизации обменных процессов в нем. Наиболее полный обзор печатных работ представлен группой авторов в [11].
2. Механическая модель. В качестве механической модели межпозвонкового диска рассмотрим объемный многогранник, ограниченный плоскими пластинами и составленный из трапециевидных ячеек, как показано на рис. 4.
Внутренние кольцевые многогранники моделируют внутрифиброзные кольца. Таким образом, границы многогранников являются хордами по отношению к моделируемым кольцам. Следует отметить, что при увеличении числа радиальных разбиений структура модели, а следовательно, и его вычислительная сетка будут в большей степени приближены к реальной структуре межпозвонкового диска по сравнению с другими моделями (см., например, рис. 5).
Центральная часть, моделирующая пульпозное ядро, исключается из деления и дальнейших расчетов. Будем считать, что внутренности ячеек заполнены несжи-
Рис. 4• Механическая модель межпозвонкового диска а — общий вид; б — радиальный фрагмент модели.
Рис. 5. Модель межпозвонкового диска [4]
маемой жидкостью. С учетом наличия постоянной нагрузки на позвоночный столб и сжатия дисков межпозвонковыми связками давление Ро внутри ячеек превосходит атмосферное. В связи с этим в исходном положении хордовые стенки ячеек выгнуты в противоположную от центра диска сторону (рис. 6).
Рис. 6. Механическая модель межпозвонкового диска (К — высота ячейки, А — ширина прогиба)
При повороте соседних позвонков в сагиттальной плоскости на некоторый угол а вокруг центральной точки О (рис. 7) произойдет перераспределение давления, которое за некоторый промежуток времени выравнивается с изменением прогиба стенок ячеек, как показано для сагиттальной плоскости. Это изменение происходит так, чтобы сохранялся объем ячеек ввиду несжимаемости жидкости.
Для упрощения определения изменяющейся во времени формы прогиба стенки воспользуемся известным методом Рэлея, в котором предлагается определять изменяемую форму дуги как произведение заданной наперед формы / (х) на множитель д(£), как показано на рис. 8.
Сразу следует заметить, что у Рэлея метод был продемонстрирован для процесса изменения во времени формы балки. Однако, согласно самому методу, его можно применить и для другого множителя. В изучаемом случае в качестве множителя была
Рис. 1. Изменение прогиба стенок ячеек после поворота на угол а в сагиттальной плоскости
Рис. 8. Иллюстрация к методу Рэлея
выбрана ширина прогиба ячейки А, которая будет различна для разных ячеек. Эти ячейки можно пронумеровать индексами г и ^ где г — номер кольца в модели, а ] — номер ячейки в кольце, т. е. А—. Тогда можно было бы записать, что форма ячейки задается произведением А— на / (х). Но в методе Рэлея расстояние между точками опоры (как видно на рис. 8) остается постоянным, в то время как в рассматриваемом случае оно является высотой ячейки к—, которая также будет неодинаковой для разных ячеек.
В качестве функции /(х) была выбрана синусоида причем величина
к— зависит еще и от угла поворота позвонков а.
3. Математическое моделирование. Алгоритм расчетов. Для выполнения дальнейших вычислений предположим, что давление внутри ячейки Р— в каждый момент времени пропорционально ширине прогиба А—.
В первую очередь необходимо рассчитать исходные объемы V— и ширину прогибов стенок ячеек А—. В связи с тем, что ширина прогиба принята пропорциональной давлению, то
А— =
г.]
К '
где К — коэффициент жесткости стенок.
Объемы ячеек вычисляются согласно формуле
V— = Б— к— (а) + V* - V-и,
в которой вц = \ Бт(ф) (гц(а)2 — (се)2) — площадь верхней грани ячейки, Гц — радиус ячейки (расстояние от центра диска до большего основания грани ячейки),
ф — шаг радиального деления, V? = агз Агз эт (х к ^ ^ (1х = ац _
объем выпуклости (рис. 9).
Модель имеет одну степень свободы — угол между гранями диска и положительным направлением горизонтальной оси. В систему в течение некоторого времени вносится возмущение в виде изменения угла наклона а от а0 = 0 до заданного предельного значения ацш по линейному закону.
Рис. 9. Выпуклая часть ячейки
Каждый шаг алгоритма можно разбить на следующие этапы:
1) пересчет угла (если время работы алгоритма не превышает Т):
. alim
а = ak+i = ак + —jT\
2) расчет изменения объемов за счет диффузии:
AVj — D (Pi-1j + Pi+1j - 2Pij) + d -i + Pij+i - 2Pj),
здесь D — коэффициент межкольцевой диффузии, d — коэффициент внутрикольце-вой диффузии;
3) пересчет объемов:
Vnew — V + AV ••
Vij — Vij + A Vij;
4) пересчет ширины прогибов стенок ячеек Aij-:
где aij — величина большего основания трапеции в верхней грани ячейки;
5) пересчет давлений:
р new _ к Anew
Pij K Aij '
Процесс продолжается до тех пор, пока разность давлений между всеми ячейками будет не больше заданного е.
4. Результаты. На рис. 10, а и б представлены графики динамики давлений внутри двух выбранных наугад ячеек, одна из которых находится в области растяжения, а вторая — в области сжатия.
Рис. 10. Динамика давления в ячейке в областях растяжения (а) и сжатия (б)
На рис. 10, о показано, что давление в ячейке, находящейся в области растяжения межпозвонкового диска, сначала падает, а потом восстанавливается до прежней величины, на рис. 10, б— что давление в ячейке, расположенной в области сжатия, — наоборот: сначала возрастает, а затем уменьшается до исходной величины.
5. Заключение. Построенная механическая модель межпозвонкового диска позволяет моделировать процесс стабилизации внутреннего давления при динамическом изменении угла меду соседними позвонками. Важное практическое значение имеет возможность определить время этой стабилизации. Ключевой особенностью модели является то, что она отражает внутреннюю структуру межпозвонкового диска.
Литература
1. Castro А. P. G., Wilson W., Huyghe J.M., Ito К., Alves J.L. Intervertebral disc creep behavior assessment through an open source finite element solver // Journal of Biomechanics. 2013. Vol. 47. P. 297301.
2. Gao X., Zhu Q., Gu W. Prediction of glycosaminoglycan synthesis in intervertebral disc under mechanical loading // Journal of Biomechanics. 2016. Vol. 49. P. 2655-2661.
3. Guo L.-X., Li R., Zhang M. Biomechanical and fluid flowing characteristics of intervertebral disc of lumbar spine predicted by poroelastic finite element method // Acta of Bioengineering and Biomechanics. 2016. Vol. 18. N 2. P. 19-29.
4. Jacobs N. Т., Cortes D. H., Peloquin J. M., Vresilovic E. J., Elliott D. M. Validation and application of an intervertebral disc finite element model utilizing independently constructed tissue-level constitutive formulations that are nonlinear, anisotropic, and time-dependent // Journal of Biomechanics. 2014. Vol. 47. P. 2540-2546.
5. Schmidt, H., Bashkuev M., Galbusera F., Wilke H.-J., Shirazi-Adl A. Finite element study of human lumbar disc nucleus replacements // Computer Methods in Biomechanics and Biomedical Engineering. 2014. Vol. 17. Iss. 16. P. 1762-1776.
6. Zanjani-Pour S., Winlove C. P., Smith C. W., Meakin J. R. Image driven subject-specific finite element models of spinal biomechanics // Journal of Biomechanics. 2016. Vol. 49. P. 919-925.
7. Gohari E., Nikkhoo M., Haghpanahi M., Parnianpour M. Analysis of different material theories used in a FE model of a lumbar segment motion // Acta of Bioengineering and Biomechanics. 2013. Vol. 15. N 2. P. 33-41.
8. Schmidt H., Heuer F., Simon U., Kettler ARohlmann AClaes L., Wilke H.-J. Application of a new calibration method for a three-dimensional finite element model of a human lumbar annulus fibrosus // Clinical Biomechanics. 2006. Vol. 21. P. 337-344.
9. Мерой E. А., Натали A. H., Паван П. Г., Скарпа К. Численный анализ механического поведения межпозвонкового диска с учетом структуры коллагеновых волокон // Рос. журн. биомеханики. 2005. Т. 9. № 1 (33). С. 36-51. "
10. Жарнов А. М., Жарнова О. А. Биомеханические процессы в межпозвонковом диске шейного отдела позвоночника при его движении // Рос. журн. биомеханики. 2013. Т. 17. № 1 (59). С. 32-40.
11. Хорошев Д. В., Ильялов О. Р., Устюжанцев Н. Е., Няшин Ю. И. Биомеханическое моделирование межпозвоночного диска поясничного отдела человека — современное состояние проблемы // Рос. журн. биомеханики. 2019. Т. 23. № 3. С. 411-422.
Статья поступила в редакцию 19 июля 2021 г.
Статья принята к печати 13 октября 2021 г.
Контактная информация:
Чубарое Дмитрий Александрович — студент; st050417@student.spbu.ru
Трегубое Владимир Петрович — д-р физ.-мат. наук, проф.; v.tregubov@spbu.ru
Егорова Надежда Константиновна — аспирант; nadezhda_ego@mail.ru
Mathematical modeling of pressure distribution during deformations of the intervertebral disc
D.A. Chubarov, V. P. Tregubov, N.K. Egorova
St. Petersburg State University, 7-9, Universitetskaya nab., St. Petersburg, 199034, Russian Federation
For citation: Chubarov D. A., Tregubov V. P., Egorova N. K. Mathematical modeling of pressure distribution during deformations of the intervertebral disc. Vestnik of Saint Petersburg University. Applied Mathematics. Computer Science. Control Processes, 2021, vol. 17, iss. 4, pp. 381-388. https://doi.org/10.21638/11701/spbul0.2021.406 (In Russian)
Based on the analysis of the mechanical models of the intervertebral disc given in the literature, it is concluded that the finite element computational grids used in them do not reflect the real structure of the intervertebral disc. In this regard, a mechanical model of the intervertebral disc was built, the structure of which is close to its real structure. The proposed model was used to determine the dynamics of the pressure distribution in the intervertebral disc when one of the vertebrae is rotated by a given angle. To determine the resulting bulges of each cell in the structure of the fibrous ring, the Rayleigh method and its modification were used. This made it possible to rationally calculate the volumes of cells. When calculating the pressure in each cell, its linear dependence on the deflection value of the cell was used. As a result of the proposed algorithm, the pressure dynamics in each cell of the intervertebral disc model was determined when the angle between the vertebrae changed.
Keywords: mechanical model, intervertebral disc, pressure distribution, mathematical modeling.
References
1. Castro A. P.O., Wilson W., Huyghe J.M., Ito K., Alves J.L. Intervertebral disc creep behavior assessment through an open source finite element solver. Journal of Biomechanics, 2013, vol. 47, pp. 297301.
2. Gao X., Zhu Q., Gu W. Prediction of glycosaminoglycan synthesis in intervertebral disc under mechanical loading. Journal of Biomechanics, 2016, vol. 49, pp. 2655-2661.
3. Guo L.-X., Li R., Zhang M. Biomechanical and fluid flowing characteristics of intervertebral disc of lumbar spine predicted by poroelastic finite element method. Acta of Bioengineering and Biomechanics, 2016, vol. 18, no. 2, pp. 19-29.
4. Jacobs N. T., Cortes D. H., Peloquin J. M., Vresilovic E. J., Elliott D. M. Validation and application of an intervertebral disc finite element model utilizing independently constructed tissue-level constitutive formulations that are nonlinear, anisotropic, and time-dependent. Journal of Biomechanics, 2014, vol. 47, pp. 2540-2546.
5. Schmidt H., Bashkuev M., Galbusera F., Wilke H.-J., Shirazi-Adl A. Finite element study of human lumbar disc nucleus replacements. Computer Methods in Biomechanics and Biomedical Engineering, 2014, vol. 17, iss. 16, pp. 1762-1776.
6. Zanjani-Pour S., Winlove C.P., Smith C.W., Meakin J.R. Image driven subject-specific finite element models of spinal biomechanics. Journal of Biomechanics, 2016, vol. 49, pp. 919-925.
7. Gohari E., Nikkhoo M., Haghpanahi M., Parnianpour M. Analysis of different material theories used in a FE model of a lumbar segment motion. Acta of Bioengineering and Biomechanics, 2013, vol. 15, no. 2, pp. 33-41.
8. Schmidt H., Heuer F., Simon U., Kettler A., Rohlmann A., Claes L., Wilke H.-J. Application of a new calibration method for a three-dimensional finite element model of a human lumbar annulus fibrosus. Clinical Biomechanics, 2006, vol. 21, pp. 337-344.
9. Meroi E. A., Natali A. N., Pavan P. G., Skarpa K. Chislennyj analiz mekhanicheskogo povedeniya mezhpozvonkovogo diska s uchetom struktury kollagenovyh volokon [Numerical analysis of the mechanical behavior of the intervertebral disc taking into account the structure]. Russian Journal of Biomechanics, 2005, vol. 9, no. 1 (33), pp. 36-51. (In Russian)
10. Zharnov A.M., Zharnova O.A. Biomekhanicheskie processy v mezhpozvonkovom diske shejnogo
otdela pozvonochnika pri ego dvizhenii [Biomechanical processes in the intervertebral disc of the cervical spine during its movement]. Russian Journal of Biomechanics, 2013, vol. 17, no. 1 (59), pp. 32-40. (In Russian)
11. Horoshev D.V., Il'yalov O.R., Ustyuzhancev N. E., Nyashin Yu.I. Biomekhanicheskoe modeli-rovanie mezhpozvonochnogo diska poyasnichnogo otdela cheloveka — sovremennoe sostoyanie problemy [Biomechanical modeling of the intervertebral disc of the human lumbar spine — current state of the problem]. Russian Journal of Biomechanics, 2019, vol. 23, no. 3, pp. 411-422. (In Russian)
Received: July 19, 2021.
Accepted: October 13, 2021.
Authors' information:
Dmitry A. Chubarov — Student; st050417@student.spbu.ru
Vladimir P. Tregubov — Dr. Sci. in Physics and Mathematics, Professor; v.tregubov@spbu.ru
Nadezhda K. Egorova — Postgraduate Student; nadezhda_ego@mail.ru
