Журнал «Травма» Том 15, №6, 2014
Вернуться к номеру
Приближенный анализ анатомии, механических характеристик и напряженно-деформированного состояния позвоночника человека
Авторы: Чуйко А.Н. - Украинский аналитически-исследовательский центр в области медицины М. Угрина, г. Львов
Рубрики: Травматология и ортопедия
Разделы: Справочник специалиста
Версия для печати
Статья опубликована на с. 100-109
Введение
Вопросы клиники и биомеханики позвоночника волнуют человечество с древнейших времен. Еще в 1929 г. Н.А. Бернштейн в статье «Клинические пути современной биомеханики» наметил основные задачи и перспективы биомеханических исследований в клинике, показав их значимость для понимания вопросов патогенеза различных нарушений опорно-двигательного аппарата. Несколько позже появилось крылатое высказывание М.И. Ситенко: «Биомеханика — философия ортопедического мышления».
Все публикации по биомеханике и оценке напряженно-деформированного состояния позвоночника условно можно разделить на три этапа.
Первый этап, с древнейших времен до 70-х годов прошлого столетия, характеризуется в основном изучением анатомии и анализом рентгенограмм [1, 2].
Второй этап связан с возникновением в середине 1970-х годов метода конечных элементов (МКЭ) и практически одновременно — компьютерной томографии (КТ). На смену обычным представлениям, подкрепленным рентгенографическими данными, приходит 3D (трехмерное) моделирование, построенное на базе компьютерной томографии, с последующим анализом по методу конечного элемента. Однако на протяжении длительного времени 3D-модели создавались на основе инструментального определения размеров, а свойства мягких и костных тканей задавались по статистике. В подтверждение этого положения сошлемся на некоторые работы харьковских ученых [3–5], подчеркнув, что работа [4] опубликована в 2011 г., а [5] — в 2012. В [4] отмечено «При моделировании функционального взаимодействия элементов опорно-двигательной системы мы сталкиваемся с тем, что механические свойства биологических тканей практически не изучены. Все известные исследования проводятся на трупном материале… Методов изучения механических свойств живых тканей пока не существует».
Третий этап основан на современных компьютерных технологиях, например комплексе MIMICS-ANSYS, который позволяет изучать анатомические и механические (прочностные и жесткостные характеристики) свойства отдельных частей позвоночника конкретного пациента на основе данных компьютерной томографии в режиме реального времени, т.е. фактически in vivo, с последующим конечно-элементным анализом [6–8].
Целью работы является демонстрация возможностей современных компьютерных технологий по изучению анатомических особенностей и механических свойств позвоночника конкретного пациента и анализу его напряженно-деформированного состояния.
Материалы и методы
Условие прочности и условие жесткости — основа биомеханического анализа
В основе биомеханического анализа позвоночника в норме, при любых патологических изменениях и при наличии элементов, используемых при реконструкции (имплантатов, аппаратов внешней фиксации и т.п.), в соответствии с [6] должно лежать условие прочности или аналогичное ему по механическому смыслу условие жесткости.
Эти условия взаимно связывают нагрузку — расчетную модель (анатомию, конструкцию) анализируемого элемента — свойства костных тканей (конструкционных материалов). Представляя эти условия в виде треугольника, можно определить любую из вершин треугольника при известных значениях двух других.
Так, например, если в рассматриваемом случае анатомия (расчетная модель) позвоночника пациента задана, т.е. является неизменной, то при известных свойствах костных тканей можно определить максимально допустимую нагрузку; либо при заданной нагрузке можно определить необходимые свойства костных и мягких тканей, которые эту нагрузку способны выдержать. Точность получаемых результатов, их приближенность к конкретному пациенту будут зависеть от точности определения основных механических характеристик: линейных размеров элементов позвоночника, его поперечных сечений и свойств мягких и костных тканей — в первую очередь модуля упругости и предела прочности.
Очевидно, что при гармоничном развитии и росте позвоночника условие прочности выполняется автоматически, соответственно основным законам эволюции, т.е. при увеличении веса человека (основной нагрузки для позвоночника) должно происходить увеличение объема костной ткани и улучшение ее механических свойств. Нарушение этой закономерности неизменно приведет к развитию различных патологий, начало которым может быть положено и в детском возрасте.
Условие прочности, которое лежит в основе «глубокого биомеханического анализа» [6], формулируется достаточно четко и просто: необходимо определить действующие напряжения Q и сравнить их с допускаемыми Qu (разрушающими, индекс «u» от английского слова ultimate — предельный), т.е. проверить справедливость неравенства Q ≤ Qu (1).
Фактически эта формула объединяет три взаимосвязанных понятия, очерченных выше (рис. 1). Основной смысл этого неравенства — действующие напряжения не должны превышать значения разрушающего (травмирующего) напряжения для каждой структурной составляющей позвоночника [6].
Нагрузка. Необходимо знать нагрузку для конкретного пациента в норме, при наличии заболевания и после реконструкции, например установления имплантатов, аппаратов внешней фиксации и т.п. Так как обычно оценка напряженно-деформированного состояния (НДС) проводится в линейной постановке, то результаты, полученные при конкретной нагрузке (весе пациента), могут быть интерполированы на другие расчетные случаи.
Расчетная модель (анатомия, конструкция) анализируемого элемента — это, собственно, и есть построение модели и ее рационализация в соответствии с поставленными исследовательскими задачами.
Следует всегда иметь в виду, что модель может дать только те результаты, которые предусмотрены в ее функционировании. Например, оценка опорных свойств позвоночника может быть проведена на сравнительно простой модели, учитывающей только упругие свойства элементов позвоночника, а для оценки амортизационной функции модель должна учитывать вязкоупругие свойства элементов, во всяком случае межпозвоночного диска. В статье [8] применительно к периодонтальной связке отмечено, что «амортизирующие свойства периодонта в статической задаче никак не проявляются. Поясним это на простом примере. Положите на весы кусочек резины с хорошими амортизирующими свойствами, а сверху поставьте гирьку. Весы покажут суммарный вес гири и резины, сила веса никуда не исчезает и не уменьшается. Совершенно по-другому будет работать эта система при динамическом нагружении». В работе [9] сделана попытка рассмотреть ударное воздействие на зубочелюстной сегмент (ЗЧС), используя основы теории колебаний. Теория удара основана на теории колебаний и ее более сложном разделе — теории волн. Соответствующая постановка задачи и для позвоночника имеет, на наш взгляд, большие перспективы.
Далее под конструкцией мы будем понимать объект исследования, т.е. это может быть и позвоночный столб целиком, и любой из отделов позвоночника, и два смежных «полупозвонка» с межпозвоночным диском, и любой из этих элементов с имплантатом любого типа.
Свойства костных тканей. Если мы не знаем количественных прочностных характеристик кости конкретного пациента, то все наши рассуждения будут носить качественный характер — нужно больше, нужно меньше. Методика определения механических характеристик костных и мягких тканей подробно рассмотрена в работах [6, 7], и этому вопросу будет посвящена значительная часть предлагаемого исследования.
Отметим еще раз, что при известной нагрузке (силе) и конструкции (в первую очередь размерах) можно определить действующие напряжения (левая часть неравенства (1)), которые следует сравнивать с допускаемыми (травмирующими) напряжениями (правая часть неравенства (1)). Обе величины, входящие в это неравенство, являются как бы двумя сторонами одной медали. С одной стороны, необходимо постоянно совершенствовать методы определения действующих напряжений — добиваться максимальной корректности расчетной схемы: геометрических размеров, механических свойств структурных составляющих (например, учета физической нелинейности — пластичности костных тканей или гиперупругости мягких тканей, методов расчета) и т.п. С другой стороны, необходимо постоянное пополнение базы данных разрушающих напряжений, которые могут быть получены в основном экспериментальным путем, как при патологоанатомических исследованиях, так и при испытаниях живых тканей в зависимости от пола, возраста, типа заболевания и пр. Кроме того, устанавливать для костных тканей показатели типа разрушающие (допускаемые) напряжения, как в металлах, не совсем информативно. Здесь более продуктивно ввести показатель травмирующие напряжения костной и мягкой ткани [6].
Условие прочности, математическая и механическая сущность, которого определена соотношением (1), справедливо только для осевой нагрузки при условии ее равномерного распределения по сечению. Для остальных случаев — сдвига, изгиба и пр. — оно будет иметь другой вид. В более сложных случаях будут возникать зоны концентрации напряжений, которые при современном уровне развития науки могут быть наиболее успешно проанализированы с помощью метода конечных элементов (МКЭ), реализуемого на ЭВМ.
Для учета одновременно всех компонентов поля напряжений — и нормальных, и касательных — существуют так называемые эквивалентные напряжения (напряжения по Мизесу). Напряжения по Мизесу, показатель, к которому мы будем часто обращаться, рассчитываются по известной формуле
и характеризуют общее напряженное состояние в точке. В формуле (2) через Q и T обозначены соответственно нормальные и касательные напряжения, а индексы при них — направления действия напряжений вдоль осей x, y и z. Современные программы, реализующие метод конечного элемента (МКЭ), рассчитывают эквивалентное напряжение в автоматическом режиме.
Отметим, что напряжения по Мизесу (эквивалентные), определяемые по формуле (2), по своему механическому смыслу предполагают сравнение полученной величины с пределом текучести материала. Поэтому здесь и далее мы фактически будем предполагать, что травмирующие напряжения и предел текучести — идентичные понятия.
Предварительный анализ анатомии и механических свойств позвоночника конкретного пациента
Для анализа анатомии и механических свойств позвоночника воспользуемся данными компьютерной томограммы пациентки К. (18.07.11) 62 лет, с онкозаболеванием, потребовавшим резекции верхней челюсти и установки протеза, проведенным плановым лечением с радиационной нагрузкой и выраженным остеопорозом костных тканей.
Визуализацию объектов позвоночника на этапе предварительного анализа проведем с использованием программы MIMICS. MIMICS (Materialises Interactive Medical Image Control System) — интерактивный программный пакет для визуализации и сегментации изображений, полученных томографией (КТ, микроКТ, МРТ и др.), и 3D-изображения объектов [10]. Пакет предоставляет пользователю широкий набор функций по преобразованию наборов изображений в 3D-объекты и подготовке этих объектов для различных областей применения.
На рис. 2. приведен общий вид двух отделов скелета позвоночника пациентки К. при оттенках серого (GV) = 1250–4095. 3D-изображение тазового отдела свидетельствует о существенном искривлении позвоночника в этом отделе.
На рис. 3, 4 представлены результаты анализа компьютерной томограммы грудного отдела пациентки при разных уровнях (порогах) фиксируемой плотности кости — GV или чисел Хаунсфилда (HU), т.е. GV = 1024 + HU. Сравнение 3D-изображений на рис. 2–4 показывает, что плотность костных тканей разных элементов позвоночника различна и принимаемое обычно по статистике [3, 5] одно значение для всех элементов позвоночника далеко от реальности.
Рассмотрим возможность определения механических характеристик мягких и костных тканей путем создания аналитических зависимостей между числами HU либо значениями GV, определяющими рентгенологическую плотность ткани в условных единицах, и физической (реальной) плотностью мягких и костных тканей и их механическими характеристиками — пределом прочности p и модулем упругости E. По предлагаемой методике вначале создается единая линейная зависимость между числами GV (HU) и физической плотностью p.
В программе MIMICS предопределены следующие пороги (уровни) плотности для разных структурных составляющих костных и мягких тканей взрослого человека (табл. 1 и гистограмма на рис. 5).
Программа MIMICS позволяет определить основные механические характеристики мягких и костных тканей фактически в режиме in vivo, по эмпирическим формулам, зависящим от типа рассматриваемой ткани. Для компактной кости Femur рекомендуются следующие формулы.
Для определения плотности p или DN (Density)
DN = –13,4 + 1017 GV, (3)
где GV (Grey Values) — значения серого на томограмме.
Значение модуля упругости E вычисляется по формуле
E = –388,8 + 5925 DN . (4)
Результаты расчетов по этим формулам во всем диапазоне значений серого приведены на рис. 5 — штрих-пунктирная линия.
При дальнейшем анализе для определения основных механических характеристик костных тканей (модуля упругости и предела прочности) в зависимости от их плотности будем использовать также эмпирические формулы [11]
E = 2195 p3 и Q = 60 p2, (5)
где p — плотность костной ткани, которая вычисляется в г/см3. В этом случае модуль упругости и напряжения имеет размерность мегапаскаль (MПa).
Результаты расчета по этим формулам и данным табл. 1 представлены на гистограмме (рис. 5) в виде графических зависимостей Q и Е. Как видим, при плотности кости p около 1,8 г/см3 эти кривые для Е пересекаются, т.е. значения по разным методикам совпадают.
Предварительный анализ показал, что диапазон GV = 1250–4095 не перекрывает значений плотности всех структур межпозвоночного диска. Поэтому проведено расширение этого диапазона до GV = 1125–4095 и существенное упрощение 3D-модели с сохранением только двух полупозвонков (Th6-Th7) с межпозвоночным диском. На рис. 6 показано фронтальное и осевое сечение этой упрощенной системы с размерами и значениями плотности кости GV (Mean), которые будут использованы для предварительных расчетов.
Представив поперечное сечение межпозвоночного диска в виде эллипса и используя обычные формулы сопротивления материалов [6], получим:
— площадь поперечного сечения
— напряжение от сжатия (1), приняв, следуя [5], силу F = 400 Н,
По данным рис. 6б, среднее значение плотности кости в сечении GVср = 1289. По формуле (3) получим:
DN = –13,4 + 1017GV = –13,4 + 1017•1289 = = 1,31•106 г/м3 или p = 1,31 г/см3.
По формуле (4) получим:
E = –388,8 + 5925DN = –388,8 + 5925•1,31•106 = = 7,76•109 Па = 7,76•103 МПа.
Осевое смещение торца системы
Эти сугубо оценочные данные будут использованы ниже для сравнения с результатами расчетов по МКЭ.
Результаты и их анализ
Определение механических характеристик костных тканей с использованием программного комплекса MIMICS-ANSYS
Рис. 2–4 и 6б показывают, что плотность костных тканей разных элементов позвоночника различна. Для дифференцированной оценки плотности костных тканей воспользуемся возможностями программного комплекса MIMICS-ANSYS [7, 10].
На рис. 7а представлена 3D-модель рассматриваемого фрагмента позвоночника. Для удовлетворения требованиям программы ANSYS проведено сглаживание средствами MIMICS, которое распространяется не только на наружную поверхность, но и, что особенно важно, на внутреннюю границу между компактной и губчатой костью (рис. 7б). Отметим, что любые программы, реализующие МКЭ, не любят острых углов и мелких элементов, так как размер сетки должен быть меньше любого характерного размера модели.
MIMICS генерирует на 3D-модели поверхностную сеть и после ее оптимизации экспортирует оптимизированный объект в ANSYS (рис. 7в). Этот шаг по оптимизации сети является ответственным и может быть достаточно трудоемким, определяющим возможность дальнейшего применения всей технологической цепочки. Оптимизированная поверхностная сеть применяется далее для генерирования объемной сети в ANSYS, где выбирается тип твердотельного элемента из библиотеки элементов ANSYS.
Далее объемная сеть, созданная в ANSYS, загружается обратно в MIMICS для назначения материала каждому элементу. MIMICS может рассчитывать свойства материалов как по значениям GV, так и по значениям HU, выбрав соответствующую опцию. На специальной гистограмме MIMICS показывает количество элементов в модели для каждого значения серого цвета (рис. 8). Можно ввести количество материалов, например 10, и весь диапазон значений серого, встречающийся в объемной сети, будет разделен на 10 равных интервалов, каждый из которых будет отображать один материал, окрашенный в определенный цвет (рис. 9).
После вычисления свойств материалов для каждого элемента конечноэлементной сети они окрашиваются в определенный цвет в соответствии с цветовой гаммой, заложенной в MIMICS, и сводятся в специальную таблицу редактором MIMICS (рис. 10). На рис. 11 приведены эти же значения, но в памяти программы ANSYS (только для групп 6–10) в размерности МПа и г/мм3.
На рис. 12 показана 3D-модель фрагмента позвоночника с элементами, окрашенными в соответствии с цветовой гаммой MIMICS, а на рис. 13 — с цветовой гаммой ANSYS.
Разрезы 3D-модели фронтальными плоскостями в ANSYS представлены на рис. 14.
Можно констатировать, что поставленная выше задача — выявить структуру межпозвоночного диска — решена не полностью, на наш взгляд, по двум причинам: 1) принятый размер сетки (рис. 7в и 13) слишком крупный; 2) явление остеопороза привело к выравниванию механических характеристик позвонков и диска. Требуется более тщательное повторное моделирование системы.
В заключение этого раздела отметим, что плотность костных тканей в разных зонах модели отличается в 2,15 раза (рис. 10), пропорционально будут изменяться и другие механические характеристики, что соответствующим образом будет сказываться при анализе работы позвоночника под нагрузкой.
Анализ напряженно-деформированного состояния позвоночника с использованием метода конечного элемента
Для анализа напряженно-деформированного состояния позвоночника с использованием метода конечного элемента рассмотрим тот же фрагмент позвоночника. На рис. 15 показан этот фрагмент позвоночника с граничными условиями — для нижнего торца перемещения всех 75 узлов равны нулю (UZ = 0 — голубые треугольники) и нагрузкой — к 20 узлам верхнего торца приложена осевая сила по FZ = 20 Н (красные стрелки), т.е. суммарная нагрузка, как и выше, F = 400 Н. Отметим, что красные стрелки распределены неравномерно (справа их больше), что позволяет моделировать нагрузку на позвонки при обычном изгибе позвоночного столба.
Вся нагрузка, приложенная к верхнему торцу, полностью реализуется на нижнем торце, что подтверждается суммарной реакцией — см. последнюю строку на рис. 16. Поэтому задача по определению напряженно-деформированного состояния позвоночника заключается в анализе распределения нагрузки между структурными составляющими фрагмента позвоночника или, говоря другими словами, выявлении зон концентрации напряжений.
На рис. 17–19 представлены результаты анализа напряженно-деформированного состояния позвоночника. На рис. 17 показано поле напряжений по Мизесу SEQ (2) для 3D-модели и его разрез, на рис. 18 — основная компонента напряжений по Мизесу, осевые напряжения SZ (1) и разрез, а на рис. 19 — поле суммарных перемещений для 3D-модели и его разрез.
Особо отметим, что сила обычно изображается отрезком со стрелкой, как на рис. 15. Однако сил в виде стрелок (так принято изображать векторные величины) в природе практически не существует. Так их изображают для удобства восприятия и анализа. В сопротивлении материалов и теории упругости широко используется принцип Сен-Венана, который гласит: «В точках тела, достаточно удаленных от мест приложения нагрузок, внутренние силы весьма мало зависят от конкретного способа приложения этих нагрузок». Этот принцип во многих случаях позволяет производить замену одной системы сил другой системой, статически эквивалентной, что может упростить расчет. Поэтому в рассматриваемом случае, исходя из принципа Сен-Венана, будем анализировать только НДС вдали от точек приложения сил на верхнем торце и возникающих реакций на нижнем торце, сравнивая их с результатами по формулам сопротивления материалов (см. выше).
Краткий анализ полученных полей напряжений и перемещений позволяет сделать следующие выводы:
1. Наибольшая концентрация напряжений возникает в межпозвоночном диске (рис. 17, 18).
2. Основной компонентой напряжений по Мизесу SEQ (2) в рассматриваемом случае являются осевые напряжения SZ (1). В соответствии с рис. 17 максимум напряжений по Мизесу SEQ лежит в диапазоне 1,11–1,66 МПа (голубое поле), а осевых напряжений SZ (рис. 18) — в диапазоне –1,66…–1,33 МПа (зеленое поле). Напряжения SZ при сжатии имеют знак «минус». Вычисленное выше по формуле сопротивления материалов значение напряжения 0,676 МПа предполагает, что сила по сечению распределена равномерно, в отличие от результатов по МКЭ, которое учитывает напряжения в каждой точке. Поэтому можно считать, что эти цифры коррелируются удовлетворительно.
3. Перемещения в обоих подходах вычисляются интегрально. Поэтому, учитывая неравномерность приложенной нагрузки и принцип Сен-Венана, перемещения, вычисленные по формуле сопротивления материалов 3,46•10–3 мм и по МКЭ (желтое поле на рис. 19) в диапазоне 0,00279–0,003263 мм, можно считать, совпадают удовлетворительно.
В заключение этого раздела отметим, что в предлагаемой работе отработан алгоритм определения действующих напряжений (левая часть неравенства (1)) — в зависимости от величины действующей нагрузки, ее распределения по элементам системы и механических свойств костных тканей. Вопрос о величине травмирующих (разрушающих) напряжений (правая часть неравенства (1)) остается открытым. Мы согласны с выводами авторов работы [4]: достоверных «методов изучения механических свойств живых тканей пока не существует». Поэтому необходимо постоянное пополнение базы данных разрушающих напряжений, которые могут быть получены в основном экспериментальным путем, как при патологоанатомических исследованиях, так и при испытаниях живых тканей в зависимости от пола, возраста, типа заболевания и пр. [6].
Заключение
Разработан алгоритм (методика) использования программного комплекса MIMICS-ANSYS для уточнения анатомических особенностей позвоночника конкретного пациента, количественного определения основных механических характеристик мягких и костных тканей, анализа параметров напряженно-деформированного состояния в режиме реального времени, т.е. фактически in vivo.
Намечены направления уточнения предлагаемой методики.
Показана необходимость постоянного пополнения базы данных разрушающих напряжений, которые могут быть получены в основном экспериментальным путем, как при патологоанатомических исследованиях, так и при испытаниях живых тканей в зависимости от пола, возраста, типа заболевания и пр.
1. Воробьев В.П. Анатомия человека. — М.: Госмедиздат, 1932. — Т. 1.
2. Михайлов С.С. Анатомия человека. — М.: Медицина, 1984. — 704 с., илл.
3. Tkachuk N.A., Veretelnyk Y.V., Tkachuk N.N. Generalized parametrical approach to research of biomechanical systems elements. The International Conference «Аdvanced information and telemedicine technologies for healh». — Minsk, 2005. — P. 63-67.
4. Тяжелов А.А., Карпинский М.Ю., Карпинская Е.Д. Проблемы биомеханики в ортопедии. «Современные проблемы математики и ее приложения в естественных науках и информационных технологиях». Международная конференция «Тараповские чтения». — Харьков, 2011. — С. 104-105.
5. Стауде В.А., Кондратьев А.В., Карпинский М.Ю. Численное моделирование и анализ напряженно-деформированного состояния крестцово-подвздошного сочленения при различных вариантах поясничного лордоза // Ортопедия, травматология и протезирование. — 2012. — № 2. — С. 50-56.
6. Чуйко А.Н., Шинчуковский И.А. Биомеханика в стоматологии: монография. — Х.: Форт. 2010. — 516 с.
7. Чуйко А.Н., Калиновский Д.К., Левандовский Р.А., Грибов Д.А. Биомеханическое сопровождение операций в челюстно-лицевой хирургии с использованием комплекса MIMICS-ANSYS // Ортопедия, травматология и протезирование. — 2012. — № 2. — С. 57-63.
8. Чуйко А.Н. Еще раз о биомеханике пародонта. Часть 2 // Пародонтология. — 2007. — № 4. — С. 45-52.
9. Чуйко А.Н., Вовк В.Е. Особенности биомеханики в стоматологии. — Х.: Прапор, 2006. — 304 с.
10. Mimics 12. Пакет обработки изображений. Базовый обучающий курс // Materialise. — 2008. — С. 81.
11. Mow C., Hayes W.C. Basic Orthopedic Biomechanics. — New York, 1991.