Нелинейное вязкоупругое поведение волокон миокарда: гистерезис кривых «Сила-длина» и релаксация деформации Текст научной статьи по специальности «Физика»

УДК 531/534:57+612.7

Российский

Журнал

Биомеханики

www. biomech. ас. ru

НЕЛИНЕЙНОЕ ВЯЗКОУПРУГОЕ ПОВЕДЕНИЕ ВОЛОКОН МИОКАРДА: ГИСТЕРЕЗИС КРИВЫХ «СИЛА-ДЛИНА» И РЕЛАКСАЦИЯ ДЕФОРМАЦИИ

А.В. Кобелев**, Р.М. Кобелева*, Ю.Л. Проценко*, И.В. Берман***

‘Екатеринбургский филиал института экологии и генетики микроорганизмов УрО РАН, ул. Первомайская, 91, 620219, Екатеринбург

‘‘Институт физики металлов УрО РАН, e-mail: kobelev@imp.uran.ru ‘‘‘Государственный университет Сан Хосе, Калифорния, США

Аннотация. Исследовано динамическое поведение образцов сердечной мышцы в пассивном состоянии при периодической и ступенчатой нагрузке. При теоретическом описании использовалась часть из идентифицированных реологических двумерных моделей, предложенных нами ранее для анализа статических кривых «сила-деформация» миокарда. Вязкие свойства миокарда моделировали введением демпферов, включенных параллельно некоторым первичным элементам модели. Сравнение теоретических кривых с экспериментальными кривыми гистерезиса «сила-деформация», полученными нами на образцах миокарда морских свинок и кошек, позволило получить оценки не только упругих модулей составляющих элементов, как в статическом случае, но и соответствующих коэффициентов вязкости. Установлено, что по форме кривой гистерезиса можно сделать заключение о том, вязкость каких элементов преобладает в сердечном волокне. Впервые получены теоретические зависимости скорости релаксации первичных элементов, обладающих различной вязкостью, при заданной нагрузке в виде периодической или ступенчатой функции времени.

Ключевые слова: биомеханика сердечной мышцы, нелинейные кривые

«деформация-напряжение», вязкость, гистерезис, релаксация деформации.

Введение

До настоящего времени при анализе зависимостей «пассивное напряжение-деформация» для мышцы, как правило, используются одномерные модели реологического состояния, другими словами, не учитываются поперечные деформации [1, 2]. В то же время, экспериментально установлена нелинейная зависимость поперечного размера деформируемого препарата от продольного растяжения [3, 4]. Впервые учет существенного вклада поперечных деформаций в происхождение нелинейной зависимости «сила-длина» для изолированных препаратов миокарда был проведен в работах [4-6], где исследованы двумерные топологические структуры, при деформации которых определяющую роль играет геометрический эффект «схлапывания». Возможность адекватного описания статической зависимости «пассивное напряжение-деформация» для папиллярных мышц была исследована с использованием классических форм в виде линейных упругих элементов,

© А.В. Кобелев, Р.М. Кобелева, Ю.Л. Проценко, И.В. Берман, 2003

организованных в структуру, подобную морфологически-функциональной единице миокардиальной ткани [4].

В рамках предложенной модели показано, что имеет достаточное основание предположение о функционировании морфологических структур миокардиальной ткани с различной линейной жесткостью в разных участках физиологического диапазона деформаций без введения дополнительных предположений об общем характере нелинейности во всем диапазоне деформаций. Так, в диапазоне начальных длин препарата при деформации сначала включаются более мягкие элементы, оказывающие наименьшее сопротивление генерации мышцей активного напряжения. С ростом деформации подключаются более жесткие элементы, препятствующие разрушению ткани. Отметим, что введение геометрического фактора ранее было успешно использовано при математическом моделировании зависимости активного напряжения от длины мышцы [3].

Одной из физиологических функций вязкости в живых тканях является предохранение их от разрывов при резких деформациях. Вязкостью обладают практически все биологические ткани, однако ее природа в мышечной ткани не ясна до сих пор. Установлено, что ответ на достаточно быструю деформацию сопровождается релаксацией мышечного напряжения в изолированном невозбужденном состоянии. В ряде экспериментов [1, 2, 8] показано, что процессы релаксации напряжения зависят не только от величины скачка деформации, но и от исходной (начальной) длины препарата. Такая зависимость не воспроизводится ни в одной из реологических моделей без введения дополнительного искусственного предположения о зависимости коэффициента вязкости от длины препарата [9].

Вязкие свойства проявляются на всех уровнях организации мышечной ткани -от клеточного до органного. Причинами вязкого поведения могут являться, по-видимому, фильтрация тканевой жидкости через поры, трение скольжения сократительных элементов и других внутриклеточных структур при их движении в цитоплазме при деформации и, наконец, собственное вязкоупругое поведение молекул белков цитоскелета, в частности, титина [9].

