Оптимизационная модель стационарной формы эритроцита в двумерном приближении Текст научной статьи по специальности «Физика»

УДК 531.2

DOI 10.18413/2075-4639-2019-51-1-135-144

ОПТИМИЗАЦИОННАЯ МОДЕЛЬ СТАЦИОНАРНОЙ ФОРМЫ ЭРИТРОЦИТА В ДВУМЕРНОМ ПРИБЛИЖЕНИИ

OPTIMIZATION MODEL OF THE ERYTHROCYTE STATIONARY SHAPE IN TWO-DIMENSIONAL APPROXIMATION

А.В. Голочалова, Э.Б. Кулумбаев A.V. Golochalova, E.B. Kulumbaev

Белгородский государственный национальный исследовательский университет, Россия, 308015, г. Белгород, ул. Победы, 85

Belgorod National Research University, 85 Pobeda St, Belgorod, 308015, Russia

E-mail: kulumbaev@bsu.edu.ru

Аннотация

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

Abstract

A model of conditional nonlinear minimization of the erythrocyte membrane's elastic energy has been developed for the calculation its stationary shapes in the two-dimensional approximation in normal conditions and in the existence of external influences. The equilibrium shape was calculated for a given volume and surface area of a normal red blood cell, which qualitatively agrees with the experimentally observed and with the results of a two-dimensional dynamic description of the erythrocyte shape relaxation. The external effects of osmotic pressure or the capillary walls during a red blood cell aspiration are simulated by imposing additional conditions in the form of constraints on the vector-argument of the objective function. The results of the calculation of deformation under the action of osmotic pressure or capillary aspiration of a red blood cell qualitatively correspond to observational data.

Ключевые слова: эритроцит, форма мембраны, упругая энергия, оптимизационная, модель. Keywords: erythrocyte, membrane shape, elastic energy, optimization model.

Введение

Эритроциты или красные кровяные тельца - это наиболее многочисленные клетки крови, концентрация которых примерно равна 5-106 mm-3, а занимаемый ими объем составляет 40-45 % от общего объема крови [Herman, 2007; Caro et al., 1978]. Поэтому динамика крови - суспензии форменных элементов (клеток) крови в плазме - в

наибольшей степени определяется механическими свойствами эритроцитов и в значительно меньшей степени зависит от свойств лейкоцитов и тромбоцитов, занимающих только 0,3 % и 0,15 % объема цельной крови соответственно. Тогда для установления адекватных реологических свойств крови, являющейся неньютоновской жидкостью [Caro et al., 1978], актуальными становятся экспериментальные и теоретические исследования движения и деформирования отдельных эритроцитов. В частности, значимым является нахождение равновесных (стационарных) форм эритроцита при различных внешних условиях.

Нормальный зрелый эритроцит человека состоит из эластичной мембраны и ограниченного ею почти насыщенного раствора гемоглобина [Herman, 2007]. Мембрана эритроцитов обеспечивает гомеостаз и функциональное состояние эритроцита, имеет сложную структуру, включая липидный бислой и подстилающую его спектриновую сеть. Толщина мембраны составляет ~10 nm, а ее массовая доля - не более 3 % от массы всей клетки [Caro et al., 1978].

Во взвешенном состоянии в плазме крови или в другом изотоническом растворе нормальный эритроцит имеет дисковидную двояковогнутую форму. «Усредненный» эритроцит имеет диаметр диска 7.8 цт, площадь поверхности 134 цт2, объем 94 цт3 [Evans, Fung, 1972]. «Целесообразность» создания природой эритроцита такой структуры, формы и размеров для системы переноса кислорода в организме аргументировано обсуждается в работе Атауллаханова и др. [2018].

Считается, что форма эритроцита в норме не деформирована и обусловлена свойствами его мембраны. Последнее подтверждается экспериментально наблюдаемым примерным восстановлением формы и размеров "тенями" эритроцитов - эритроцитов, первоначально лишенных гемоглобина в результате гемолиза и вновь заполненных подходящим изотоническим раствором [Caro et al., 1978]. Дисковидная двояковогнутая форма обеспечивает способность эритроцитов сильно деформироваться без изменения объема и площади поверхности в микрососудах (например, при диаметре 8 цт они могут проходить через капилляры диаметром 3 цт и длиной 12 цт [Caro et al., 1978]) и выполнять свою основную функцию - транспорт газов в системе кровообращения. Неизменность объема и площади поверхности эритроцита обусловлены малыми сжимаемостью жидкого клеточного содержимого и растяжимостью мембраны соответственно.

