CC BY
19
DOI: 10.15593/RZhBiomeh/2019.1.10 УДК 531/534:[57+61], 539.3
ПОСТРОЕНИЕ МЕХАНИЧЕСКОЙ МОДЕЛИ ЭЛЕМЕНТОВ ГЕТЕРОГЕННОЙ СРЕДЫ НА ОСНОВЕ ЧИСЛЕННО-ЦИФРОВОГО АЛГОРИТМА ОБРАБОТКИ ДАННЫХ КОМПЬЮТЕРНОЙ ТОМОГРАФИИ
О.В. Герасимов1, Д.В. Бережной1, П.В. Большаков1, Е.О. Стаценко2, О.А. Саченков1,3
1 Кафедра теоретической механики Института математики и механики им. Н.И. Лобачевского Казанского (Приволжского) федерального университета, Россия, 420008, Казань, ул. Кремлевская, 18, e-mail: [email protected]
2 Научно-исследовательская лаборатория «Современные геоинформационные и геофизические технологии» Института геологии и нефтегазовых технологий Казанского (Приволжского) федерального университета, Россия, 420008, Казань, ул. Кремлевская, 4/5
Кафедра основ конструирования Казанского национального исследовательского технического университета им. А.Н. Туполева, Россия, 420008, Казань, ул. Карла Маркса, 10
Аннотация. В работе рассматривается один из возможных численных подходов к описанию процессов деформирования гетерогенных сред под действием внешних нагрузок, основанный на моделировании структуры расчетной области с учетом данных ее компьютерной томографии. Предлагаемый подход позволяет моделировать поведение пористой среды с учетом ее структурных свойств на основе методов неразрушающего контроля. На первом этапе исследования методом компьютерной томографии проводится сканирование расчетного образца, далее полученные данные оцифровываются, после чего цифровой прототип структуры образца с соответствующими весами вносится в конечно-элементный алгоритм расчета на этапе формирования локальных матриц жесткости расчетной области, интегрирование которых проводится методом центральных прямоугольников. В работе проведена оценка влияния данных томографического исследования на сходимость предложенной численной методики, решены тестовые задачи. Приведены результаты решения модельной задачи для участка диафиза бедренной кости. Полученные численные результаты отображают влияние точности аппроксимации геометрии образца, а также иллюстрируют зависимость поля перемещений от структуры материала.
Ключевые слова: математическое моделирование, негомогенные среды, компьютерная томография, метод конечных элементов, цифровой прототип.
Введение
В настоящее время наиболее перспективным направлением в моделировании поведения гетерогенных сред является использование данных компьютерной томографии [15, 16]. Этот подход позволяет получать информацию о структуре сложных анизотропных, мелкозернистых или пористых материалов, а также иных тел, характеризующихся неоднородностью. Эта задача особенно актуальна в ортопедической клинической практике. Так, дефекты костной ткани, связанные с всевозможными патологиями, могут оказать большое влияние на качество проводимого лечения, потому на этапе диагностики необходимо получить максимально возможную информацию о прочности и жесткости костной ткани [6, 12, 18, 20, 21].
© Герасимов О.В., Бережной Д.В., Большаков П.В., Стаценко Е.О., Саченков О.А., 2019 Герасимов Олег Владимирович, лаборант-исследователь кафедры теоретической механики, Казань Бережной Дмитрий Валерьевич, к.ф.-м.н., доцент, доцент кафедры теоретической механики, Казань Большаков Павел Владиславович, лаборант-исследователь кафедры теоретической механики, Казань Стаценко Евгений Олегович, младший научный сотрудник научно-исследовательской лаборатории, Казань
Существует несколько подходов к построению моделей элементов пористых сред, деформирующихся под действием внешних нагрузок. В первую очередь к ним можно отнести аппроксимацию распределения неоднородности методом средних пересечений линий [1, 2, 10]: в этом случае формулируются соотношения, связывающие компоненты тензора упругих констант и тензора структуры [1, 2, 7, 9, 10], характеризующего осреднённое направление пор. Другой подход - сведение анизотропии материала к ортотропии путём определения констант из численных экспериментов [5, 8, 13, 14, 17, 22, 23, 24]. В данной работе предлагается подход, основанный на учёте особенностей структуры пористого материала, выявленных по данным компьютерной томографии при построении численной модели.
Проведение компьютерной томографии элемента пористой среды предполагает создание его цифрового прототипа, представляющего собой трехмерный целочисленный массив определенной структуры [15, 16]. Элементами этого массива являются единицы или нули, которые характеризуют наличие или отсутствие вещества в определенном микроэлементе объема и отражают рентгеновскую плотность согласно шкале Хаунсфилда. Таким образом, цифровой прототип представляет структуру элемента пористой среды в виде совокупности элементарных микрообъемов (кубиков), про каждый из которых известно - содержит он вещество (кости) или нет. По этим данным на основе какого-либо приближенного метода можно построить дискретную механическую модель элемента пористой среды. Наиболее удобным с алгоритмической точки зрения в этом случае методом дискретизации (а в дальнейшем и расчета) будет метод конечных элементов.
Наивысшая точность расчета будет достигнута в случае моделирования каждого микрообъема трехмерным конечным элементом сплошной среды [15, 25, 26, 27]. Но в этом случае затраты ресурсов ЭВМ на создание дискретной модели, постпроцессорную обработку результатов и на этапе процессорных вычислений превысят все разумные пределы. Поэтому представляется целесообразным увеличить размеры конечных элементов, а каждый определенный в цифровом прототипе микрообъем считать окрестностью точки интегрирования при численном вычислении матрицы жесткости конечного элемента [19, 28]. Остается неясным вопрос определения числа квадратур в каждом из направлений при интегрировании локальной матрицы жесткости. Если число точек интегрирования по каждой координате мало, то точность решения задач на основе такого элемента может быть невысокой, так как в самом простом случае при интегрировании придется применять метод прямоугольников. Однако использование конечных элементов с очень большим числом квадратур хоть и повышает точность интегрирования в пределах элемента, но может увеличить жесткость расчетной области из-за возможно малого количества элементов.
Поэтому целью работы представляется реализация методики статического расчета элементов трехмерных пористых объектов на основе трехмерного изопараметрического линейного конечного элемента сплошной среды, построенного на основе выявленного по данным компьютерной томографии исследуемой области ее цифрового прототипа.
Материалы и методы
Создание цифрового прототипа
Для некоторого элемента расчетной области с помощью методов компьютерной томографии определяется соответствующий ему набор данных, представляющих собой трехмерную структуру. Интенсивность окрашивания микрообъемов структуры, слагающих расчетную область, отражает их рентгеновскую плотность согласно
шкале Хаунсфилда. Процесс создания цифрового прототипа исследуемого объекта предполагает его разделение на большое число виртуальных микрокубиков (вокселей) размером Ах х Ау х & с координатами центра вокселя хк, ук, , где линейные размеры Ах, Ау и Д1 обычно бывают равны друг другу и определяются разрешающей способностью компьютерного томографа.
Значения данных компьютерной томографии, соответствующие вокселю, бинаризируются по заданному порогу, определяемому, например, методом Отсу, отделяя плотную костную структуру от вещества в порах: если значение компьютерной томографии в текущем вокселе больше порогового значения, то этот воксель характеризует костную ткань, в противном случае - пору (рис. 1).
Рис. 1. Значения данных компьютерной томографии после бинаризации для одного
слоя в локальной области
Численное интегрирование матрицы жесткости с учетом данных компьютерной томографии
Рассмотрим методику построения конечного элемента пористой сплошной среды на основе широко известной методики [3, 4, 11, 29] построения восьмиузлового трехмерного изопараметрического конечного элемента сплошной среды с линейной аппроксимацией геометрии и поля перемещений в локальных координатах. В рамках изопараметрической концепции для аппроксимации геометрии (радиуса-вектора точки {г}) и исходных перемещений (вектора перемещений точки 9) используется одинаковая система функций:
ix] 8 'xn '
Г = = <j У = Z yn
z n=1 zn ,
u 8 un
{e} = - v = Z vn
w n=1 Wn .
N (UZ),
N (UZ),
(1)
(2)
где Ып ( ) = _ (1 + )(1 + П«л)(1 + С «С ) - известные линейные функции формы, 8
Ъ,п, , С« - локальные координаты узлов элемента, и, V, w - проекции вектора перемещений на орты глобальной декартовой системы координат х, у, z . Соотношения (2) могут быть переписаны в матричной форме
{9} = [N]{9е}, (3)
где [ N ] - матрица аппроксимирующих функций, {9е} - вектор узловых перемещений.
Деформацию среды описывают компоненты линейных деформаций 8^ ,8^ ,8^ и деформации сдвига уху ,ууг , у. , которые представляются в виде приведенного вектора {8} и выражаются через перемещения (2) известными соотношениями Коши
du
dv
dw
du dv
dv dw
du dw
^xx ^ , ^ yy ^ , ^ zz ^ , Y xy ^ ^ ^ , Y yz ^ ^ ^ , Y zx ^ ^ ^ :
dx dy dz dy dx dz dy dz dx
(4)
которые, в свою очередь, также могут быть записаны в матричной форме
{8}=[ ¿]{е). (5)
Напряженное состояние представляет тензор напряжений Коши в виде компонент нормальных а.., а, а^ и касательных тху , т^ , т2Х напряжений, который,
как и тензор малых деформаций, можно записать в форме приведенного вектора {а}.
Закон Гука, связывающий приведенные векторы напряжений и деформаций, представляется в форме
{а} = [^(г)]{в}, (6)
где [Б (г )] - матрица упругих постоянных кусочно-однородного тела, которую можно представить в виде [5]
X + 2ц X X X + 2ц X X 0 0 00 00
[ D ( r )] = [ D\ ш(г>
X 0 0 0
X 0 0 0
X + 2Ц 0 0 0
0 Ц 0 0
0 0 Ц 0
0 0 0 ц
• ю(Г).
(7)
Здесь [Б] - тензор упругости для изотропного материала, X и ц - постоянные
Ламе, определяемые через модуль упругости Юнга Е и коэффициент Пуассона и в виде
Ц =
E
2 (1 + и)'
X =
2ци 1 - 2и
(8)
ю(г) - некая весовая функция, определяемая точкой в пространстве по данным компьютерной томографии.
Для вычисления матрицы жесткости конечного элемента используются известные соотношения [4, 29] вида
[Г шB(r)Т [D(r)][B(r)]dV ,
(9)
где [В (г)] - матрица, связывающая приведенный вектор малых деформаций и вектор узловых перемещений [3, 4, 29]:
{в } = [В (г )]{е*}. (10)
При переходе к локальным координатам соотношение (9) можно переписать в
виде
[ к ] = { } } [B (л, С )! [ Г> (л, С)][B (л, С )]х (11)
-1 -1 -1 (11)
х| J (п, С)| ®(5, п, м^.
где J (5,П,С)| - детерминант матрицы Якоби преобразования координат.
В отличие от построения локальной матрицы жесткости изотропной сплошной среды, когда для вычисления интеграла (11) используется метод Гаусса, для пористых сред приходится использовать метод центральных прямоугольников. Конечный элемент пористой сплошной среды, как и в случае изотропии, представляет собой выпуклый шестигранник с линейчатыми четырехузловыми боковыми поверхностями. Квадратурными точками в таком элементе являются координаты центров вокселей из цифрового прототипа модели, только их, как и размеры вокселей, необходимо пересчитать в локальных координатах для каждого элемента. Тогда
1 3 К г ~.т
и ;=1 к=1 (12)
[ D] [ B ( г ( ^, n,, Z * ))] | J ( ,Пj ,Z * )| A^AZ,
где Ъ>1 ,nj ,Zк - локальные координаты точек интегрирования; А^, An, AZ - величина шага в трёх направлениях в локальных координатах; CT ("t>i,nj ,Zk ) - веса квадратурной
формулы, определяемые значением данных компьютерной томографии в точке интегрирования; I, J, K - число квадратурных точек по каждой локальной координате в элементе.
В ходе работы была реализована программа на языке программирования С+ + для расчётов с использованием восьми узлового конечного элемента с тремя степенями свободы в каждом узле: перемещения вдоль соответствующих декартовых координатных осей. Решение задач осуществлялось на компьютере следующей комплектации: процессор - AMD Ryzen 7 1700 с 8 физическими и 16 логическими ядрами с частотой 3,7 ГГц, оперативная память - G.Skill Aegis 16GB DDR4 16GISB K2 3000 C16 с частотой 2933 МГц, материнская плата - MSI B350M MORTAR. Ввиду больших трудозатрат при интегрировании локальной матрицы жесткости в программе было реализовано распараллеливание алгоритма с использованием технологии OpenMP, в рамках которого расчет локальной матрицы жесткости для каждого конечного элемента производился на разных потоках. Полученные результаты визуализировались в программе ParaView.
Результаты и обсуждение
Для оценки влияния плотности заполнения конечного элемента вокселями были решены тестовые задачи, результаты которых сравнивались с результатами решения в пакете Ansys 14.5.
х
Тестовые задачи
Для оценки влияния плотности заполнения конечного элемента вокселями были решены тестовые задачи для одного конечного элемента в условиях различного нагружения. Рассмотрим полученные результаты для нагружения сжимающей силой: к верхней грани прикладывалась нагрузка P = -400Н, равномерно распределённая по четырём узлам; нижние узлы фиксировались в перемещениях по трём направлениям, модуль Юнга составлял 2 ГПа, коэффициент Пуассона - 0,3, размер куба -10 х 10 х10 мм.
Варьируя величину вокселей и устанавливая в каждом из них значение из данных компьютерной томографии, равное единице, произвели серию численных экспериментов для сплошного материала. В табл. 1 отображены количество вокселей и их соответствующие размеры, используемые при расчётах.
На рис. 2 приведены результаты, полученные численным интегрированием методом центральных прямоугольников для сплошного изотропного материала с использованием разного количества вокселей внутри конечного элемента.
Таким образом, при увеличении количества вокселей решение стремится к результатам, полученным в программном комплексе Ansys. Исходя из зависимости времени выполнения программы от количества вокселей внутри конечного элемента, можно сделать вывод о том, что для наиболее эффективного решения (оптимального по точности и времени выполнения) предпочтительно использовать порядка 50 вокселей на стороне куба. Результаты для других случаев нагружения конечного элемента аналогичны.
С целью оценки влияния густоты сетки при постоянной плотности компьютерной томографии были решены тестовые задачи для одноосного сжатия и изгиба. Задача ставилась для геометрии в виде параллелепипеда с одним характерным размером. В этом случае на одном торце фиксировались все перемещения в узлах, на противоположном - прикладывалась распределённая по узлам нагрузка (см. рис. 3).
Рассматривался образец со следующими параметрами данных компьютерной томографии: количество вокселей в направлении осей Ox, Оу, и 02 - 100х 500х 100 соответственно, размер вокселя - 0,2 х 0,2 х 0,2 мм. Свойства материала использовались из предыдущей тестовой задачи. Параметры конечно-элементной сетки приведены в табл. 2. Исследовались отдельно два вида нагружения: влияние сжимающей нагрузки вдоль оси Оу и действие изгибающей силы в направлении оси 02 (решение для силы, действующей вдоль оси Ох, аналогично в силу симметрии).
На рис. 4 отображены результаты, полученные для трёх видов конечно-элементной сетки. Следует отметить, что относительная погрешность решения несущественно зависит от количества вокселей внутри одного конечного элемента, что позволяет в допустимых пределах изменять размер и количество конечных элементов с целью лучшей аппроксимации геометрии.
Таблица 1
Количество вокселей на ребре куба и их соответствующие размеры
N 1 yvox 10 16 20 25 32 40 50 64 80 100 125 160 200 -
Svox, мм 2 1,25 1 0,8 0,625 0,5 0,4 0,3125 0,25 0,2 0,16 0,125 0,1 -
а
50
100
150
200
250
1000 2000 3000 4000 5000 6000 7000
б
Рис. 2. Результаты численных расчётов: а - относительная погрешность значений перемещений для узлов с приложенной нагрузкой; б - время выполнения программы, с; ось абсцисс - количество вокселей по одной координате
A
P,
, A
A-A ,л z
x
P?
b
Рис. 3. Расчётная схема: Р1, Р2 - прикладываемая сила; Ь, Ь и h - параметры геометрии в направлении осей Ох, Оу и 02 соответственно
0
h
Таблица 2
Параметры конечно-элементной сетки
NbK3 Количество вокселей в одном конечном элементе
1 5 1 1003
2 10 2 503
3 20 4 253
40 35 30 25 ^ 20 15 10 5 0
1003
503
1003 503 253
0
50
100
150
200
250
300
350
Рис. 4. Относительная погрешность значений перемещений для узлов с приложенной нагрузкой: красная линия - результаты для модели в условиях действия сжимающей нагрузки, синяя линия - под действием изгибающей силы; ось абсцисс - количество конечных элементов
Модельная задача
В модельной задаче использовались данные компьютерной томографии для диафиза бедренной кости крысы (рис. 5, а). Сканирование выполнялось с применением микро- и нанофокусной системы рентгеновского контроля для компьютерной томографии и 2.0-инспекции Phoenix V\tome\X S240 в лаборатории рентгеновской компьютерной томографии Института геологии и нефтегазовых технологий Казанского (Приволжского) федерального университета. Система оснащена двумя рентгеновскими трубками: микрофокусной с максимальным ускоряющем напряжением 240 кВ мощностью 320 Вт и нанофокусной с максимальным ускоряющем напряжением 180 кВ мощностью 15 Вт. Для первичной обработки данных и создания объёмной (воксельной) модели образца на базе рентгеновских снимков (проекций) использовалось программное обеспечение datos\x reconstruction. Зафиксированный в держателе образец помещался на вращающийся столик камеры рентгеновского компьютерного томографа на оптимальном расстоянии от источника рентгеновского излучения. Съёмка проводилась при ускоряющем напряжении 90-100 кВ и токе 140-150 мА. Размер исследуемой области - 4,21х 3,46 х 5,89 мм, количество вокселей в направлении соответствующих координатных осей Ox, Oy и Oz - 624 х 512 х 874, размер вокселя -6,747 х 6,747 х 6,747 мкм.
•/■Sa
а
Рис. 5. Исследуемый образец диафиза бедренной кости крысы: а - представление данных компьютерной томографии; б - сетка 1; в - сетка 3
В табл. 3. отображены параметры расчетных конечно-элементных сеток, а на рис. 5, б и в приведены конечные элементы для сетки 1 и 3 соответственно. Узлы нижнего торца фиксировались в перемещениях по трём направлениям, на верхнем торце прикладывалась равномерно распределённая по узлам продольная сжимающая сила, равнодействующая которой равна Р = —400 Н. Модуль упругости Юнга принимался равным 30 000 МПа, коэффициент Пуассона - 0,3.
На рис. 6, а и б приведены распределения осевых перемещений, полученных с использованием двух типов конечно-элементных сеток.
На сетке с меньшим количеством конечных элементов (первая сетка) значения перемещений усредняются одним конечным элементом по толщине, при дискретизации по толщине двумя конечными элементами (для второй и третьей сеток) порядок аппроксимации распределения перемещений на торце повышается. На рис. 6, в показаны доверительные интервалы средних перемещений узлов верхнего торца объекта для конечно-элементных сеток различной густоты, увеличение стандартного отклонения объясняется увеличением количества узлов на исследуемом торце и, как следствие, изменением формы его деформированного состояния. Время выполнения расчетов составило: сетка 1 - 3,07 ч, сетка 2 - 3,16ч, сетка 3 - 3,53ч. Видно, что затраченное время реализации алгоритма решения задачи на ЭВМ для различных сеток сопоставимо. Это объясняется тем, что формируемые матрицы жесткости конечного элемента для разных сеток вычисляются в одних и тех же квадратурных точках. Однако эта ресурсозатратная операция для конкретной расчетной области проводится только один раз, а при дальнейших расчетах сформированная ранее матрица жесткости просто считывается с носителя информации.
Таблица 3
Параметры разбиения на конечно-элементную сетку
Сетка 1 Сетка 2 Сетка 3
По толщине 1 2 2
По высоте 4 10 10
В окружном направлении 11 30 60
б
в
-0,14
0
44 600 1200
в
Рис. 6. Распределение осевых перемещений, мм, в рассматриваемом образце в условиях действия равномерно распределённой по узлам сжимающей нагрузки: а - сетка 1; б - сетка 3, в - среднее продольное перемещение верхнего торца со стандартным отклонением при различном разбиении (ось абсцисс - количество
конечных элементов)
Заключение
1. В работе представлен один из возможных численных подходов к описанию процессов деформирования гетерогенных сред под действием внешних нагрузок, основанный на моделировании структуры расчетной области с учетом данных ее компьютерной томографии.
2. Для оценки влияния данных компьютерной томографии на сходимость предложенной численной методики были решены тестовые задачи. Полученные результаты применялись для определения подходящих параметров (числа вокселей в конечном элементе) модельной задачи.
3. Полученные численные результаты отображают влияние точности аппроксимации геометрии образца, а также иллюстрируют зависимость поля перемещений и напряженно-деформированного состояния от структуры материала.
4. Данный подход позволяет исследовать поведение пористого материала под действием внешних нагрузок на основе оптической плотности материала. Исследования могут быть расширены с использованием другого вида обработки данных компьютерной томографии, путем усовершенствования формулы численного интегрирования локальной матрицы жесткости, а также с применением иных типов конечных элементов. Вопрос построения полей напряжений и деформаций остался вне рамок исследований данной работы, поскольку прямое построение напряженно-деформированного состояния в рамках предложенной модели метода конечных элементов проводится легко и однозначно, но возникает вопрос достоверности полученных результатов. Вполне возможно, что наилучшую точность можно будет получить путем осреднения напряженно-деформированного состояния, вычисленного в вокселях, но этот дискуссионный вопрос требует дополнительных экспериментальных исследований.
Благодарности
Публикация осуществлена при финансовой поддержке РФФИ и Правительства Республики Татарстан в рамках научных проектов № 18-41-160025 и № 18-41-160018.
Авторы выражают отдельную благодарность сети независимых диагностических центров «Пикассо».
Список литературы
1. Киченко А.А., Тверье В.М., Няшин Ю.И., Заборских А.А. Экспериментальное определение тензора структуры трабекулярной костной ткани // Российский журнал биомеханики. - 2011. - Т. 15, № 4. -С. 78-93.
2. Киченко А.А., Тверье В.М., Няшин Ю.И., Симановская Е.Ю., Еловикова А.Н. Становление и развитие классической теории описания структуры костной ткани // Российский журнал биомеханики. - 2008. - Т. 12, № 1. - С. 69-89.
3. Крылов О.В. Метод конечных элементов и его применение в инженерных расчетах: учеб. пособие для вузов. - М.: Радио и связь, 2002. - 104 с.
4. Сагдатуллин М.К., Бережной Д.В. Постановка задачи численного моделирования конечных деформаций // Прикладные математические науки. - 2014. - Т. 8, № 35. - С. 1731-1738. DOI: 12988/ams.2014.4283
5. Саченков О.А., Герасимов О.В., Королева Е.В., Мухин Д.А., Яикова В.В., Ахтямов И.Ф., Шакирова Ф.В., Коробейникова Д.А., Хань Х.Ч. Построение негомогенной конечно-элементной модели по данным компьютерной томографии // Российский журнал биомеханики. - 2018. - Т. 22, № 3. - С. 332-344.
6. Саченков О.А., Хасанов Р.Ф., Андреев П.С., Коноплев Ю.Г. Численное исследование НДС тазобедренного сустава при ротационной остеотомии проксимального участка бедренной кости // Российский журнал биомеханики. - 2016. - T. 20, №3. - С. 257-271.
7. Тверье В.М., Няшин Ю.И., Никитин В.Н., Оборин Л.Ф. Механическое давление как основа биомеханического моделирования зубочелюстной системы человека // Российский журнал биомеханики. - 2014. - Т. 18, № 1. - С. 24-35.
8. Харин Н.В., Воробьев О.В., Бережной Д.В., Саченков О.А. Методика построения репрезентативной модели по данным компьютерной томографии // Вестник Пермского национального исследовательского политехнического университета. Механика. - 2018. - № 3. - С. 95-102.
9. Чикова Т.Н., Киченко А.А., Тверье В.М., Няшин Ю.И. Биомеханическое моделирование трабекулярной костной ткани в состоянии равновесия // Российский журнал биомеханики. - 2018. -Т. 22, № 3. - С. 282-291.
10. Экспериментальные методы в биомеханике / под ред. Ю.И. Няшина, Р.М. Подгайца. - Пермь: Изд-во Перм. гос. техн. ун-та, 2008. - 400 с.
11. Abdrakhmanova A.I., Gariffulin I.R., Davydov R.L., Sultanov L.U., Fakhrutdinov L.R. Investigation of strain of solids for incompressible materials // Applied Mathematical Sciences. - 2015. - Vol. 9. -P. 5907-5914.
12. Baltina T., Ahmetov N., Sachenkov O., Fedyanin A., Lavrov I. The influence of hindlimb unloading on bone and muscle tissues in rat model // BioNanoScience. - 2017. - Vol. 76 № 1. - P. 67-69.
13. Carniel T.A., Klahr B., Fancello E.A. On multiscale boundary conditions in the computational homogenization of an RVE of tendon fascicles // Journal of the Mechanical Behavior of Biomedical Materials. - 2019. - Vol. 91. - P. 131-138.
14. Carniel T.A., Klahr B., Fancello E.A. A multiscale numerical approach for the finite strains analysis of materials reinforced with helical fibers // Mechanics of Materials. - 2018. - Vol. 126. - P. 75-85.
15. Eggermont F., Derikx L.C., Free J., van Leeuwen R., van der Linden Y.M., Verdonschot N., Tanck E. Effect of different CT scanners and settings on femoral failure loads calculated by finite element models // Journal of Orthopaedic Research. - 2018. - Vol. 36, iss. 8. - P. 2288-2295.
16. Eggermont F., Derikx L.C., Verdonschot N., Van Der Geest I.C.M., De Jong M.A.A., Snyers A., Van Der Linden Y.M., Tanck E. Can patient-specific finite element models better predict fractures in metastatic bone disease than experienced clinicians // Bone and Joint Research. - 2018. - Vol. 7, iss. 6. -P. 430-439.
17. Gabidullin M.G., Kayumov R.A., Rakhimov R.Z., Temlyakov A.V. Inter-relation between structures and heat-transfer properties of porous ceramic materials // Stroitelnye Materialy. - 2005. - Vol. 9. - P. 62-66.
18. Gerasimov O., Shigapova F., Konoplev Y., Sachenkov O. Evaluation of the stress-strain state of a one-dimensional heterogeneous porous structure // IOP Conference Series: Materials Science and Engineering. -
2016. - Vol. 158. - 012036. DOI: 10.1088/1757-899X/158/1/012036
19. Giovannelli L., Rodenas J.J., Navarro-Jimenez J.M., Tur M. Direct medical image-based Finite Element modelling for patient-specific simulation of future implants // Finite Elements in Analysis and Design. -
2017. - Vol. 136. - P. 37-57.
20. Greatrex F., Montefiori E., Grupp T., Kozak J., Mazza C. Reliability of an integrated ultrasound and stereophotogrammetric system for lower limb anatomical characterisation // Applied Bionics and Biomechanics. - 2017. - Vol. 2017. - 4370649. DOI: 10.1155/2017/4370649
21. Hettich G., Schierjott R.A., Ramm H., Graichen H., Jansson V.b, Rudert M., Traina F., Grupp T.M. Method for quantitative assessment of acetabular bone defects // Journal of Orthopaedic Research. - 2018. DOI: 10.1002/jor.24165
22. Kayumov R.A., Muhamedova I.Z., Tazyukov B.F., Shakirzjanov F.R. Parameter determination of hereditary models of deformation of composite materials based on identification method // Journal of Physics: Conference Series. - 2018. - Vol. 973, iss. 1. - 012006.
23. Kayumov R.A., Muhamedova I.Z., Tazyukov B.F., Shakirzjanov F.R. Parameter determination of hereditary models of deformation of composite materials based on identification method // Journal of Physics: Conference Series. - 2018. - Vol. 973, iss. 1. - 012006.
24. Kayumov R.A. Structure of nonlinear elastic relationships for the highly anisotropic layer of a nonthin shell // Mechanics of Composite Materials. - 1999. - Vol. 35, № 5. - P. 409-418.
25. Marcian P., Florian Z., Horackova L., Kaiser J., Borak L. Microstructural finite-element analysis of influence of bone density and histomorphometric parameters on mechanical behavior of mandibular cancellous bone structure // Solid State Phenomena. - 2017. - Vol. 258. - P. 362-365.
26. Marcian P., Wolff J., Horackova L., Kaiser J., Zikmund T., Borak L. Micro finite element analysis of dental implants under different loading conditions // Computers in Biology and Medicine. - 2018. - Vol. 96. -P. 157-165.
27. Ridwan-Pramana A., Marcian P., Borak L., Narra N., Forouzanfar T., Wolff J. Finite element analysis of 6 large PMMA skull reconstructions: A multi-criteria evaluation approach // PLoS ONE. - 2017. - Vol. 12. - e0179325. DOI: 10.1371/journal.pone.0179325
28. Semenova E., Gerasimov O., Koroleva E., Ahmetov N., Baltina T., Sachenkov O. Automatic processing and analysis of the quality healing of derma injury // Advances in Intelligent Systems and Computing. - 2019. -Vol. 831. - P. 107-113.
29. Sultanov L.U. Analysis of finite elasto-plastic strains: integration algorithm and numerical examples // Lobachevskii Journal of Mathematics. - 2018. - Vol. 39, iss. 9. - P. 1478-1483. DOI: 10.1134/S1995080218090056
MECHANICAL MODEL OF A HETEROGENEOUS CONTINUUM BASED ON NUMERICAL-DIGITAL ALGORITHM PROCESSING COMPUTER TOMOGRAPHY DATA
O.V. Gerasimov, D.V. Berezhnoi, P.V. Bolshakov, E.O. Statsenko, O.A. Sachenkov (Kazan, Russia)
The paper considers one of the possible numerical approach for describing mechanical behavior of heterogeneous media based on modeling the structure using CT data is presented. The proposed approach allows to model the mechanical behavior of a porous medium taking into account the structure properties. At the first stage of the study sample is scanned by computed tomography then the data is digitized and digital prototype is built. This prototype with corresponding weights is used in the forming the local stiffness matrices for finite element method. The integration of stiffness matrices carried out by the central rectangles method. The effect of tomographic data density on the convergence of the proposed numerical method was investigated. The results of solving a typical problem for the diaphysis of the femur are given. The obtained numerical results reflect the influence of the accuracy of approximation of the sample geometry, and also illustrate the dependence of the displacement field on the structure of the material.
Key words: mathematical modelling, inhomogeneous media, computed tomography, finite element method, digital prototype.
Получено 5 января 2019