При динамических нагружениях геометрическая структура сказывается не только на упругих, но и на вязких свойствах миокарда. Следовательно, при нагрузке, зависящей от времени и от начальных условий, наряду с упругой нелинейностью, могут проявляться и нелинейные вязкие свойства миокарда.

Цель исследования

Исследование ставит целью описать упругие и вязкие свойства пассивного миокарда при периодической и ступенчатой нагрузке в области физиологических деформаций (а ~ 0,3) на основе предложенных в [4] простых двумерных структур, отражающих функциональные особенности морфологической единицы миокарда, составленных из линейных гуковских элементов с включенными параллельно им вязкими элементами (ньютоновскими демпферами). На основе сравнения с экспериментально полученными кривыми гистерезиса «сила-длина» оценить адекватность той или иной модели без введения зависимости коэффициентов упругости и вязкости от длины и определить значения вязких постоянных составляющих элементов. Исследовать особенности релаксации деформации при заданном периодическом или ступенчатом изменении нагрузки в зависимости от величины параметров вязкости первичных элементов.

Методика физиологического эксперимента Приготовление препаратов

Эксперименты выполнены на папиллярных мышцах правого желудочка крыс, кроликов, кошек и морских свинок. Перед вскрытием животного проводилась инъекция гепарина для уменьшения свертываемости крови. Доза гепарина для кроликов составляла в среднем 0,5 мл/кг, для крыс - 0,2 мл/кг. Спустя 5-15 минут после инъекции гепарина проводилась эвтаназия в соответствии с международными правилами [10].

Использовали раствор Кребса-Рингера, содержащий (в миллилитрах): 98 №С1,

4,7 КС1, 2,4 М§Б04*7Н20, 25 КаНС03, 1,25 СаС12*2Н20, глюкоза 4,5, при температуре 300С, барбатируемый 100%-м кислородом. После приготовления раствора измеряли величину pH при температуре раствора 300С. Буферные компоненты раствора (ТЯНМЛ, 81&та-ЛШпсИ) обеспечивали уровень pH в интервале от 7,30 до 7,35.

После усыпления животного производили разрез по средней грудинной линии, затем сердце быстро извлекали и помещали в препаровальную термостатируемую кювету объемом 50 мл с физиологическим раствором. Препарирование вели под бинокулярным микроскопом МБС-10. Правый желудочек сердца вскрывали по нисходящей коронарной артерии вдоль межжелудочковой перегородки, иссекали папиллярную мышцу и помещали в ванночку экспериментальной установки. При этом отбирали преимущественно тонкие и длинные препараты. Концы препарата фиксировали на штоках датчика силы и мотора длины наложением ниточных лигатур на оба конца препарата (рис. 1).

На начальном этапе подготовки к эксперименту осуществляли отбор и оценку сократительных свойств препаратов. Если сразу после закрепления препарата развиваемая им активная сила не превышала 50 мг, то препарат считался поврежденным и дальнейшие эксперименты над ним не проводились. Если же развиваемая препаратом сила составляла около 100 мг для крыс, а для препаратов кролика около 400 мг, то после наложения дополнительных лигатур препарат изометрически сокращался в течение 30-60 минут (период «врабатывания» препарата). В течение этого периода времени при непрерывном питании физиологическим раствором при температуре 300С, перфузии кислородом и включенной стимуляции активная сила, развиваемая препаратом, возрастала в 2-3 раза, а пассивное напряжение

падало. По истечении периода врабатывания препарата, при условии восстановления его сократительных свойств, реализовывали экспериментальный протокол.

Протокол эксперимента

В качестве тестирующих механических воздействий использовали периодические деформации синусоидальной, прямоугольной, треугольной, “пилообразной” формы с изменяющейся частотой. С помощью такого типа воздействий нами измерены динамические реологические характеристики папиллярных мышц (компоненты динамического модуля упругости и диссипативные потери энергии при деформации). Деформация задавалась в пределах 10% от рабочей длины с целью выявить наличие нелинейности. Пилообразное изменение длины папиллярной мышцы задавали с амплитудой от 0,05 до 0,3 1К со скоростью от 0,001 до 0,1 длины мышцы в секунду.

Препарат разделяли на ряд сегментов с помощью шипов кактуса. Изменения длины сегментов регистрировали с помощью компьютерной системы технического зрения [11]. Для выявления продольного перемещения метки по оси препарата использовали среднее значение положений метки с одного и другого края препарата. Этот прием позволил оценить характер чисто продольного перемещения. Несмотря на изометрический режим, справедливый для препарата в целом, наблюдали изменения длины сегмента в продольном направлении, которые были наиболее выражены в фазе расслабления и достигали нескольких процентов от длины сегмента.

Образцы препаратов

ЛВ.013: морская свинка, 10 = 3,88 мм, = 4,82 мм, диаметр к = 0,42 мм, четыре

метки, амплитуда продольного растяжения в форме пилы = 940 мкм, период Т = 1000 с и

