Содержание к диссертации
Введение
1 Математическая модель 14
1.1 Построение математической модели механики черепно-мозговой травмы 14
1.1.1 Анатомическое строение головы человека 14
1.1.2 Двумерные модели механической реакции головы человека на ударные нагрузки 19
1.2 Уравнения механики деформируемого твердого тела 22
1.2.1 Уравнения движения и реологические соотношения 22
1.2.2 Матричная форма уравнений 23
1.2.3 Гиперболические свойства систем уравнений линейной упругости 24
2 Численные методы 25
2.1 Численные методы для решения системы уравнений механики деформируемого твердого тела 25
Оглавление З
2.1.1 Конечно-разностные методы для одномерных гиперболических систем уравнений 25
2.1.2 Разностные схемы для двумерных гиперболических уравнений 29
2.1.3 Разложение матриц для двумерных уравнений линейной упругости 30
2.1.4 Разностные схемы для неструктурированных треугольных сеток 34
2.2 Расчет граничных узлов и постановка граничных условий 38
2.2.1 Расчет граничных узлов в одномерном случае 38
2.2.2 Расчет граничных узлов в двумерном случае 40
2.2.3 Примеры граничных условий 43
2.3 Решение задачи контактного разрыва 46
2.3.1 Решение задачи контактного разрыва в одномерном случае . 46
2.3.2 Решение задачи контактного разрыва в двумерном случае . 47
2.3.3 Алгоритм обнаружения столкновений (collision detection) . 49
2.3.4 Примеры контактных условий 55
2.4 Перестройка расчетных сеток 63
2.4.1 Структурированные четырехугольные сетки 63
2.4.2 Неструктурированные треугольные сетки 67
3 Результаты моделирования 69
3.1 Интегральная оценка механического воздействия 69
3.2 Результаты расчетов 72
3.2.1 Качественное описание механической реакции головы человека на ударную нагрузку 72
3.2.2 Сравнение с клиническими данными 77
Заключение 81
Основные результаты и выводы диссертации 81
Литература
- Двумерные модели механической реакции головы человека на ударные нагрузки
- Матричная форма уравнений
- Разностные схемы для двумерных гиперболических уравнений
- Качественное описание механической реакции головы человека на ударную нагрузку
Введение к работе
Постановка задачи о моделировании механических факторов черепно-мозговой травмы
Одной из активно развивающихся областей приложения вычислительной механики деформируемого твердого тела является моделирование процессов, происходящих в биомеханических системах. Изучение механических аспектов функционирования различных биосистем, в частности, органов человека, позволяет получить новые качественные и количественные характеристики систем, выработать новые принципы диагностики на ранних стадиях различных заболеваний, уточнить существующие критерии безопасности в различных областях человеческой деятельности.
В диссертационной работе рассматривается задача о численном моделировании механической реакции системы череп-мозг на ударные воздействия. Эта задача является актуальной с точки зрения выяснения механизмов повреждаемости тканей мозга при различных типах внешней нагрузки.
Математическое моделирование как метод исследования в данном случае обладает как рядом очевидных достоинств (малая стоимость вычислительного эксперимента, легкость и широта диапазона изменения основных параметров, полнота получаемой в результате картины динамики скоростей, напряжений и деформаций во всем объеме), так и рядом очевидных недостатков (существенное упрощение и сужение спектра рассматриваемых процессов, большая стоимость разработки необходимого программного обеспечения, сложность или принципиальная невозможность количественного сопоставления данных вычислительного эксперимента и натурных испытаний).
Оглавление
Тем не менее, механико-математические модели реакции головы человека на ударные нагрузки широко используются как инструмент для исследования качественных зависимостей влияния геометрических параметров (например, размера головы, места приложения удара, возрастной атрофии мозга) на степень риска при черепно-мозговых травмах.
Исследованию этой задачи посвящен ряд зарубежных работ, в которых основным инструментом является метод конечных элементов (МКЭ), позволяющий решать системы дифференциальных уравнений в двумерных и трехмерных областях произвольной формы. В нашей стране, по результатам анализа литературы и электронных публикаций, подобные работы не ведутся.
Научные результаты диссертации опубликованы в работах [1—11].
Цель работы Целью работы является получение качественной картины распределения механических нагрузок на мозг человека при различных ударных нагрузках на основе двумерных механико-математических моделей системы череп-мозг с применением сеточно-характеристических численных методов на регулярных и неструктурированных сетках.
Для достижения поставленной цели возникла необходимость в разработке и реализации адаптированных сеточно-характеристических численных методов, пригодных для решения двумерных динамических задач при наличии нескольких взаимодействующих деформируемых тел произвольной формы с неоднородной структурой, претерпевающих конечные деформации.
Научная новизна. Научная новизна работы заключается в следующем:
Реализован метод обобщения известных сеточно-характеристических схем на неструктурированные треугольные сетки путем восстановления регулярного сеточного шаблона интерполяцией значений в вершинах треугольной сетки.
Разработан алгоритм определения контактных границ между взаимодействующими деформируемыми телами, учитывающий непрерывность их движения.
Алгоритмы перестройки регулярной сетки для адаптации к деформирующимся границам области интегрирования обобщены на случай неструктури-
Оглавление
рованных треугольных сеток.
4. Реализованные сеточно-характеристические методы на неструктурированных треугольных сетках использованы для моделирования реакции головы человека на ударные нагрузки.
Практическая ценность. Построенная математическая модель воздействия ударных нагрузок на голову человека может быть использована как инструмент исследования различных механизмов повреждения мозга при черепно-мозговых травмах. Численный эксперимент может служить источником данных для построения экспертных систем диагностики и прогнозирования состояния больных с черепно-мозговой травмой. Полное численное моделирование механической реакции головы человека на ударную нагрузку может использоваться как вспомогательное средство в судебной медицине для проверки гипотез о природе тех или иных телесных повреждений.
Разработанный в рамках исследования комплекс программ может применяться для решения целого ряда задач механики деформируемого твердого тела, характеризующихся ударным характером механических нагрузок, сложной геометрией областей интегрирования, наличием динамических контактных взаимодействий между деформируемыми телами.
Оглавление
Анализ имеющихся работ по моделированию ударного воздействия на голову человека
По результатам анализа работ, опубликованных в печатных и электронных изданиях, можно сделать вывод о том, что подавляющее количество исследовательских работ в области моделирования ЧМТ было проведено за рубежом с применением различных вариаций методов конечных элементов (МКЭ).
Несмотря на то, что моделирование головного отдела человека с помощью МКЭ в последние десять лет интенсивно развивается, оно еще далеко от того, чтобы иметь возможность объяснить механизмы повреждения мозга и предсказать последствия черепно-мозговой травмы. Такие вопросы, как реология биоматериалов, составляющих голову человека, механическое взаимодействие поверхностей мозга и черепа, воздействие желудочков и пустот на распределение напряжений, рассмотрение многослойного строения оболочек мозга с воспроизведением его точной геометрии требуют дальнейшего исследования.
Для большинства случаев механического воздействия может быль адекватна лишь трехмерная модель черепно-мозгового отдела. Тем не менее, двумерные модели полезны в качестве начального приближения. Они также упрощают введение геометрических неоднородностей, которые могут моделироваться отдельно.
В работе [12] был предложен ряд двумерных моделей для трансверсально-го и сагиттального сечений для анализа влияния различий в механических свойствах белого и серого вещества на распределение максимальных сдвиговых напряжений. Автор показал, что различие в механических свойствах серого и белого веществ, а также включение в модель желудочков необходимо для соответствия вычисленного распределения максимальных сдвиговых напряжений наблюдаемым аксональным повреждениям в белом веществе. В данной модели желудочки описывались жидкостными элементами, а модуль Юнга белого вещества был на 60% больше, чем серого.
В исследовании [13] использовалась двумерная модель вертикального сечения головы человека. Контакт череп-мозг моделировался двумя способами: 1) жесткое сцепление, и 2) скольжение с возможностью отрыва. Кроме того, боль-
Оглавление 8
Таблица 1. Сводная таблица двумерных конечно-элементных моделей головы человека
шое отверстие также моделировалось тремя способами: 1) жесткое сцепление, 2) скольжение с возможностью отрыва, и 3) незакрепленное. Сравнение с экспериментом [14] показало наилучшее соответствие величины давления в месте удара в случае скользящего контакта, в то время как модели с жестким сцеплением показали результат, наиболее далекий от действительности. Ни одна из моделей не показала соответствия давлений в месте противоудара, что может объясняться неспособностью контактного алгоритма противостоять растягивающим нагрузкам. Авторы заключили, что величины давления как в месте удара, так и в месте противоудара намного более чувствительны к способу моделирования контакта череп-мозг, чем к наличию или отсутствию свободного большого отверстия.
Оглавление
В работе [21] исследуются сдвиговые характеристики мозгового вещества на основе анализа экспериментальных данных, полученных с помощью физической модели. Показано, что вязкоупругая реология, использованная для моделирования мозгового вещества в большинстве работ [12,22,23], справедлива лишь для деформаций, не превышающих 1%. Предложена нелинейная вязко-упругая модель на основе уравнений гиперупругости Муни-Ривлина. Нелинейная модель мозгового вещества применена для конечно-элементной модели из [22], и было показано улучшение соответствия линейных смещений мозгового вещества при фронтальном ударе с экспериментальными данными [14].
Смещение мозга относительно черепа описаны более пятидесяти лет назад. Это вызывало особый интерес к численному моделированию этого явления. Моделирование внутричерепной жидкости линейно-упругими элементами с малым сдвиговым модулем [12,24] показало, что этот подход вносит дополнительные касательные напряжения между черепом и мозгом, а также снижают устойчивость вычислений вследствие вырождения элементов. Другим подходом является введение контактных условий на границе череп-мозг. Конкретные условия различаются у разных авторов, от жесткого сцепления до скольжения с трением. Было проведено некоторое количество сравнительных исследований с использованием различных контактных условий [22]. Все они приводят к заключению, что эффект ударного воздействия на голову человека чувствителен к способу моделирования взаимодействия череп-мозг.
Более точная трехмерная геометрическая модель черепа была использована в [23] для оценки существующих критериев повреждаемости мозга. Модель включает в себя мембраны твердой оболочки, сагиттальные синусы и эмиссар-ные вен. Численное исследование проводилось с помощью пакета нелинейного конечно-элементного анализа LS-DYNA3D. Автором предложены способы коррекции критерия НІС (Head Injury Criteria), применяющегося в машиностроении, в зависимости от размера головы, размера мозга и направления удара.
Основные трудности, возникающие в большинстве современных МКЭ моделей ударных воздействий на черепно-мозговой отдел человека, возникают в связи с высокой чувствительностью результата к аккуратному моделированию контактных взаимодействий черепа и мозга, а также к учету различных структурных
Оглавление
неоднородностей в мозге. Это дает основание полагать, что применение иных классов численных методов может служить эффективным средством повышения точности существующих моделей. Одним из таких классов является класс сеточно-характеристических методов, известных более аккуратной формулировкой граничных и контактных условий, а также способностью более адекватно описывать сложные волновые картины распространения возмущений в сильно гетерогенной среде. В то же время наличие сложной геометрии свидетельствует о необходимости обобщения и распространения существующих сеточно-характеристических методов на случай нерегулярных сеток.
Оглавление
Содержание работы
В диссертационной работе рассматриваются постановка задачи о моделировании механических факторов черепно-мозговой травмы, программная реализация сеточно-характеристических методов на регулярных и неструктурированных треугольных сетках, использованная для ее решения, проведение тестовых расчетов для отладки программной реализации и анализ результатов, полученных при расчетах.
Во введении обсуждается актуальность работы, приводится обзор работ на тему моделирования ударного воздействия на голову человека, выполненных с помощью методов конечных элементов в двумерной и трехмерной постановке. Оцениваются перспективы применения для решения этих задач сеточно-характеристических методов.
В главе 1 рассматривается анатомическое строение головы человека, уравнения механики деформируемого твердого тела и основные реологические модели: линейно-упругая и вязкоупругая, применяемые для построения математической модели черепно-мозгового отдела, гиперболические свойства систем уравнений, описывающих волновые процессы. Описаны двумерные модели механического воздействия на голову человека, использованные в данной работе.
В главе 2 обсуждаются реализованные подходы к численному решению уравнений гиперболического типа, приведены примеры конечно-разностных методов, относящихся к классу сеточно-характеристических. Приводятся результаты сравнительного анализа нескольких разностных схем для решения двумерных нестационарных задач с разрывными начальными условиями. Рассматривается важный вопрос об аппроксимации граничных и контактных условий, технике расчета граничных узлов расчетной сетки и контактных границ между двумя деформируемыми телами.
Рассматривается предложенный автором и реализованный программно метод обобщения ряда известных схем для регулярных четырехугольных расчетных сеток на случай неструктурированных треугольных сеток. На примере численного решения тестовых задач показывается его эффективность. Обсуждаются недостатки и достоинства применения неструктурированных сеток к численному моделированию черепно-мозговой травмы.
Оглавление
Обсуждаются методы построения адаптивных расчетных сеток и их перестройки в процессе расчета.
В главе 3 приведены результаты расчетов в зависимости от структуры модели, типа и места приложения механической нагрузки и контактных условий на границе череп-мозг, а также итоги их сравнения с клиническими данными. На основании полученных численно результатов делаются выводы относительно характера распределения механических нагрузок на мозг при ударном воздействии.
В заключении приводятся основные результаты и выводы, полученные автором в ходе выполнения работы по теме диссертации.
Оглавление
Положения, выносимые на защиту
На защиту выносятся:
Математическая модель механической реакции системы череп-мозг на ударное воздействие для численного расчета пространственного распределения механических нагрузок и деформаций тканей мозга при черепно-мозговой травме.
Разработанный в рамках исследования комплекс программ с использованием сеточно-характеристических методов на двумерных регулярных и неструктурированных расчетных сетках для численного решения задач механики деформируемого твердого тела с областями интегрирования произвольной формы, наличием динамических контактных взаимодействий и конечными деформациями.
Реализованный алгоритм динамического определения и расчета контактных границ между деформируемыми телами произвольной формы, учитывающий непрерывность их движения
Обобщение метода перестройки регулярных сеток для адаптации к деформирующимся границам на неструктурированные треугольные сетки.
Благодарности
Автор выражает благодарность академику Олегу Михайловичу Белоцер-ковскому за постоянное внимание, проявленное к данной работе, коллективу сотрудников Главного военного клинического госпиталя им. Н. Н. Бурденко и лично Игорю Александровичу Климову за ценные консультации и поддержку, и научному руководителю Игорю Борисовичу Петрову за постоянное вдохновение и совместную научную работу на протяжении нескольких лет.
Двумерные модели механической реакции головы человека на ударные нагрузки
В данной работе использовались двумерные модели системы череп-мозг, построенные для трансверсального сечения. Простейшей используемой моделью являлась двухкомпонентная модель (рис. 1.6, а). Ткани кости и мозга описывались однородными изотропными материалами, имеющими усредненные механические свойства. Внешняя нагрузка задавалась в виде соударения системы череп-мозг с абсолютно жесткой неподвижной преградой с заданной начальной скоростью (1-3 м/с).
Согласно результатам [12,22,23], большой вклад в формирование волновой картины распространения возмущений в упругой среде вносят неоднородно сти. Поэтому двухкомпонентная модель получила дальнейшее развитие в виде трехкомпонентной модели (рис. 1.6, б) , включающей, помимо кости и мозга, еще желудочки.
Мембраны твердой мозговой оболочки (dura mater) оказывают сдерживающее влияние на перемещения мозга внутри черепной коробки. Поэтому впоследствии было принято решение ввести в модель falx cerebri (рис. 1.6, в), вертикальную мембрану, разделяющую полушария в теменной области.
Реологические свойства биоматериалов также подвергались вариации. Так, реология мозгового вещества менялась от идеальноупругой до вязкоупругой. Кость моделировалась идеальноупрутим материалом со средними свойствами пластинчатой и губчатой кости. Желудочки представлялись квази-жидкостью - вязко-упругим материалом с модулем сдвига, близким к нулю.
Моделирование взаимодействия между черепом и мозгом является сложной задачей ввиду того, что в действительности мозг имеет большое количество различных по механическим свойствам (многие из которых до сих пор неизвестны) оболочек, складчатых структур, врастающих друг в друга, с полостями, заполненными ликвором (цереброспинальной жидкостью). В виду невозможности имитировать детально структуру контактной поверхности, общепринятый подход состоит в замене ее на геометрически простую контактную поверхность со сложными контактными условиями, которые подбираются экспериментально для получения качественного или количественного соответствия между расчетными и экспериментальными данными.
В данной работе применялся метод явного выделения контактного разрыва с контактными условиями, которые варьировались от полного слипания до скольжения с возможностью отслоения.
Механические свойства материалов и их реология, использованные в описанных двумерных моделях, приводятся в таблице 1.1. Осредненные упругие константы для костей черепа и мозга заимствованы из [25].
Для математического моделирования волновых процессов в деформируемом твердом теле использовалась система динамических уравнений [26,27] в виде рщ = Vj оц + fi (уравнения движения), Уц — Qijki Ski + Fij (реологические соотношения).
Здесь р - ПЛОТНОСТЬ СреДЫ, Vi - КОМПОНенТЫ СКОРОСТИ СМещеНИЯ, Jij, Sij компоненты тензоров напряжений и деформаций, Vj - ковариантная производная по j -й координате, /j - массовые силы, действующие на единицу объема, F -добавочная правая часть.
В случае малых деформаций тензор скоростей деформаций е - = ёу выражается через компоненты скорости смещения линейным образом:
Вид компонент тензора 4-го порядка q i определяется реологией среды. Для линейно-упругого тела Qijki = 5кі + p (SikSji + SuSjk), (1.2) Fij = 0. В этом соотношении, которое обобщает закон Гука, Аид- параметры Ляме, 5у- -символ Кронекера. Для вязкоупругого тела правая часть Fy- определяет временное затухание напряжений в среде: т?. __ У rij — г , где г - время релаксации. Для замыкания системы уравнений (1.1) ее необходимо дополнить уравнением состояния, определяющим зависимость плотности от напряжений: ж. р = р0ек, где р — — 2 akk - давление, К = А + /І - коэффициент всестороннего сжатия. Матричная форма уравнений Уравнения (1.1, 1.2) допускают запись в матричной форме: (1.3)
Здесь й = {vx, vy, тхх, аху, ауу, czz] - вектор искомых функций, / - вектор правых частей той же размерности, х, у - независимые пространственные переменные, t -время, Л-х — При использовании в численных расчетах подвижных регулярных расчетных сеток, например, лагранжевых, удобно записывать уравнения (1.3) в системе координат (, 7?), связанной с узлами расчетной сетки. Ее связь с фиксированной декартовой системой координат (х, у) фактически задает уравнения движения узлов сетки: t t х{, v, t) = Хо(, г))+ сх(, г], t)dt, у(, г), t) = Y0(, 77)+/ Су(, V, t)dt, о о где функции {Xo(,rj),Yo(,r))} определяют начальную конфигурацию расчетной сетки при t = О, а с(, 77) = {сх, Су} - скорость движения узлов сетки в декартовой системе координат (для лагранжевой сетки с = {vx,vy}).
Матричная форма уравнений
В случае двумерных гиперболических уравнений разностные схемы на регулярной расчетной сетке {xi = lh, ут = mh, tn = пт} могут быть получены из перечисленных одномерных разностных схем следующим образом. Каждую явную схему на трехточечном шаблоне можно представить в виде действия оператора
Тогда соответствующая двумерная схема на пятиточечном шаблоне будет иметь вид аК1 = Я + ос д(±Ах, т, /i)K_1 m, дцт, йг+1,т}+ + (і _ а) g{ Ay,T,h){ m_ m m+l}. В случае линейного оператора Q(A, т, К) получаемая разностная схема не зависит от коэффициента а. Например, двумерная схема Куранта-Изаксона-Риса на пятиточечном шаблоне имеет вид: l,m V?,m Т х Лх і(иГ+і,т — Uttm) Тр х Лг ( — "Г-1,т) П П (2.5) - А-Пу( т+1 - Щ - - Й ). Максимальный шаг по времени для схемы (2.5) достигается при maxlAJ а = maxAz -ЬтахІЛуІ и равен h т max АХ + max \ХУ\
Видно, что для регулярной прямоугольной сетки шаг интегрирования уменьшается вдвое по сравнению с одномерным случаем.
Расщепление по пространственным направлениям. Идея другого метода заключается в замене двумерной системы (1.3) тремя одномерными системами
Решение их производится с помощью некоторой одномерной разностной схемы поэтапно: сначала во всей области интегрирования решается система (2.6), затем результат используется для решения системы (2.7), и, наконец, конечное решение получается в результате прибавления правой части (2.8).
Максимальный шаг интегрирования для каждого этапа при использовании данного подхода совпадает с соответствующим шагом для одномерной схемы.
Для сравнения двух методов построения двумерных разностных схем была численно решена тестовая задача о «точечном взрыве», т.е. о распространении первоначального возмущения, заданного в круговой области.
Ниже приведены изолинии модуля скорости, полученные численно с применением перечисленных разностных схем: без расщепления (рис. 2.3) и с расщеплением (рис. 2.4). Вычисления проводились на прямоугольной сетке размером 200 х 200 узлов.
Как видно из рисунка, для всех схем, кроме схемы Лакса, метод с расщеплением дает лучший результат: численное решение обладает меньшей степенью размытия волнового фронта и большей изотропностью (меньшей разницей между решениями вдоль сеточных линий и вдоль диагональных направлений), чем при использовании метода без расщепления. Отличительным признаком всех схем с расщеплением является отсутствие полной симметрии между х- и у-направлениями, что, однако, не сказывается на общем качестве результата. Анализ свойств различных двумерных разностных схем для гиперболических уравнений приведен в [4].
Разложение матриц для двумерных уравнений линейной упругости Для построения перечисленных выше сеточно-характеристических схем требуется знать разложение матриц системы (1.3)
При расчетах на подвижных расчетных сетках (см. раздел 1.2.2) необходимо иметь разложение матрицы более общего вида, соответствующей переходу в произвольную криволинейную систему координат, движущуюся относительно - компоненты скорости узла расчетной сетки относительно точек сплошной среды, выраженные в системе координат { ,г]}. Для лагранжевой сетки сх = vx, Су = vy, поэтому AV — Avv = 0.
Набор собственных векторов не зависит от первого слагаемого в (2.9) и, следовательно, от Jt. Матрица левых собственных векторов П( Jx, Jy) равна
Трудности с построением разностных схем на неструктурированных сетках связаны с тем, что у расчетного узла в такой сетке не имеется фиксированных соседей по координатным направлениям. Для аппроксимации пространственных производных используются значения в произвольно расположенных соседних узлах. Степень произвольности зависит от топологической «структуры», описывающей неструктурированную сетку. Например, в случае триангуляции Делоне, использовавшейся в данной работе, на соседство узлов накладываются определенные ограничения, гарантирующие для любого треугольника отсутствие других узлов внутри описанной вокруг него окружности. Таким образом, с каждым расчетным узлом связан набор треугольников, содержащих его в качестве вершины, которые можно рассматривать как «неструктурированный шаблон».
Разностные схемы для двумерных гиперболических уравнений
Схемы без расщепления. Описанный выше алгоритм расчета граничных узлов легко обобщается на сеточно-характеристические схемы для двумерных уравнений, построенные без расщепления по пространственным направлениям. В таких схемах приращение искомой функции состоит из двух слагаемых, каждое из которых вычисляется по своему направлению независимо от другого. Например, двумерную схему Куранта-Изаксона-Риса (2.5) можно переписать в виде гС1 = tC + 8& + «, где
Заметим, что при использовании нелагранжевой расчетной сетки количество выводящих характеристик для внешнего направления может отличаться от двух: это зависит от знака и величины компонент Дг и Avv скорости движения узла сетки относительно точек сплошной среды (см. раздел 2.1.3).
Допустим, {, г]} является локальной системой координат граничного узла нелагранжевой расчетной сетки (рис. 2.9). Для того, чтобы в данном узле 5чй, можно было вычислить с привлечением двух граничных условий, мы должны потребовать Avr, = (vx - сх)— + (vy - су) — Ov дх = (Vx - CX)(-QC) + (Vy - Су) — т = {v — с) п = 0.
Вектор (v — с) есть скорость точек сплошной среды относительно узла сетки, an- направление нормали к области интегрирования в данной точке.
Таким образом, требование наличия только двух выводящих характеристик эквивалентно вполне естественному условию на движение узлов сетки: граничные узлы должны двигаться относительно точек сплошной среды по касательной к границе области интегрирования.
Рассмотрим случай регулярной четырехугольной сетки. Построение алгоритма расчета для бокового узла (см. рис. 2.10, а) не вызывает трудностей: приращение по оси вычисляется так же, как и для внутреннего узла, с помощью значений в узлах (/ — 1, М), (I, М), (1 + 1, М), а приращение по оси ц вычисляется с помощью алгоритма, аналогичного описанному в предыдущем разделе, с привлечением значений в узлах (I, М — 1), (ly М) и двух граничных условий для верхней границы области интегрирования.
Более сложный случай представляет собой угловой узел (см. рис. 2.10, б): использование описанного алгоритма для обоих направлений и ц потребовало бы постановки четырех граничных условий, что не является физически корректным. Даже в случае совместности граничных условий на верхней и правой границе области интегрирования, подобный алгоритм является численно неустойчивым.
Подход, применяемый в данной работе, заключается в том, что угловая точка сетки рассматривается как «гладкая» граничная точка области интегрирования (см. рис. 2.10, в): ось f направляется по касательной к границе тела, а ось т) - через диагональный узел (L — 1, М — 1). При этом соседними узлами по оси считаются узлы (L - 1, М) и (L,M — 1), а единственным соседом по оси ц -узел (L - 1,М - 1). Далее угловой узел рассчитывается так же, как боковой, с привлечением двух граничных условий для объединенной «гладкой границы.
Описанный подход может рассматриваться как способ регуляризации численного решения задач упругости, имеющих особенности в области абсолютно острых углов. Допущение о «гладкости» границы области интегрирования как правило имеют место в большинстве практических задач.
Схемы с расщеплением. Для использования двумерных сеточно-характерис-тических схем, построенных методом расщепления по пространственным направлениям, в данной работе был реализован комбинированный алгоритм, в котором часть узлов - там, где это возможно, - рассчитывается по схеме с расщеплением, а другая часть - по соответствующе]"] двумерной схеме без расщепления.
Рассмотрим опять случай регулярной четырехугольной сетки (рис. 2.11, а). Комбинированный алгоритм в каждом узле хранит флаг, соответствующий этапу расщепления.
Как видно из рисунка, в случае регулярной сетки на 3 этапе рассчитываются только граничные узлы. В случае неструктурированной используется аналогичный способ маркировки узлов, только условие на соседние узлы распространяется на вершины всех треугольников, по которыми производится интерполяция значений для виртуального шаблона. Поэтому конфигурация узлов может быть более сложной, и на 3 этапе могут рассчитываться с помощью двумерной схемы также и внутренние узлы.
При использовании комбинированного алгоритма подавляющее большинство узлов рассчитываются по схемам с расщеплением по пространственным направлениям (так как являются внутренними), что позволяет сохранить все полезные качества этих схем (меньшее размазывание волнового фронта и большую изотропность численного решения).
Качественное описание механической реакции головы человека на ударную нагрузку
В этом разделе приведены некоторые результаты численного решения задачи моделирования механической реакции головы человека на ударную нагрузку.
На рис. 3.2 приведены изолинии давления для бокового соударения головы человека, движущейся с начальной скоростью VQ = 1 м/с с неподвижной плоскостью слева. Продолжительность соударения порядка 10 мс.
На рис. 3.3, 3.4 и 3.5 показаны зависимости интегральных параметров механического воздействия - максимального сжатия, максимального растяжения и максимальных сдвиговых нагрузок, - от направления удара.
По результатам анализа пространственных распределений пиковых механических воздействий на ткани мозга были сделаны следующие выводы:
1. Численно получены области концентрации растяжения и сжатия в местах удара и противоудара, вне зависимости от варианта модели. Максимальные сжатия в области удара и максимальные растяжения в области противоудара в 2 - 5 раз превосходят соответствующие величины на противоположной стороне. Сопоставляя эти данные с описанным в литературе противоударным механизмом ушибов мозга, можно заключить, что растягивающие напряжения являются более сильным фактором поражения мозговой ткани, чем сжимающие.
2. Области локализации растяжения и сжатия имеют примерно эллиптическую форму с одним главным максимумом, лежащим близко или на поверхности мозга. Сдвиговые напряжения имеют множество фокусов концентрации, которые более локальны и располагаются внутри объема мозга. Поэтому можно заключить, что поражения внутренних областей мозга, не прилежащих к его поверхности, обусловливаются преимущественно сдвиговыми нагрузками.
3. Характер распределения максимальных сжатий и растяжений относительно оси удара практически не изменяется для моделей без мембраны.
4. Характер распределения сдвиговых нагрузок, равно как и их максимальные значения, при изменении направления удара меняются нерегулярным образом. На основании этого можно предположить, что их пространственное распределение существенным образом определяется геометрией черепа (т.е. отличием его формы от сферической), что затрудняет предсказание конкретных мест их концентрации при помощи двумерной модели.
5. Для модели с мембраной сжатия и растяжения в противоударной области при боковом ударе уменьшаются примерно вдвое по сравнению с моделями без мембраны. Это можно объяснить экранирующим действием мембраны, отражающей продольные упругие волны и сдерживающей смещения мозгового вещества.
6. Сдвиговые нагрузки для модели с мембраной увеличиваются в 1.25-1.5 раз по сравнению с моделями без мембраны. Это можно объяснить более сложной картиной деформации мозгового вещества вокруг отростков мембраны при их «обтекании».
7. Характер контактного взаимодействия на границе череп-мозг оказывает существенное влияние на распределение сдвиговых нагрузок (см. рис. 3.6, а-б).
8. Наличие в модели желудочков практически не оказывает влияния на распределение в пространстве областей сжатия/растяжения, но существенно сказывается на распределении сдвиговых напряжений (см. рис. 3.6, б-в). Сдвиговые напряжения обтекают область желудочков, поэтому количество фокусов их концентрации увеличивается.
Для верификации модели было проведено сравнение полученных результатов с клиническими данными, предоставленными коллективом сотрудников отделения неврологии Главного военного клинического госпиталя им. Н. Н. Бурденко (руководитель коллектива И. А. Климов) в рамках совместного проекта при поддержке Российско-Индийского Центра компьютерных исследований (РИЦКИ).
В таблице 3.1 приведены некоторые данные по пациентам с черепно-мозговыми травмами различной тяжести. На снимках компьютерной томографии более свет
лые пятна соответствуют повышенной плотности вещества в местах кровоизлияний, более темные - участкам водянистой раздробленной ткани. Кроме КТ-снимков, по каждому пациенту имелись сведения о характере полученной травмы, направлении основного удара, наличии повреждений тканей мозга того или иного типа.
По результатам сравнения расчетных данных с клиническими данными по 18 пациентам, получившим черепно-мозговые травмы различной тяжести, а также с имеющимися в медицинской литературе качественными описаниями биомеханики черепно-мозговой травмы, сделан вывод об удовлетворительной описательной способности модели. На рис. 3.7 приведен пример качественного совпадения рассчитанных и фактических зон поражения мозга.
При уменьшении площади соударения (аналог удара острым предметом) наблюдается область концентрации сдвиговых напряжений, совпадающая с очагом гематомы в 3 случаях. Зависимость между наличием повреждений в области противоудара и концентрацией отрицательных давлений не наблюдается.