В настоящее время численное моделирование формы и деформирования эритроцита проводится обычно на основе континуального и дискретного подходов [Ju et al., 2015]. В рамках первого подхода мембрана считается двумерной сплошной средой, равновесная конфигурация которой при изотропной нагрузке поддерживается за счет внутренних натяжений, изгибающих моментов и поперечных сил, и рассчитывается на основе различных вариантов нелинейной теории тонкостенных оболочек [Pozrikidis, 2001]. Второй подход использует дискретное приближение мембраны системой упруго-связанных мезоскопических частиц, движение и равновесие которых рассчитывается на основе ньютоновской механики с учетом сил упругих деформаций растяжения, изгиба, диссипативной силы вязкого трения и штрафных сил, стремящихся упругим образом обеспечить заданные значения площади поверхности мембраны и ограниченного ею жидкого объема [Pivkin, Karniadakis, 2008].

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

минимуму упругой энергии изгиба мембраны при заданных площади ее поверхности и ограничивающем ею объеме [Canham, 1970]. Однако из-за ненадлежаще проведенного сравнения теории с экспериментальными данными [Marchenko, Podolyak, 2015] эта модель в дальнейшем неоправданно усложнялась учетом спонтанной кривизны [Deuling, Не1Й1с^ 1976], нелокального взаимодействия, сдвиговой упругости мембраны и цитоскелета [Mukhopadhyay, 2002]. В работе Марченко и Подоляк [Marchenko, Podolyak, 2015] показано, что принцип [Canham, 1970] количественно описывает наблюдаемую форму эритроцита без введения новых материальных параметров мембраны, неизвестных из опыта и поэтому являющихся по сути подгоночными. Следует отметить, что нелинейность соответствующих вариационных задач с неизбежностью обусловливает необходимость использования численных методов. Например, краевая задача для вариационного уравнения решается численным интегрированием методом стрельбы[Marchenko, Podolyak, 2015]. С этой точки зрения не проигрышным будет нахождение формы эритроцита в исходной постановке из минимизации упругой энергии изгиба мембраны при заданных площади поверхности и объеме эритроцита без сведения экстремальной задачи к вариационному уравнению. Такая постановка соответствует задаче условной нелинейной оптимизации, реализованной хорошо развитыми, апробированными и эффективными численными методами, например, в пакете Matlab [Messac, 2015]. Однако более существенным и значимым доводом для реализации такого подхода является оценка применимости принципа [Canham, 1970] для эритроцита, деформированного в результате внешнего воздействия.

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

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

В отсутствие внешних воздействий равновесная форма двумерного эритроцита -это замкнутая плоская упругая кривая Ь, которая минимизирует целевую функцию -упругую энергию изгиба этой кривой:

где Еь - упругая энергия изгиба, В - изгибная жесткость, к - кривизна элемента <1 плоской кривой.

Ограничения на искомую кривую ставятся в виде условий заданной площади 5"0, ею ограниченной, и заданного ее периметра Ь0:

которые соответствуют в двумерном приближении наблюдаемым свойствам неизменности объема и площади поверхности эритроцита.

Структура (1) и исключительно геометрический характер ограничений (2) обусловливают независимость искомой упругой кривой от значения В ее изгибной жесткости; значением В определяется только абсолютная величина Еь упругой энергии ее минимизирующей формы.

Метод исследования

(1)

(2)

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

В данной работе искомая непрерывная кривая L аппроксимируется дискретным набором из N точек с радиус-векторами ri (здесь и везде далее i —1,2,..., N). Тогда аргумент целевой функции - это вектор длины 2N, составленный из всех ri . Геометрическая интерпретация такой аппроксимации соответствует замене плоской кривой замкнутой ломаной (многоугольником), составленным из замкнутой цепочки

векторов l — r{1-^/N))i+1 - ri с фиксированной длиной | lt | —l0 — L0 / N (рис. 1). Здесь |_i/Nj -

\i +1, for i < N

неполное частное, поэтому