Л2038: кот, 10 = 3,65 мм, 1К = 4,54 мм, диаметр к = 0,9 мм, две метки, амплитуда продольного растяжения в форме пилы =890 мкм, период Т = 333 с.

Результаты эксперимента

При линейном растяжении препарата пассивное напряжение возрастало приблизительно по экспоненциальному закону. Характер изменения длины центральных сегментов при этом был самый разнообразный - от линейного удлинения, через немонотонный, до укорочения сегмента при растяжении всего препарата. Последнее наблюдение при нанесении меток тушью зарегистрировано и в работе [12].

При “пилообразной” деформации препарата в координатах “напряжение-деформация” нами зарегистрирована характерная петля гистерезиса [1]. При этом показано, что движение по петле гистерезиса в исследованной вязкоупругой структуре пассивного миокарда всегда направлено по часовой стрелке, что характерно для диссипативных систем. Вид петли гистерезиса свидетельствует о нелинейности реологических характеристик препарата миокарда. В данном случае под “линейностью” понимается пропорциональность амплитуд напряжений и деформаций. Нелинейность петли гистерезиса проявляется при деформациях, превышающих 1-2 % от текущей длины препарата.

Обнаружено, что нелинейность связи “деформация-напряжение” при динамических испытаниях проявляется тем больше, чем больше рабочая длина папиллярной мышцы, т.е. длина, на которую растянули препарат перед испытаниями.

Такое поведение характерно как для всего препарата, так и для его центрального сегмента.

Значения комплексного модуля упругости и энергии диссипации препарата и его центрального сегмента зависят от скорости деформации. При увеличении скорости деформации препарата на порядок единиц длин мышцы в секунду величина модуля возрастает на десять процентов. Характер зависимости в целом препарате не отличается от вида зависимости для сегмента. Однако принципиально важным является то обстоятельство, что абсолютные значения комплексного модуля упругости при всех исследованных скоростях деформации для целого препарата выше, чем для сегмента. В группе “мягких” препаратов (Рра5 / Рс < 0,5) тангенциальный модуль Юнга по

начальной части гистерезиса в целом препарате составляет 15+2 кН/м2, а в группе жестких препаратов - на порядок больше. Наличие более податливого центрального участка папиллярной мышцы в эксперименте проявляется в том, что при любой деформации растяжения на концах препарата его центральный сегмент деформируется сильнее. Модуль упругости составляет для центральных сегментов “мягких” препаратов 11+2 кН/м . Таким образом, модуль упругости как целого препарата, так и его сегмента возрастает с увеличением скорости деформации и с увеличением исходной длины.

Величины значений тангенциального модуля упругости, полученные в наших экспериментах, совпадают со значениями модуля пассивного напряжения миокарда, полученными другими авторами: 14 кН/м2 - для папиллярной мышцы кошки [13],

15,7 кН/м2 - для папиллярной мышцы крысы [14]. Различная величина деформации и значений модуля упругости целого препарата и его сегмента свидетельствует о выраженной продольной неоднородности реологических характеристик миокарда.

Кроме того, показано, что при задании ступенчатой деформации вязкоупругое поведение сердечной мышцы проявляется в релаксации напряжения на концах препарата, которое не сопровождается изменением длины центрального сегмента. В сопоставлении с данными об изменении поперечных размеров, полученными в других экспериментах, можно предположить, что релаксация сопровождается изменениями размеров в различных сечениях образца.

Выбор и описание моделей

При моделировании вязкоупругих динамических свойств сердечной мышцы (без учета ее массы, что справедливо при достаточно низкой частоте периодической нагрузки) мы будем опираться на предложенные в работе [4] двумерные структуры, составленные из гуковских линейных упругих элементов, соединенных шарнирами без трения в точках типа а и Ь (рис. 2а). Вязкость учитывается введением ньютоновских демпферов, включенных параллельно упругим элементам. При этом часть из моделей, исследованных в [2], мы опустим, учитывая их неадекватные статические характеристики. Таким образом, в данной работе мы остановимся на моделях ромбической и параллельной геометрии 2-6 и 9 в обозначениях работы [4].

Рассмотрены два варианта постановки задачи исследования динамики. В первом из них считается заданной зависимость от времени деформации какого-либо первичного элемента модели в виде периодической функции синусоидальной или пилообразной формы. При этом предполагается, что этот элемент обладает вязкостью. Например, при деформации поперечного элемента а3 (V) вязкая сила пропорциональна

г]ъ йг3>1 Л, где г]ъ - коэффициент вязкости. Зависимости от времени деформации всех

Рис. 2. Схемы моделей RLH1 (a), RLH2 (b), RLH3 (c) ромбической (слева) и параллельной PH1 (d), PH2 (e), PH3 (f) (справа) геометрии. Первичные элементы, обладающие вязкостью, заштрихованы серым цветом. Кроме представленных, рассмотрены ромбические модели без продольного элемента (RH1, RH2, RH3), а также варианты всех моделей при h=0

