CC BY
42
УДК 519.6: 616.832-02
Р. Л. Седов, С. В. Орлов, К. С. Латышев
ОПТИМИЗАЦИЯ ФИКСИРУЮЩИХ УСТРОЙСТВ ПОЗВОНОЧНИКА НА ОСНОВЕ МАТЕМАТИЧЕСКОЙ МОДЕЛИ ТРЕХПОЗВОНКОВОГО КОМПЛЕКСА ЧЕЛОВЕКА
Проведен анализ математической модели стабильности трехпозвонкового комплекса человека. Оптимизацией механической конструкции с целевой функцией смещения поврежденного позвонка в сторону спинного мозга с сохранением подвижности всей системы получен интервал жесткости стабилизирующей пластины 609,4 <Cct1 <6403,6 н/мм. Доказана целесообразность применения полуригидных конструкций при лечении подобных травм позвоночника.
For backbone complex of the human's spine with alternative of a cuneiform strain mean backbone the stabilizing rigid lamina (for front and back abutment complexes) holds optimization of stiffness of stabilizing of a lamina СсЬ on an objective function of offset 2 backbone, at its sufficient mobility.
Ключевые слова: математическое моделирование, обыкновенные дифференциальные уравнения, численные методы, оптимизация, биомеханика.
Keywords: mathematical modelling, the ordinary differential equations, numerical methods, optimization, biomechanics.
Введение
В работе [1] описывается математическая модель стабильности трехпозвонкового комплекса человека. В основу модели положено математическое описание дифференциальными уравнениями Лагранжа II рода динамических процессов трехпозвонкового комплекса — дискретных сосредоточенных масс, связанных упругодемпфирующими элементами и обладающих определенными геометрическими параметрами. За основу модели принят принцип стабильности позвоночного столба, изложенный Rene [6], согласно которому стабильность позвоночника представлена в вертикальной, горизонтальной и аксиальной плоскостях (ротация), что обеспечивается телами позвонков с дугоотростчатыми суставами, которые связаны между собой упругодемпфирующими элементами (межпозвоночные диски, мышечно-связочный аппарат).
Учитывались следующие параметры трехпозвонкового комплекса:
1. Механическая система является диссипативной.
2. Распределение нагрузок соответствует трехстолбовой концепции Denis [3].
Вестник Российского государственного университета им. И. Канта. 2009. Вып. 10. С. 42 — 48.
3. Предел прочности тел позвонков и упругодемпфирующих элементов, а также их упругая деформация и плотность считались условно установленными по данным работы [2].
4. Изменение геометрических характеристик трехпозвонкового комплекса соответствовало типичным статико-динамическим нарушениям стабильности позвоночника [5; 7].
Расчетная схема фрагмента позвоночника человека, изображающая три позвонка с клиновидным средним позвонком и стабилизирующими конструкциями, показана на рисунке 1.
На расчетной схеме (рис. 1) третий позвонок связывается посредством жестких элементов Сор1 и Сор2 с опорой по оси Х, а первый — по оси У через Су [1]. Для фиксации вариантов нестабильности позвоночника предусмотрено применение условных жестких плоскостных конструкций с коэффициентами жесткости Сст1 и Сеть что позволяет моделировать как жесткие ригидные металлические системы, так и полури-гидные пружинные элементы.
Математическая модель позволяет на основе вычисления внутренних нагрузок опорных комплексов каждого позвонка трехпозвонкового комплекса рассчитывать варианты переломов и нестабильности позвонков в различных зонах при их патологии. Кроме этого, возможен расчет смещения позвонков по оси У под воздействием силы Q2y, что чаще всего является причиной стеноза позвоночного канала и может приводить к сдавливанию дуального мешка. Выбранная динамическая модель трехпозвонкового комплекса человека (рис. 1) является механической системой, для которой уравнение Лагранжа II рода имеет вид
Рис. 1. Расчетная схема трехпозвонкового комплекса с патологией среднего позвонка и ее двухсторонней стабилизацией
д , дТ , дП дФ
д" (д—) + д— + д— = ; к = '
ді дх к дх к дх к
., 7,
где Т, П — кинетическая и потенциальная энергия системы; Ф — диссипативная функция, определяемая спинными мышцами; Qk — внешние воздействия.
В качестве обобщенных координат Хк принимаются координаты хі;
Х2; хз; Х4(яі = Хі; а2 = хз);Х5; Хб(а3 = Хб^ Х5 ); у; Оі = йі + йі; і = і, 2, 3.
Оі
Оі
Оз
В работе [і] выписана система уравнений
С=
й2Х _ йХ _ М =- + В + С йі2 йі X = д,
тгс 1 т10 0 00 0 0 "
т^с! т2 0 00 0 0
0 0 тзс 3 т20 0 0 0
М = 0 0 т20С 3 т4 0 0 0 ;
0 0 0 0 т5 т30 0
0 0 0 0 тз0 т6 0
0 0 0 00 0 М! + м 2 + М 3 _
С1с1 +Ссщ 0 -С1С3 0 Ссщ 0 0'
0 С2 +Ссщ 0 -С 0 8- с С - 0
-С1с1 0 (С1 +Сі К 0 -С3 0 0
0 -С2 0 С2 +С4 0 -С4 0
-Ссщ 0 -С3 0 Сор +Ссі +С3 0 0
0 - -Ссщ 0 С4 0 С4 +Сор2 +Ссі 0
_ -5, 0 53 0 0 0 С У _
.(31 <22 0 0 0 0 0ву ]Г А N X Х 2 Х3 Х 4 Х 5 Х 6 У
а=
ББг = В^т^)^^); ББз = Бз*з8т(Рз)со8(Рз).
Умножая обе части векторного уравнения (1) на обратную матрицу М"1, получим следующее приведенное векторное уравнение:
й х + М-1 • В •—+М-1 • С • X = М-1
йі
2
йі
д,
2
или —X + W • — + А • X = Б, А = М-1 С; W = М-1В; Б = М-10. йі2 йі
(2)
Особенностью численного решения системы с большим числом обыкновенных дифференциальных уравнений первого порядка с квазилинейными коэффициентами, состоящей из двух векторных уравнений, является то, что она решается относительно производных от обобщенных координат траекторий движений материальных точек, являющихся тангенсами углов наклонов касательных к точкам этих траекторий [10]. Даже незначительные погрешности вычисления указанных производных по разностным схемам первого порядка точности О(Ь) приводят к значительным погрешностям численных значений функций обобщенных координат. А с учетом большого числа уравнений (к = 14) и неточностью задания начальных условий для них делают эти разностные схемы неприемлемыми для решения обыкновенных дифференциальных уравнений первого порядка с квазилинейными коэффициентами (2). Поэтому используются разностные схемы, повышающие порядок точности аппроксимации дифференциальных уравнений, основанные на введении полуцелых временных слоев ] + 0,5. Причем система дифференциальных уравнений 1-го порядка (4) решается по разностным схемам второго порядка точности по времени О(Ь2).
Оптимизация
При помоши построенной модели решим следующие задачи.
1. Оптимизация стабилизирующей конструкции трехпозвонкового комплекса.
2. Разработка алгоритма принятия решения об оперативном вмешательстве при различных степенях клиновидности позвонка.
Первая задача. Цель оптимизации — рассчитать параметры стабилизирующей конструкции со следующими свойствами: 1) фиксировать средний клиновидный позвонок; 2) минимизировать его смещение, сохранять физиологическую подвижность трехпозвонкового комплекса; з) максимально увеличить свободу перемещения не более чем на 1 мм [8; 9]. Критерием оптимальности будет являться система требований
где Д; = ^, Ccti, Х1, Х3, р1, р3) - Pinorm , { = 5, 6, — компонентах двумер-
ного вектора Д целевых функций. В силу положительности целевых функций тривиальным решением будет нулевое решение. В этом случае нагрузки Pi тождественно приблизились к константам Pnormi, поэтому вырожденное решение — физиологически допустимое. Целевые функции также зависят от параметров С1 и Сз — коэффициентов жесткости межпозвоночных дисков, но эти параметры выбираются для каждой отдельной практической задачи частным образом. Напомним, что параметры математической модели в разных отделах позвоночника принимают разные значения. Поэтому входные параметры (жесткости, геометрические параметры, углы деформации, масса и др.) выбраны условно, но дальнейшие исследования позволят построить таблицы параметров модели для каждого отдела позвоночника.
Задавались следующие исходные данные (см. рис. 1): С1 = 3,26 • 103н/мм; С2 = 0,92 • 103 н/мм; С3 = 0,46 • 103 н/мм; С4 = 3,26 • 103 н/мм; С5= 0,92 • 103 н/мм; Сб = 0,46 • 103 н/мм; Су = 5 • 103 н/мм; Copl = Cop2 = 3,26 • 103 н/мм; О = 800 Н (внешняя сила приложена к центру тяжести позвонка при у = 21 мм); Р1 = 10°; Р3 = 15° (деформация позвонка); Сст2 = 0 (стабилизирующие пластины справа от позвонков отсутствуют (рис. 1), число временных слоев ] = 1000, а шаг интегрирования по времени т = 10 — 2с.
Вторая задача. Прежде чем решать задачу, необходимо выбрать параметры вариации. Тяжесть патологии зависит от углов клиновидного позвонка. В качестве целевой функции задается функция смещения 2-го позвонка у2^ вдоль оси у, как и в первой задаче,— формула (3).
При варьировании одного из углов клиновидности Р1 в диапазоне 2—40° требуется, чтобы у2^ ^ тт. Угол клина среднего позвонка выражается через углы р1 и р3, участвующих в уравнениях как сумма.
Вариация начинается от суммарного угла в диапазоне 2—55° параллельно с вариацией остаточного угла в диапазоне 2 — 30°. Вне этих интервалов ЭВМ дает слишком большие смещения у2от, что позволяет
получить циклические численные результаты. Их интерпретация сложнее, но решение более практично, чем предыдущее.
Интерпретация полученных результатов
Показаны зависимости смещения 2-го позвонка по оси Y и от коэффициента жесткости стабилизирующей пластины СсЬ (рис. 2, 3). Как видно из графиков, наименьшее смещение поврежденного позвонка произойдет при коэффициенте жесткости пластины СсЬ = 609,4 н/мм, если приложение внешней нагрузки приходится на 2і-й мм (центр тяжести). При такой фиксации подвижность системы сохраняется.
Сс1г Сс1г
Рис. 2. Зависимость Y и Y2cm от жесткости пластины Сс^ в приложении внешней нагрузки к центру тяжести (в прямом положении позвоночника)
Рис. 3. Зависимость Y и Y2cm от жесткости пластины Сс^ в приложении внешней нагрузки к точке у = 0 мм (при наклоне позвоночника)
У2сш.
Сс1
Сс1
Из рисунка 4 следует, что уже начиная с р1 = 5° смещение клиновидного позвонка не превысит физиологической ригидности системы позвоночника. При усугублении повреждения позвонка увеличивается смещение, а значит, возрастает угроза ухудшения патологии. В этом случае требуется оперативное вмешательство. Если же внешняя нагрузка действует непосредственно на центр тяжести, то стабилизация комплекса требует более жесткой пластины СеЬ= 6403,6 н/ мм. Подвижность конструкции также сохраняется (рис. 5).
Рис. 4. Зависимость Y и Y2cm от угла р1 верхней клиновидности в приложении внешней нагрузки к точке у = 21 мм (при наклоне позвоночника)
у = 21 поп
У 21 мши
Рис. 5. Зависимость Y и Y2cm от угла Р$ общей клиновидности в приложении внешней нагрузки к точке у = 21 мм (при наклоне позвоночника)
Выводы
В работе [1] изучалось поведение конструкции при значениях жесткостей пластины 10 и 5 тыс. н/мм. В 1-м случае нагружения среднего позвонка не происходило, но сохранялась неподвижность позвоночника. Во 2-м случае нагрузка увеличивается, соответственно, увеличивается и смещение поврежденного позвонка, усиливающее патологию.
Результаты наших расчетов дают оптимальный интервал жесткости стабилизирующей пластины 609,4 < СсЬ < 6403,6 н/мм.
Таким образом, применение полуригидных конструкций является предпочтительней по отношению к использованию ригидных конструкций, так как сохраненный до некоторой степени объем движений в позвоночно-двигательном сегменте больше приближается к физиологической норме биомеханических свобод перемещений. Вместе с тем для более объективной оценки стабильности патологических типов комплекса необходимо некоторое усложнение концепции и введение дополнительных параметров в систему уравнений, что позволит оценивать не только фиксирующие элементы типа пластин, но и более сложные (такие, как транспедикулярные фиксаторы или протезы межпозвоночных дисков).
По данным модели применение полуригидных конструкций (в том числе и имитирующих работу межпозвоночного диска) ограничено некоторыми функциями, влияющими на изменения биомеханических характеристик комплекса при стандартных нагрузках. Это такие огра-
ничения, как степень клиновидности среднего позвонка в комплексе, степень жесткости передней конструкции.
Вышеизложенное позволяет сделать вывод, что вычисленный интервал жесткости позволит продолжить изыскания надежного фиксатора нового поколения для лечения нестабильных состояний позвоночника, обусловленных клиновидной деформаций среднего позвонка или разрушением его передних колонн.
Работа выполнена при поддержке Российского фонда фундаментальных исследований, проект № 06-01-00396.
Список литературы
1. Орлов С. В., Бобарыкин Н. Д., Латышев К. С. Математическая модель стабильности трехпозвонкового комплекса // Мат. моделирование. 2006. Т. 18, № 10.
2. Громов А. П. Биомеханика травмы. М., 1979. С. 179—210.
3. Denis F. Spinal instability as defined by the three column spine concept in acute spinal trauma / / Clin. Orthop. 1984. 189:65.
4. Denis F. The three column spine and significance in the classification of acute thoracolumbar spinal injuries // Spine. 1983. N. 8. C. 817—831
5. Fergusson R., Tencer A., Woodard P., Allen A. Biomechanical comparison of spinal fracture models and the stabilizing effects of posterior instrumentations / / Spine. 1988. 13:453.
6. Reno Louis. Surgery of the Spine. Berlin; Heidelberg, 1983.
7. Haher Th. R., Felmly W. T., O'Brien M. Thoracic and Lumbar Fractures: Diagnosis and Management / / Spinal Surgery. 1991. Vol. 2. P. 857 — 910.
8. Седов Р. Л. Оптимизация фиксирующих устройств на основе математической модели // Фундаментальная и клиническая медицина: матер. XI Всерос. мед-биол. конф. молодых исследователей «Человек и здоровье». СПб., 2008. С. 335.
9. Седов Р. Л., Орлов С. В., Бобарыкин Н. Д. О математическом моделировании физических свойств стабилизирующих конструкций при лечении травм позвоночника человека // Высокие технологии, фундаментальные исследования, промышленность: сб. тр. Шестой междунар. науч.-практ. конф. «Исследование, разработка и применение высоких технологий в промышленности». СПб., 2008. С. 161.
10. Sedov R., Orlov S., Bobarykin N. Modeling of human's three-vertebra system and optimisation of fixation constructions // The International Conference Mathematical modeling and computational physics. Dubna, 2009.
Об авторах
Р. Л. Седов — ассист., КГТУ, e-mail: [email protected].
С. В. Орлов — канд. мед. наук, директор Института биомеханики позвоночника и крупных суставов, Калининград, e-mail: [email protected].
К. С. Латышев — д-р физ.-мат. наук, проф., РГУ им. И. Канта.
Authors
R. L. Sedov — assistant, KSTU, e-mail: [email protected].
Dr S. V. Orlov — director BSLAI, Kaliningrad, e-mail: [email protected].
Professor K. S. Latyshev — IKSUR.