I I [i +1, for i < N

(1 "Li/Nj) i +1fori = N

что учитывает периодичность

индекса /. Тогда из «скоростного» смысла кривизны как отношения угла поворота орта касательной при смещении его начала на ё1 вдоль кривой к dl в используемом приближении кривой (рис. 1) следует, что

(3)

к =е/ /0, cos3 = ¡t_i ■ i / /2

а ограничения (2) на вектор-аргумент представляются равенствами:

1 N

1 z e ■ (-(

(1—Li / N j) i+1

+ ï)Jex ■ / ) = 5о, | /, | = /о,

(4)

i=1

где ex, ey - орты декартовой системы координат.

Рис. 1. Аппроксимация двумерной формы эритроцита N -угольником Fig. 1. Approximation of the two-dimensional shape of the erythrocyte by a polygon

Окончательно оптимизационная задача (1-2) с учетом дискретизации (3-4) записывается в координатном виде как минимизация целевой функции:

1 N

Еь(x1,xN, y^ yN) = -ЕЫ)Zee

2 1=1

где cos6> = [(х, - х,

(1-[г / N J) i+1

X) + (у, - У,-1)(y

41-Li / V J),+1

y,)]/12,

Eb0 = B /l0

характерное значение упругой энергии изгиба с ограничениями

1 V VV

(y(1-Li/nJ)i+1 + yi)(x(1-Li/NJ)i+1 - xi) = S0; (x(1-Li/VJ),+1 - xi) + (y(1-Li/nJ),+1 - yi) = l0 . (6)

2

i=1

Здесь xi, yi - абсцисса и ордината i - й вершины N -угольника.

Оптимизация осуществляется численно в пакете Matlab.

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

Равновесная форма эритроцита доставляется решением оптимизационной задачи (5 - 6), значения параметров в которой после задания числа N вершин многоугольника выбираются равными: l0 = L0/N - длина стороны N -угольника; L0 = 2kR0 - периметр двумерного эритроцита, равный длине окружности радиуса R = 3 цт; S0 = kR -площадь эритроцита, которая составляет только часть ( s» = 0.46 ) максимально возможной площади, ограниченной замкнутой кривой с периметром L0 - круга радиуса R. Такой набор значений устанавливается из соответствия L0 и S0 двумерного эритроцита аналогичным характеристикам осевого сечения усредненного осесимметричного эритроцита (рис. 2) - экспериментально наблюдаемой формы, полученной в работе [Evans, Fung, 1972] путем усреднения измерений геометрических размеров нормального эритроцита и представленной в виде формулы, образующей форму линии в цилиндрической системе координат (r, р, z) :

z(r ) = - r2 (c0 + cr 2 + c2r4), 0 < r < 1,

где c0 = 0.207161; c = 2.002558; c2 = -1.122762 - эмпирические константы усреднения измерений; r = r / R, z = z / R - безразмерные цилиндрические координаты, нормированные на радиус R = 3.91 prn эритроцита; z(r) - это половина образующей, которая дополняется ее зеркальным отражением относительно радиальной оси.

Рис. 2. Равновесная форма эритроцита: наблюдаемая (сплошная кривая) по данным [Evans,

Fung, 1972] и расчетная (тонкая кривая с N = 76 узлами) по модели (5-6) при s = 0.46. Fig. 2. The equilibrium shape of the erythrocyte: observed (solid curve) [Evans, Fung, 1972] and predicted (thin curve with N = 76 nodes) by the model (5-6) for s = 0.46.

Результаты расчета для N — 76 (таким часто задается число мезоскопических частиц при динамическом описании двумерной формы эритроцита [Barns et al., 2017]) представлены на рис. 2. Видно, что равновесная форма качественно соответствует экспериментально наблюдаемой, а количественно - рассчитанные радиус и толщина в центре вогнутости эритроцита завышены не более чем на 10 и 30 % соответственно. Аналогичного качества результаты получаются и при динамическом описании двумерной формы эритроцита [Nayanajith et al., 2016].

Влияние числа вершин аппроксимирующего форму эритроцита многоугольника, которое оценивается тестовыми расчетами с увеличенным и уменьшенным в два раза значением N — 76 , оказывается несущественным, так как результаты расчета изменяются мало настолько, что визуально практически неразличимы на рис. 2.

При динамическом описании двумерной формы эритроцита упругая энергия изгиба N одинаковых упруго-связанных мезоскопических частиц рассчитывается по формуле [Barns et al., 2017]:

1 N ( й \

Е = 1К X tg2 1, (7)

2 ¿=1 V2 У

где кь - упругая константа деформации изгиба.

Результаты тестового расчета с целевой функцией (7) и релаксации упруго связанной системы мезоскопических частиц при динамическом описании мембраны [Голочалова, 2017] согласовываются между собой и с равновесной формой, приведенной на рис. 2, для 5* = 0.46 , до визуальной неразличимости в масштабе рисунка.

Формоизменение эритроцита под действием осмотического давления. Изменение формы эритроцита в зависимости от параметра 5 = {0.3; 0.7; 0.9} иллюстрируется на рис. 3. Видно, что радиус эритроцита уменьшается, а его толщина в центре вогнутости увеличивается с ростом 5 , причем сама вогнутость реализуется в диапазоне значений 5 , ограниченном сверху значением между 0,7 и 0,9. Поскольку 5 задает фактический «недобор» площади двумерного эритроцита от максимально возможной площади, ограниченной замкнутой кривой с тем же периметром, то влияние увеличения (уменьшения) 5 на результаты расчета интерпретируются как формоизменение эритроцита при его осмотическом набухании (сжимании). Аналогичные результаты получаются и при динамическом описании мембраны эритроцита [Nayanajith et а1., 2016].

Формоизменение эритроцита при капиллярной аспирации. Рассматривается формоизменение эритроцита при его «засасывании» в капиллярный канал. В оптимизационной модели аспирация и капилляр формируются дополнительными к (6) ограничениями:

Уп-1 = ^ - =0; Уп+1 = ^ (8а)

у—к — с < 0, (8б)

где п - номер вершины N -угольника, совпадающей с центром вогнутости формы эритроцита в исходном положении, d - глубина проникновения в капилляр соседних с ней (п — 1) и (п +1 )-й вершин при аспирации, к, с - параметры геометрии стенок капилляра, которые задают фактическую область определения аппроксимирующего эритроцит N -угольника.

Результаты расчета формоизменения эритроцита при его «засасывании» в капилляры с параметрами к = 5, с = 5 и с = 5 -103 приводятся на рис. 4 и 5 соответственно.

з

0.9

Е

0

si

-2

-1

2

-3

-4 -3 -2 -1 0

х, fj, m

2

3

4

Рис. 3. Стационарная форма эритроцита при s» = 0,3; 0,46; 0,7; 0,9 Fig. 3. The stationary shapes of the erythrocyte for s = 0,3; 0,46; 0,7; 0,9

Перед аспирацией оптимизационная модель (5-6, 8б) доставляет равновесную, недеформированную, форму эритроцита в нижней полуплоскости. Эритроцит только касается образующих капилляр горизонтальных стенок (рис. 4 а, 5 а). Это положение эритроцита считается исходным перед его проникновением в капилляр.

При аспирации оптимизационная модель (5-6, 8а, 8б) определяет форму эритроцита для каждой заданной глубины его проникновения в капилляр (рис. 4 Ь-Г 5 Ь-1). Видно, что изменение формы зависит от размера капилляра. В узком капилляре (с шириной около 2 цш) вначале, при й < 3.3 цт, формируется звездчатая (трехконечная) форма эритроцита (рис. 4 Ь-ё). При дальнейшем углублении на 0.1 цш форма эритроцита становится снова двояковогнутой, как и в исходном положении в недеформированном состоянии перед аспирацией (рис. 4а), но повернутой на п/2 (рис. 4е). Далее, при й > 3.4 цт, эритроцит «углубляется» торцом и форма его не изменяется (рис. 41).

В более широком капилляре (с шириной около 4 цш) с самого начала реализуется куполообразная форма (или форма парашюта) эритроцита (рис. 5Ь). С увеличением й она деформируется к подковообразной (рис. 5 с-1), а далее качественной перестройки формы при аспирации эритроцита в такой капилляр не наблюдается.

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

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

Заключение

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

5 4 3

5 4

3 2

1

Е

а- о si

-1 -2 -3 -4 -5

Е

а.

4 3 2 1

Е

а- 0

si

-1

-2 -3 -4

(a) J

0 x, fim 5

■ (c) ■ ■ J V

6 ■

-5 0 x, (¡m

2

1

Е

а- о

Ss

-1 -2 -3 -4 -5

5 4 3 2 1

Е

а. о г;

-1 -2 -3 -Л -5

5

- (e) 4 (f)

3

2

J V 1 E

: a- 0 si -1 -2

-3

w -4

(b)

О

х, рт

■ (d)

х, дт

0

х, рт

Рис. 4. Форма эритроцита перед (а) и при аспирации в зависимости от глубины проникновения d

= 1 (b); 2 (c); 3 (d); 4 (e); 5 (f) цт в капилляр (8б) с параметрами к — 5, c — 5 . Fig. 4. The shape of the erythrocyte before (a) and during aspiration on the depth of penetration d = 1 (b); 2 (c); 3 (d); 4 (e); 5 (f) цт in the capillary (8б) with parameters к — 5, c — 5 .

Рис. 5. Форма эритроцита перед (а) и при аспирации в зависимости от глубины проникновения d

= 1 (b); 2 (c); 3 (d); 4 (e); 5 (f) цт в капилляр (8б) с параметрами k = 5, c = 5 ■ 103. Fig. 5. The shape of the erythrocyte before (a) and during aspiration on the depth of penetration d = 1 (b); 2 (c); 3 (d); 4 (e); 5 (f) цт in the capillary (86) with parameters k = 5, c = 5 ■ 103.

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

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

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

1. Атауллаханов Ф.И., Борсакова Д.В., Протасов Е.С. и др. 2018. Эритроцит: мешок с гемоглобином или живая, активная клетка? Вопросы гематологии/онкологии и иммунопаталогии в педиатрии. 17(1): 108 - 116.

Ataullakhanov F.I., Borsakova D.V., Protasov E.S. i dr. 2018. Eritrotsit: meshok s gemoglobinom ili zhivaya, aktivnaya kletka? Voprosy gematologii/onkologii i immunopatalogii v pediatrii. 17(1): 108 - 116.

2. Вольмир А.С. 1967. Устойчивость деформируемых систем. - М.: Наука, 984.

Vol'mir A.S. 1967. Ustoychivost' deformiruemykh sistem. - M.: Nauka, 984.

3. Голочалова А.В. 2017. Особенности расчета двумерной формы эритроцита человека методом частиц. Естественнонаучные, инженерные и экономические исследования в технике, промышленности, медицине и сельском хозяйстве. Материалы I молодежной научно-практической конференции с международным участием. - Белгород ИД «Белгород» НИУ «БелГУ», 409 - 412.

Golochalova A.V. 2017. Osobennosti rascheta dvumernoy formy eritrotsita cheloveka metodom chastits. Estestvennonauchnye, inzhenernye i ekonomicheskie issledovaniya v tekhnike, promyshlennosti, meditsine i sel'skom khozyaystve. Materialy I molodezhnoy nauchno-prakticheskoy konferentsii s mezhdunarodnym uchastiem. - Belgorod ID «Belgorod» NIU «BelGU», 409 - 412.

4. Barns S. et al. 2017. Investigation of red blood cell mechanical properties using AFM indentation and coarse-grained particle method. BioMedical Engineering OnLine. 16:140.

5. Canham P.B. 1970. The minimum energy of bending as a possible explanation of the biconcave shape of the human red blood cell. J. Theor. Biol. 26(1): 61-81.

6. Caro C.G., Pedley T.J., Schroter R.C., Seed W.A. 1978. The mechanics of the circulation. Oxford: Oxford University Press, 527.

7. Deuling H.J., Helfrich W. 1976. Red blood cell shapes as explained on the basis of curvature elasticity. Biophysical Journal. 16(8): 861-868.

8. Evans E., Fung Y.-C. 1972. Improved measurements of the erythrocyte geometry. Microvascular Research. 4(4): 335 - 347.

9. Herman I.P. 2007. Physics of the Human Body. Springer - Verlag Berlin Heidelberg, 992.

10. Ju M. et al. 2015. A review of numerical methods for red blood cell flow simulation. Computer Methods in Biomechanics and Biomedical Engineering. 18(2): 130-140.

11. Marchenko V.I., Podolyak E.R. 2015. On the equilibrium shape of erythrocytes. JETP. 120(4): 751 - 752.

12. Messac A. 2015. Optimization in practice with Matlab. - Cambridge University Press, 465.

13. Mukhopadhyay R., Lim H.W.G., Wortis M. 2002. Echinocyte Shapes: Bending, Stretching, and Shear Determine Spicule Shape and Spacing. Biophysical Journal. 82(4): 1756-1772.

14. Nayanajith P.G.H. et al. 2016. A coupled SPH-DEM approach to model the interactions between multiple red blood cells in motion in capillaries. International Journal of Mechanics and Materials in Design. 12(4): 477-494.

15. Pivkin I.V, Karniadakis G.E. 2008. Accurate coarse-grained modeling of red blood cells. Phys. Rev. Lett. 101(11): 118105.

16. Pozrikidis C. 2001. Effect of membrane bending stiffness on the deformation of capsules in simple shear flow. J Fluid Mech. 440: 269 - 291.

17. Skalak R., Branemark P.-I. 1969. Deformation of Red Blood Cells in Capillaries. Science 164(3880):717-719.