(RL1, RL2, RL3, P1, P2, P3, R1, R2, R3)

других элементов, общая продольная деформация модели, а также результирующая сила Е(1) находятся в результате решения системы нелинейных уравнений с учетом вязкой силы, пропорциональной скорости деформации. В такой постановке удается получить решение задачи, задавая деформацию того элемента, который обладает вязкостью. Таким образом, рассмотрены варианты ромбической и параллельной геометрии, в которых последовательно учитывается вязкость наклонного цх с^ / Ж, продольного т]2 с^2/С и поперечного 7]3 С£3/С элементов (см. рис. 2).

Во втором варианте постановки задачи считается заданной зависимость от времени внешней нагрузки Е(/) в виде периодической или ступенчатой функции.

Скорости деформации первичных элементов с^/С, с^2/С и С£3/С, а также сами

деформации ех, £2 и е3 определяются при решении системы нелинейных

дифференциальных уравнений первого порядка с нулевыми начальными условиями. Полученные решения описывают релаксацию деформации различных элементов и структуры в целом в форме переходных процессов с выходом на квазистационарный предел.

В данной работе мы не рассматриваем вопросы, связанные с решением обратной задачи, т.е. с релаксацией напряжения в сердечной мышце при ступенчатой

деформации. Это требует использования моделей, в которых упругий и вязкий первичные элементы соединены последовательно (аналогичных модели Максвелла), что выходит за рамки настоящей работы и составит предмет отдельной публикации.

Уравнения для ромбической модели и модели параллелограмма

Приведем основные уравнения, описывающие реологические характеристики миокарда при учете вязких свойств в отдельных упругих элементах модели. Решения, описывающие кривые гистерезиса силы, получены при задаваемой деформации для моделей, наиболее адекватно описывающих статические упругие свойства [4]. Используем обозначения:

а - угол наклонного элемента;

И - поперечник муфты закрепления (заштриховано на рис. 2);

Е - сила нагрузки на муфту;

1, 10 - продольная длина, начальная продольная длина;

А, 12, 13 - длины наклонного, продольного и поперечного элементов;

Ао, 12о, 13о - соответствующие начальные длины;

вх, е2, £3 - однородные относительные деформации;

в1 = Св^/С, в2 = Св2/ С, в3 = Св3/С - скорости деформаций;

К, К2, К3 - коэффициенты жесткости элементов;

, Л2, Л3 - коэффициенты вязкости.

Вязкость наклонного элемента

Считая деформацию наклонного элемента е = в^) и ее производную в] = Св]/С заданными функциями времени, в модели ромба (рис. 2а) без продольного элемента КИ1 из условия равенства упругих и вязких сил в точках а и Ь получим:

2sina

f^L , Гі є Л

к k 1

уКі «2 У

є3

+ — = 0,

2cosa

к3

г Л 3 (1>

- + г є V1 ■

V К1 к 2 У

где кі = ^ , к2 = К2/20 5 к3 = К313о , / = ^/к2 5 К1 = к2 /к1 5 ^3 = к2 /к3 •

Используя геометрическое соотношение

БІпа = ——— 5 (2>

Уі

где У1 = Л (і + є1) 5 У3 = Лд (і + Є3 ) 5 Л = /10//0 5 Л3 = /30//0 5 — = Й//0 5 из первого уравнения (1) для деформации поперечного элемента получим:

є3 = -

1 - у

1 + _ у/2к3^3

(3)

A3 ^ ^ +^1/^ у

Выражение (3) позволяет получить полную деформацию:

s = -1 + VУ2 - - г)2 (4)

и в соответствии со вторым уравнением (1) и соотношением cos a = (1 + s')/y1 решение задачи определения силы f (t) при заданной деформации.

При решении этой задачи для модели ромба с продольными элементами КЬН1 с учетом £ = £2 справедливо второе уравнение (1), в левую часть которого добавляется слагаемое £2. Отметим, что для модели параллелограмма выполняется формула (4) с добавлением слагаемого у2 в правую часть и при учете соотношения /(V) = 2£2 (V) .

Вязкость продольного элемента

Для модели с вязким продольным элементом КЬН2 соответствующая система уравнений имеет вид:

Учитывая (2), из первого уравнения (5) найдем соотношение, обратное (3) и выражающее деформацию наклонного элемента через деформацию поперечного:

Зависимость от времени продольной деформации £2 и ее производной £2 считается заданной. При этом зависимость £3(£2) находится из решения уравнения

с учетом (6) и далее из (6) находим зависимость £(£2). Зависимость силы от времени / (V) находим далее из второго уравнения (5).

Для модели параллелограмма РН2 в системе уравнений (5) исключается слагаемое £2, и дополнительно появляется уравнение

Из последнего уравнения, учитывая (6), находим зависимости £3(£2) и £(£2). Для силы из (5) и (8) имеем

Вязкость поперечного элемента

Считаем заданной функцией времени деформацию поперечного элемента є3 и соответствующую скорость деформации Є3. Система уравнений для модели ромба КЬН3 в этом случае

*1

позволяет получить соотношение, аналогичное (6), для деформации наклонного элемента

_ . Єї Є3

2віпа — + — = 0,

*1 К3

(5)

Є, 772 .

2 СОБ^ —1 + є2 +— є2 = /.

К1 ^2

V КЛЄ3 у

(6)

(7)

(8)

(9)

и аналогично (4) для полной деформации

(10)

*1 *з к 2

_ 81 г

2соб а---------Ъ8 2 = у

(11)

Є1 = -

1 +

Уз - У

КА(Єз/Кз + 7з/к2 Є3)у

(12)

При этом справедлива формула (4) для продольной деформации £ = £2 , которая совместно со вторым уравнением (11) позволяет получить решение задачи определения

I (V ).

Решение задачи для параллелограмма РН3 дает аналогичную (11) систему уравнений без слагаемого £2 и с дополнительным условием / = 2£2 .

Уравнения релаксации деформации при заданной нагрузке

Другой подход к исследованию вязкоупругой динамики в рамках предложенных двумерных моделей состоит в решении системы дифференциальных уравнений для деформаций элементов при заданной зависимости от времени нагрузки / (V). Приведем систему уравнений для модели ромба без продольного элемента Ю3 при И=0.

Л

2б1п«

1 71 -

— + — Є

\

V К1 ( „

к

2соб«

2 У

1

К

к

(13)

V К1

к

= У (г).

2

Преобразуем систему (13) к виду, удобному для применения стандартных методов решения систем нелинейных дифференциальных уравнений:

771Є1 = -—Є1 +1 Т1(Є1,Єз)ї(t),

к

К

2

7 Єз = - — Єз +1 Тз(Є1,ЄзХ/ (г X

(14)

к

К

2

1 -

{ Л2

Уз - У

У1

и Тз = -2

1-

У1

V Уз- у у

Для модели ромба с продольным элементом КЬ123 при И=0 получим систему трех нелинейных дифференциальных уравнений

2

^ 1 71 - ^

—Є1 + -1 Є1 V К1 к2 У

ТГ1 + Є2 + 7- Є2 = у (г)

ко

1 73 ^

—Є3 +-Є3 + 2

К к

(

71

—Є1 +-^ Є1

V К1

к

2 У

Уз = 0, У1

(15)

Є2 (1 + Є2 ) = ЄАУ1 - Є3^зУз.

Последнее уравнение (15) получено дифференцированием соотношения (7) при

у = 0.

1

2

2

2

1

/, мм

/, мм

/, мм

Рис. 3. Кривые гистерезиса, полученные в модели ромба с вязким продольным элементом ЯЪН2 в области линейной статической кривой сила-длина при синусоидальной (а) и кубической пилообразной (Ь и с) нагрузке. Начальная продольная деформация 0,3, максимальная 0,6. На вставках приведены зависимости от времени для задаваемой и полной деформации. Значения параметров: 10=2 мм, /о1=1,07 мм, /о3=1,25 мм, ^=1,0 г/мм, К2=К3=0,01 г/мм,

Н2=0,8 ст/м (а), Н2=0,5 ст/мм (Ь, с)

Результаты моделирования

1. Гистерезис

Влияние вида задаваемой деформации

Для формы задаваемой деформации мы использовали синусоидальную и пилообразную функцию времени (кусочно-линейного, квадратичного и кубического вида) с расположением зубцов вверх и вниз. На рис. 3а представлена кривая гистерезиса для модели ромба при синусоидальной зависимости от времени задаваемой деформации продольного элемента в области начального линейного участка стационарной кривой «сила-деформация». Кривая напряжения имеет симметричную форму в соответствии с тем, что производная деформации, определяющая вязкие свойства, является непрерывной и плавной функцией. В координатах «напряжение-деформация» в этом случае получена линейная симметричная гистерезисная зависимость эллипсовидной формы. Кривые гистерезиса в случаях с входным сигналом пилообразной функции времени квадратичного и кубического вида имеют разрывы в верхней и нижней частях, соответственно расположению зубцов (рис. 3 Ь и с). Такого же рода особенности кривых гистерезиса вблизи точек скачков производной деформации характерны для других исследованных моделей.

Рис. 4. Кривые гистерезиса, полученные в моделях ромба с вязким мягким элементом КЬЮ при синусоидальной (а) и кубической пилообразной ^, :£) нагрузках и в моделях ромба с вязким жестким элементом КЬШ, соответственно (Ь, с и е). Геометрические параметры и параметры жесткости элементов указаны в подписи к рис. 3, Н3=0.05 (а), Н3=0.02 ^ и 1) Н1=0.5 (Ь),

Н1=0.3 ст/мм (с и е)

Влияние вязкости различных первичных элементов

Кривые гистерезиса, полученные для моделей с поперечным или наклонным вязким элементом, во всем диапазоне деформаций, включающем нелинейную область, имеют асимметрию либо при малых, либо при больших деформациях. Показано, что этот эффект определяется тем, что мягкий поперечный элемент работает в основном при малых деформациях и кривая гистерезиса имеет расширение при малых деформациях (рис. 4а). В моделях с вязким наклонным элементом при достижении больших величин деформаций, когда начинает работать жесткий элемент, кривая гистерезиса имеет расширение при больших деформациях (рис. 4Ь). Аналогичный эффект наблюдается и для пилообразной зависимости от времени задаваемой деформации (рис. 4 с-!}).

Рис. 5. Сопоставление кривых гистерезиса, полученных в различных моделях, с данными эксперимента ЛЯ013. Значения геометрических параметров и параметров упругости одинаковы во всех вариантах: /01=2,276 мм, /о2=3,88 мм, /о3=2,8 мм, ^=0,42 мм, ^=0,9 г/мм, К2=К3=0,02 г/мм. Синусоидальная нагрузка: (а) модель с вязким поперечным элементом ЯЬН3, Н3=1,5 ст/мм, максимальная поперечная деформация £3тах =-0,64, максимальная полная деформация £тах = 0,249; (Ь) модель с вязким наклонным элементом КЬШ при Н^16 ст/мм, £1тах = 0,07, £тах = 0,246. Пилообразная нагрузка: (с) Н3=2 ст/мм, £3тах =-0,61, £тах = 0,246,

(ф Н1=20 ст/мм, £1тах = 0,07 , £тах = 0,246

Сравнение с экспериментом

В натуральных экспериментах задавалось смещение закрепляющих концов препарата в виде пилообразной функции времени в связи с тем, что при линейном нарастании и спаде нагрузки упрощается анализ ответной реакции нелинейного объекта. Введение меток позволило исследовать продольную неоднородность реологических характеристик. Однако оказалось практически невозможно выделить структуры, локально ответственные за те или иные реологические характеристики препарата в целом. Однородность деформации по всей длине препарата не была подтверждена.

В численных экспериментах задавалась пилообразная и синусоидальная форма однородной деформации какого-либо из составляющих элементов модели. При этом форма полной результирующей продольной деформации была существенно различной и могла в той или иной мере моделировать реально задаваемую линейную форму

нагрузки. Кроме того, в численных экспериментах в предположении, что продольные линейные смещения концов образца могут трансформироваться в неоднородные продольные локальные деформации, рассмотрены нелинейные синусоидальные, квадратичные и кубические по форме задаваемые деформации. Неоднородность продольной деформации может быть исследована в моделях параллелограмма РЬН (рис. 2 ё, е, !), в которых деформации продольного элемента отличаются от деформации структур вблизи муфты закрепления.

Для сравнения экспериментальных данных с результатами моделирования использованы два эксперимента ЛВ.013 и Л^038, проведенные на папиллярных мышцах с различными геометрическими размерами и начальной жесткостью. В качестве первого шага мы ограничились ромбическими моделями (рис. 2 а, с).

Кривая гистерезиса, полученная в ромбической модели с продольным элементом и муфтой КЬН3, при учете вязкости мягкого поперечного элемента для синусоидальной задаваемой деформации, достаточно хорошо в целом описывает эксперимент ЛВ.013 (рис. 5а), хотя имеются отклонения как при малых, так и при максимальных деформациях. Теоретическая кривая в модели с вязким жестким наклонным элементом КЬН1 (рис. 5Ь) хорошо описывает максимальные деформации, но не работает при малых деформациях. На рис. 5 с и ё представлены соответствующие экспериментальные и модельные кривые гистерезиса для задаваемой деформации кусочно-линейной формы.

Следует обратить внимание, что значения параметров жесткости и вязкости, полученные в результате подгонки, несколько изменяются при изменении формы зависимости от времени задаваемой деформации элемента. Для более детального сравнения нужны численные эксперименты при разной крутизне (скорости изменения) пилообразной нагрузки и при различных начальных длинах препарата при одной и той же скорости изменения деформации в условиях нелинейности стационарной кривой «сила-длина».

Для эксперимента Л2038 характерна нерегулярная форма опытной кривой гистерезиса. В этом случае результаты моделирования при задаваемой деформации синусоидального вида (рис. 6 а, Ь) по средним квадратичным отклонениям хуже, чем при использовании пилообразной формы кривой (рис. 6 с, ё). Однако, и в этом случае согласие экспериментальных данных с теоретическими кривыми наблюдается либо при малых, либо при максимальных деформациях.

Для оценки величины упругих и вязких параметров, аналогично тому, как это сделано в работе [2], предположим, что отношение поперечника первичного элемента к его длине равно 0,2. При этом для поперечного элемента с длиной 2,8 мм и коэффициентом упругости К3 = 0,02 г/мм получим значение модуля Юнга Е3 = 0,23 г/мм2=2,2 кН/м2. Аналогично для наклонного элемента с длиной 2,28 тт и коэффициентом упругости К1 = 0,9 г/мм получим Е1 = 12,6 г/мм2=123 кН/м2. Предполагая, что площадь боковой поверхности элемента представляет рабочую поверхность возникновения ньютоновской вязкости, при значении параметра вязкости Н1 = 16 с-г/мм, получим для коэффициента вязкости ^ = 11 с-кН/м2 и, соответственно, для Н3 = 1,5 с-г/мм = 0,8 с-кН/м2. Для сравнения приведем значения модуля Юнга эластина 600 кН/м2 и коэффициента вязкости вазелина 0,5 с-кН/м2.

Рис. 6. Сопоставление кривых гистерезиса, полученных в различных моделях, с данными эксперимента Л2038. Значения геометрических параметров и параметров упругости одинаковы во всех вариантах: /о1=2,157 мм, /о2=3,65 мм, /о3=3,2 мм, ^=0,9 мм. Синусоидальная нагрузка: (а) модель с вязким поперечным элементом ККН3 при К1=0,95 г/мм, К2=К3=0,027 г/мм, Н3=0,5 ст/мм, £3тах = -0,5, £тах = 0,243 ; (Ь) модель с вязким наклонным элементом ЯЬШ при , К1=1,05 г/мм, К2=К3=0,02 г/мм, Н1=5 ст/мм, £1тах = 0,06, £тах = 0,243 . Пилообразная нагрузка ЯЬШ: (с) К1=0,85 г/мм, К2=К3=0,03 г/мм, Н3=0,8 с-г/мм, £3тах =-0,455, £тах = 0,237, (ё) Н1=5 ст/мм,

£1тах = 0,06 , £тах = 0,243

2. Релаксации деформации

На рис. 7 а и Ь представлены кривые релаксации деформации наклонного, продольного и поперечного элементов в ромбической модели при заданной внешней силовой возрастающей и убывающей ступенчатой нагрузке. Нетрудно заметить, что времена релаксации первичных упругих элементов пропорциональны соответствующим коэффициентам вязкости, и в течение периода действия нагрузки релаксация менее вязких элементов практически завершается, тогда как более вязких продолжается. Это отчетливо видно при синусоидальной и пилообразной внешней нагрузке (рис. 7 с, ё).

Сравнение с экспериментом в этом случае не было проведено, поскольку имеющиеся опытные данные относятся к релаксации напряжения при ступенчатоизменяемой продольной деформации образца. Помимо рассмотрения вопроса об

Рис. 7. Кривые релаксации деформации различных элементов модели RLH13. Геометрические параметры и параметры жесткости элементов указаны в подписи к рис. 3, ^=1,0 ст/мм; Н3=0,01 ст/мм. Ступенчатая восходящая (а), нисходящая (Ь), синусоидальная (с) и пилообразная (ё) зависимости от времени задаваемой нагрузки У(/) (на вставках)

ударной жесткости, который возникает при анализе резких скачков деформации, исследование релаксации напряжения требует, на наш взгляд, перехода к моделям с включением демпфера последовательно упругим элементам и составит предмет отдельной публикации.

Обсуждение

Полученные в результате математического моделирования кривые гистерезиса соответствуют статической кривой сила-длина с учетом ее нелинейности в определенных условиях. Таким образом, представлено объяснение нелинейной формы петли гистерезиса, наблюдаемой на эксперименте без дополнительных предположений об экспоненциальной зависимости напряжения от деформации.

Раздвижка ветвей широкой части петли гистерезиса в разных диапазонах деформаций соответствует величине вязкости элемента, который вносит основной вклад в той или иной области деформации. В случае, когда параметры жесткости и вязкости первичных элементов существенно различаются (первые, как видно из подписей к рис. 3 и 4 на два порядка, тогда как вторые на порядок), кривая гистерезиса может иметь несимметричный вид с расширением в области или малых, или максимальных деформаций. При этом, на наш взгляд, при описании формы кривой гистерезиса нет необходимости в предположении зависимости параметра вязкости от длины препарата, как это сделано, например, в работе [2].

Установлено, что благодаря наличию параметров модели, описывающих упругость и вязкость элементов, имеется возможность варьирования формы кривой гистерезиса. Поскольку величины этих параметров могут быть оценены при сопоставлении кривых с экспериментальными данными [15], по изменению формы кривой гистерезиса можно судить о соотношении параметров, характерных для тех или иных структур, в реальных образцах мышц и об их изменении в различных условиях.

Неполное соответствие экспериментальных и теоретических данных по форме петли гистерезиса во всем диапазоне деформаций предполагает необходимость введения демпферов в поперечные и наклонные элементы структуры одновременно, а также и в последовательном их соединении. Если предположить, что физическим источником вязкости является фильтрация жидкости, то на основании этого вывода можно сделать предположение о присутствии механизмов фильтрации в миокарде как через мембрану (поперечные элементы), так и межклеточной жидкости (наклонные элементы).

Для большей определенности качественных и количественных выводов при сравнении теории и опыта, необходимо повышение точности эксперимента, которая пока еще недостаточна для оказания предпочтения той или иной модели.

Таким образом, показано, что простейший вариант учета вязкости в рассмотренных ранее двумерных моделях с помощью введения демпферов параллельно первичным упругим элементам позволяет получить нелинейные кривые гистерезиса, наблюдаемые в экспериментах, и кривые релаксации деформации.

Список литературы

1. Brady A.J., Tan S.T., Ricchiuti N.V. Perturbation measurements of papillary muscle elasticity // American Journal of Physiology. 1981. V. 241. P. H155-H173.

2. Anderson J., Li Z., Goubel F. Models of skeletal muscle to explain the increase in passive stiffness in desmin knockout muscle // Journal of Biomechanics. 2002. V. 35. P. 1315-1324.

3. Kaufman K.R., An K.N., Chao E.Y. Incorporation of muscle architecture into the muscle length-tension relationship // Journal of Biomechanics. 1989. V. 22. № 8-9. P. 943-948.

4. Protsenko Yu.L, Kobelev A.V., Kobeleva R.M., Routkevich S.M. The steady-state property ‘force-deformation’ for passive myocardium // Russian Journal of Biomechanics. 2001. V. 5. № 3. P. 30-40.

5. Kobeleva R.M., Kobelev A.V., Protsenko Yu.L. Mathematical details for passive stress-strain modeling of papillary muscle // International Symposium "Biological Motility: New Trends in Research". Pushchino, 2001. P. 68-69.

6. Protsenko Yu.L, Kobelev A.V., Kobeleva R.M., Routkevich S.M. Passive force-deformation non-linear stationary dependencies for myocardium // International Symposium "Biological Motility: New Trends in Research". Pushchino. 2001. P. 69-71.

7. Изаков В.Я., Мархасин В.С., Ясников Г.П., Белоусов В.С., Проценко Ю.Л. Введение в биомеханику пассивного миокарда. М.: Наука, 2000.

8. Loeffler L., Sagawa K. A one-dimensional viscoelastic model of cat heart muscle studied by small length perturbations during isometric contraction // Circulation. Research. 1975. V. 36. № 4. P. 498-512.

9. Noble M.M. The diastolic viscous properties of cat papillary muscle // Circulation Research. 1977. V. 40. № 4. P. 288-292.

10. Guide for the care and use of laboratory animals. Washington D.C.: Natural Academy Press, 1996.

11. Protsenko Yu.L., Routkevitch S.M., Zolotushnikov V.V., Abazov V.A. The experimental technique of myocardial rheological characteristics investigation // Russian Journal of Physiology. 1998. V. 84. № 1-2. P. 133-138.

12. Pinto J.G., Win R. Non-uniform strain distribution in papillary muscles // American Journal of Physiology. 1977. V. 233. № 3. P. H410-H416.

13. Fish D., Orenstein J., Bloom Sh. Passive stiffness of isolated cardiac and skeletal myocytes in the hamster // Circulation Research. 1984. V. 54. P. 267-276.

14. Weisfeldt M.L., Loeven W.A., Shock N.W. Resting and active mechanical properties of trabeculae carnea from aged male rats // American Journal of Physiology. 1971. V. 220. P. 1921-1927.

15. Kobelev A., Protsenko Yu., Kobeleva R. Modeling of non-linear viscous-elastic properties of myocardium specimens // Acta of bioengineering and biomechanics. 2002. V. 4. № 1. P. 498-499.

A.V. Kobelev, R.M. Kobeleva, Yu.L. Protsenko, I.V. Berman (Yekaterinburg, Russia; San Xose, USA)

NONLINEAR VISCOELASTIC BEHAVIOUR OF MYOCARDIUM FILAMENTS: ‘FORCE-LENGTH’ HYSTERESIS AND RELAXATION OF DEFORMATION

Dynamical behaviour of the myocardium preparations in passive state under periodic and step-wise loading is studied. Theoretical description uses part of identified rheological 2D models proposed earlier by us to analyze static “force-deformation” relationship in the myocardium. Viscous properties of the myocardium are modeled with the help of damps connected in parallel to some of primary elements of the model. The values of both elastic modules (as in the static case) and also viscous coefficients are estimated by comparison of theoretical ‘force-deformation’ hysteresis curves with experimental data obtained on cats’ and guinea pigs’ myocardium preparations. It is stated that the conclusion on what elements’ viscosity prevails in the heart filament can be made by the analysis of hysteresis curve shape. The relaxation rates of deformation for primary elements with different viscosity at given loading of step-wise, or periodical time dependence are obtained theoretically for the first time.

Key words: myocardium biomechanics, nonlinear ‘stress-strain’ curves, viscosity, hysteresis, deformation relaxation.

Получено 5 февраля 2003 года