Exoskeleton Dynamics Simulation with the System of Three Variable-Length Links of Adjustable Stiffness

1. Введение. Несмотря на значительно увеличивающееся количество публикаций и патентов по теме “экзоскелет” [1], остаются актуальными проблемы их использования. Одной из причин недостаточной комфортабельности является наличие массивных абсолютно твердых звеньев, оказывающих существенное влияние на динамику человеко-машинной системы. В реальном опорно-двигательном аппарате человека расстояния между осями вращения в суставах изменяются [2–4], что должно учитываться при проектировании экзоскелетов. Помимо отсутствия изменения длины звена, сама масса звена экзоскелета, сопоставимая по порядку величины с массой человека, оказывает значительное влияние на движения человека в экзсокелете [5–9]. Поэтому требуется адаптация и обучение человека движениям в экзоскелете. Оценке пользы от применения современных моделей экзоскелетов, затрат энергии у человека и ее экономии от использования экзоскелетов посвящены работы [10–23]. В статье [24] представлен процесс создания структур для отображения необходимых исследований, которые помогут разработчикам экзоскелетов. В данной статье развивается идея применения в конструкции экзоскелетов звеньев переменной длины [1, 5–9]. Предлагается управление длиной участка переменной длины с использованием магнитно-реологической жидкости, путем изменения ее вязкости под действием прилагаемого внешнего магнитного поля. Научных публикаций с описанием применения магнитно-реологических сред в мехатронных антропоморфных устройствах не так много [25–30]. Имеется ряд работ, в которых описываются свойства магнитно-реологических жидкостей и их различные применения [31–45]. В данном списке работы размещены в хронологическом порядке и видно, что первые подходящие модели, которые впоследствии были использованы при описании магнитно-реологических жидкостей, созданы достаточно давно. Константы, необходимые для конструирования модели, рассматриваемой в данной работе, взяты из приведенного списка публикаций. На основе этих данных получена зависимость силы сопротивления магнитно-реологической жидкости от внешнего магнитного поля, необходимая для управления динамикой звена переменной длины с регулируемой жесткостью.

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

3. Модель экзоскелета, состоящая из трех звеньев переменной длины с регулируемой жесткостью. Рассмотрим 3D-модель экзоскелета, состоящего из трех звеньев переменной длины с регулируемой жесткостью в одноопорной фазе движения. Пусть первое, нижнее звено А₀А₁ моделирует голень опорной ноги. Второе, среднее звено А₁А₂ – бедро опорной ноги. Третье, верхнее звено А₂А₃ – корпус экзоскелета. В точках механизма А₂ и А₃ имеются сосредоточенные массы m_P, и m_H, приближенно моделирующие переносимую ногу и голову совместно с руками соответственно. Следует отметить, что если представление голени и бедра опорной ноги экзоскелета тонкими звеньями, является достаточно адекватной аппроксимацией реальных звеньев опорно-двигательного аппарата человека, то моделирование корпуса тонким стержнем, а переносимой ноги, головы и рук точечными массами является весьма грубым приближением, позволяющим провести лишь качественный анализ динамических свойств предлагаемой конструкции экзоскелета. При составлении уравнений движения используем метод, основанный на применении локальных систем координат, предложенный в работе [46]. Введем неподвижную декартову систему координат А₀x₀y₀z₀ (рис. 1), c началом в точке А₀, где закреплен шарнир.

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

В модели используются различные модели шарниров и их комбинации. Так, в шарнире А₀ имеется комбинация двух цилиндрических шарниров с взаимно ортогональными осями вращения. Такая конструкция соответствует движениям в голеностопном суставе человека, в котором при ходьбе практически не наблюдается поворота голени вокруг своей продольной оси. Поэтому в модели голеностопного сустава не используется сферический шарнир. В точке А₁ находится цилиндрический шарнир, моделирующий коленный сустав человека. Тазобедренный сустав человека является сферическим шарниром, поэтому в точке А₂ расположен сферический шарнир, моделирующий все имеющиеся реальные степени подвижности. При составлении системы дифференциальных уравнений применим подвижные локальные системы координат А₀x₁y₁z₁, А₁x₂y₂z₂ и А₂x₃y₃z₃, жестко связанные с шарнирами А₀, А₁ и А₂, в которых реализуются повороты звеньев. Такой способ описания использует углы между звеньями, соответствующие режиму работы реальных приводов роботизированного экзоскелета, посредством которых изменяются относительные углы. Оси z₁, z₂ и z₃ направим вдоль подвижных звеньев А₀А₁, А₁А₂ и А₂А₃. На (рис. 1) показаны углы поворотов, предусмотренные в модели.

Звенья одинаковой конструкции А₀А₁, А₁А₂ и А₂А₃ состоят из абсолютно твердых штоков А₁D₁, А₂D₂ и А₃D₃ с поршнями D₁, D₂ и D₃, находящимися внутри корпусов B₁C₁, B₂C₂ и B₃C₃. Корпуса являются тонкостенными цилиндрами с намотанными на них электромагнитными катушками. Цилиндры внутри заполнены магнитно-реологической жидкостью. Снизу к цилиндрам жестко прикреплены абсолютно твердые стержни А₀B₁, А₁B₂ и А₂B₃, на концах которых расположены шарниры. В шарнирах предполагается возможность создания необходимых управляющих моментов.

Каждое звено имеет по два весомых абсолютно жестких стержня: A₀B₁ = l₁₁, D₁A₁ = l₁₂, A₁B₂ = l₂₁, D₂A₂ = l₂₂, A₂B₃ = l₃₁, D₃A₃ = l₃₂, где первый индекс равен номеру звена, второй – номеру участка на звене. Длины цилиндров с магнитно-реологической жидкостью: B₁C₁ = l₁₃, B₂C₂ = l₂₃, B₃C₃ = l₃₃. Их диаметры будем считать пренебрежимо малыми в сравнении с длиной, т.е. допустим, что цилиндр с магнитно-реологической жидкостью моделируется стержнем. Массы отдельных элементов звеньев равны: m₁₁, m₁₂, m₁₃ соответственно. Моменты инерции стержней относительно осей вращений соответственно равны I₁₁, I₁₂, I₁₃. В первом приближении допустим, что сами поршни D₁, D₂ и D₃ имеют пренебрежимо малую массу и моменты инерции по сравнению со стержнями D₁A₁ , D₂A₂ и D₃A₃. Вычисление тензора инерции цилиндра с поршнем, заполненного магнитно-реологической жидкостью, требует отдельного исследования, в данной работе ограничимся существенно приближенным моделированием в целом цилиндра стержнем. Длины звеньев изменяются вследствие относительного движения поршней со штоками D₁A₁, D₂A₂ и D₃A₃ внутри цилиндров с магнитно-реологической жидкостью B₁C₁, B₂C₂ и B₃C₃ вдоль направлений звеньев А₀A₁, А₁A₂ и А₂A₃. Магнитно-реологическая жидкость реализует силы сопротивления F₁, F₂ и F₃, действующие вдоль штоков с поршнями. В результате действия внешнего магнитного поля магнитные частицы частично ориентируются вдоль силовых линий магнитного поля, причем с увеличением напряженности магнитного поля, доля ориентированных частиц увеличивается, тем самым увеличивается сила сопротивления магнитно-реологической жидкости движению поршня. Следовательно, обеспечивается управляемое изменение длины звена в соответствии с фазой ходьбы и нагрузкой от соседних звеньев.

Положения звеньев экзоскелета однозначно определяются углами α₁(t), β₁(t), α₂(t), α₃(t), β₃(t), γ₃(t) и переменными длинами участков звеньев между поршнем и дном цилиндра B₁D₁ = ξ₁(t), B₂D₂ = ξ₂(t) и B₃D₃ = ξ₃(t) (рис. 1). Следовательно, модель имеет девять степеней свободы. Управляющие моменты, развиваемые в шарнирах А₀, A₁ и A₂ с идеальными связями, обозначим M_1α, M_1β, M_2α, M_3α, M_3β, M_3γ и будем считать мгновенно развивающими требуемое усилие.

Кинетическая энергия рассматриваемого механизма вычисляется путем интегрирования по всем весомым участкам A₀B₁, B₁C₁, D₁A₁, A₁B₂, B₂C₂, D₂A₂, A₂B₃, B₃C₃, D₃A₃.

$T=\frac{1}{2}\left(m_P{V}_{A_2}^2+m_H{V}_{A_3}^2+{\displaystyle \sum _{i=1}^3{\displaystyle \sum _{j=1}^3{\displaystyle \underset{0}{\overset{l_{ij}}{\int }}{\rho }_{ij}{V}_{ij}^{2}d{\zeta }_{ij}}}}\right)$. (3.1)

Здесь: i – номер звена (i = 1,2,3); j – номер участка на звене, причем j = 1 соответствует нижним участкам A₀B₁, A₁B₂ и A₂B₃, j = 2 соответствует штокам с поршнями D₁A₁, D₂A₂ и D₃A₃, j = 3 соответствует цилиндрам B₁C₁, B₂C₂ и B₃C₃; ζ_ij – координата бесконечно малой частицы j-того инерционного участка i-того звена, ρ_ij – плотность j-того участка i-того звена, при этом m_ij = ρ_ijl_ij, I_ij = ρ_ijl_ij³/3, l_ij, m_ij, I_ij – длина, масса, момент инерции j-того участка i-того звена, m_P – масса переносимой ноги, расположенная в точке А₂, $V_{A_2}^2$ – квадрат скорости точки А₂, моделирующей тазобедренный сустав, m_H – масса головы и рук, расположенная в точке А₃, $V_{A_3}^2$ – квадрат скорости точки А₃, моделирующей голову человека, $V_{ij}^2$ – квадрат скорости бесконечно малой частицы j-того участка i-того звена.

Кинетическая энергия (3.1) предложенной модели громоздка, вследствие того, что звенья совершают сложное движение вокруг шарниров, поэтому приведем формулу кинетической энергии первого звена

$\begin{array}{c}{T}_1=\frac{1}{2}\left[\left(I_{11}+I_{12}+I_{13}+\left(l_{11}^2+l_{11}{l}_{12}+2l_{11}{\xi }_1+l_{12}{\xi }_1+{\xi }_1^2\right)m_{12}+\left(l_{11}^2+l_{11}{l}_{13}\right)m_{13}\right)\times \right.\\ \left.\times \left({\dot{\alpha }}_1^2{\cos }^2{\beta }_1+{\dot{\beta }}_1^2\right)+m_{12}{\dot{\xi }}_1^2\right].\end{array}$

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

Рассмотрим поведение магнитно-реологической жидкости. Вследствие приложения сил со стороны соседних стержней, действия реакции со стороны неподвижной опоры и силы тяжести происходит относительное движение поршней D₁, D₂ и D₃ в цилиндрах B₁C₁, B₂C₂ и B₃C₃, содержащих магнитно-реологическую жидкость, которая, вследствие изменения вязкости оказывает различное сопротивление движению поршня. В отсутствие внешнего магнитного поля поведение частиц магнитно-реологической жидкости хаотичное. При этом магнитно-реологическая жидкость не оказывает значительного сопротивления движению поршня. На (рис. 1) показаны ориентированные магнитные частицы в случае приложения внешнего магнитного поля. Изменением вязкости магнитно-реологической жидкости управляет внешнее магнитное поле, реализуемое при помощи намотанных на цилиндр катушек. При приложении внешнего магнитного поля жидкость оказывает сопротивление движению поршня, которое тем больше, чем больше ориентированных частиц в ней содержится. Становится возможным регулировать жесткость звеньев А₀А₁, А₁А₂ и А₂А₃ при изменении их длины. Управление изменением вязкости магнитно-реологической жидкости реализуется с использованием внешнего магнитного поля, что приводит к изменению силы сопротивления движению поршня внутри цилиндра. Зависимость силы сопротивления, реализуемой магнитно-реологической жидкостью от напряженности магнитного поля, имеет следующий вид [35, 37]:

$\begin{array}{c}F=\pi r^2\left\{\left[1+{\varphi }_L\left({\alpha }_n+\frac{\left({\zeta }_n+{\beta }_n{\lambda }_n+{\beta }_n\cos 2\theta \right)}{2}+\left({\chi }_n-2{\beta }_n{\lambda }_n\right)\frac{{\sin }^{2}2\theta }{4}\right)\right]{\eta }_f\dot{\gamma }+\right.\\ \left.+\frac{9{\varphi }_L{\mu }_0{H}^2\left(n-1\right)\left({\chi }_f+1\right)\sin 2\theta }{16n_c}\right\}.\end{array}$ (3.2)

Здесь: ${\varphi }_L$ – объемная доля микронных частиц, равная ${\varphi }_L$ = 0.0127, а n_с – максимальное число частиц в цепочке n_с = 50, ${\alpha }_n$, …, ${\lambda }_n$ – кинетические коэффициенты, $v_L$ – объем микронной частицы, $\dot{\gamma }$ – скорость сдвига, ${\eta }_f$ – вязкость магнитно-реологической жидкости ${\eta }_f$ = 1.36 Па с, $\theta$ – угол между вектором, направленным вдоль прямой, соединяющей центры двух микронных частиц и вектором напряженности внешнего магнитного поля $\theta =0.37$ рад [35, 37].

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

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

${\alpha }_1\left(t\right)=a_1\pi +b_1{c}_1\sin \left[d_1-\pi \left(1-\cos \omega t\right)/2\right]$ (4.1)

${\alpha }_2\left(t\right)=a_2\pi +b_2{c}_2\sin \left[d_2-\pi \left(1-\cos \omega t\right)/2\right]$ (4.2)

${\alpha }_3\left(t\right)=a_3\pi +b_3{c}_3\sin \left[d_3-\pi \left(1-\cos \omega t\right)/2\right]$ (4.3)

${\beta }_1\left(t\right)=a_4\pi +b_4{c}_4\sin \left[d_4-\pi \left(1-\cos \omega t\right)/2\right]$ (4.4)

${\beta }_3\left(t\right)=a_5\pi +b_5{c}_5\sin \left[d_5-\pi \left(1-\cos \omega t\right)/2\right]$ (4.5)

${\gamma }_3\left(t\right)=a_6\pi +b_6{c}_6\sin \left[d_6-\pi \left(1-\cos \omega t\right)/2\right]$ (4.6)

${\xi }_1\left(t\right)=e_1{l}_1-l_{11}-l_{12}+l_1{l}_{\xi }\cos 2\omega t$ (4.7)

${\xi }_2\left(t\right)=e_2{l}_2-l_{21}-l_{22}+l_2{l}_{\xi }\cos 2\omega t$ (4.8)

${\xi }_3\left(t\right)=e_3{l}_3-l_{31}-l_{32}+l_3{l}_{\xi }\cos 2\omega t$. (4.9)

Параметры механизма и ходьбы приняты равными: $t_k=0.36\text{ c}$ – время одноопорной фазы движения, $\omega =\pi /t_k\approx 8.73$ рад/с – угловая скорость, параметры ходьбы a₁ = a₄ = a₅= = a₆ = 1, a₂ = 0.9, a₃ = 0.5, b₁ = b₆ = 0.1, b₂ = 0.12, b₃ = 0.2, b₄ =1.0, b₅ = 0.5, c₁ = c₃ = c₄ = 0.25, c₂ = c₅ = c₆ = 0.279, d_i = π/2 (i = 1,2,…,6) выбирались из условия синтеза движения звеньев экзоскелета, близкого к антропоморфному; e₁ = e₂ = 1.3, e₃ = 1.2, $l_{\xi }=0.05$ – коэффициенты изменения длины звена, l₁ = 0.385 м, l₂ = 0.477 м, l₃ = 0.771 м – длины звеньев: голени, бедра, корпуса, соответствующие реальным данным человека [47]. Длины на звене распределялись следующим образом: l₁₁ = l₁₂ = l₁/2, l₂₁ = l₂₂ = l₂/2, l₃₁ = l₃₂ = l₃/2.

Параметры механизма, используемые в уравнениях движения следующие [47]: массы звеньев $m_1=2.91$ кг, $m_2=8.93$ кг, $m_3=23.86$ кг распределялись поровну между двумя весомыми абсолютно твердыми участками звеньев, т.е. $m_{i1}=m_{i2}=m_i/2$ (i = 1,2,3). Массы цилиндров с магнитно-реологической жидкостью считаем одинаковыми и равными: $m_{i3}=1.6$ кг (i = 1,2,3). Точечную массу, расположенную в шарнире А₂ и моделирующую переносимую ногу с учетом стопы примем равной m_P = 18.62 кг, а расположенную в точке А₃ массу головы и рук m_H = 10,62 кг. Моменты инерции весомых участков звеньев относительно осей, проходящих через их нижние точки равны: I₁₁ = 0.018 кг м², I₂₁ = 0.085 кг м², I₃₁ = 0.591 кг м², I_i1 = I_i2 (i = 1,2,3). Длины цилиндров с магнитно-реологической жидкостью для расчета момента инерции принимались равными l_i3 = 0.15 м, а моменты инерции I_i3 = 0.012 кг м² (i = 1,2,3). Ускорение свободного падения g = 9.81 м/с².

Графики зависимостей первых производных от времени, задаваемых формулами (4.1)–(4.9) представлены на рис. 3.

Рис. 2. Зависимости задаваемых формулами (4.1)–(4.9) обобщенных координат от времени.

Рис. 3. Зависимости обобщенных скоростей от времени.

Зависимости вторых производных от времени, задаваемые формулами (4.1)–(4.9) представлены на рис. 4.

Рис. 4. Зависимости обобщенных ускорений от времени.

В результате алгебраического решения относительно управляющих моментов и продольных сил системы уравнений движения для рассматриваемой модели (рис. 1) при задании движения выражениями (4.1)–(4.9), найдены зависимости управляющих моментов от времени, представленные на рис. 5. Момент M_3γ = 0 во все время движения, т.к. при моделировании тонким стержнем и отсутствии звена сверху не возникает усилий при вращении стержня вокруг продольной оси.

Рис. 5. Зависимости управляющих моментов в шарнирах экзоскелета от времени.

Зависимости управляющих продольных сил от времени представлены на рис. 6.

Рис. 6. Зависимости управляющих продольных сил в звеньях экзоскелета от времени.

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

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

Рис. 7. Зависимость управляющей жесткостью звена напряженности внешнего магнитного поля.

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

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

$A={\displaystyle \sum _{i=1}^k{\displaystyle {\int }_0^T\left|M_i{\dot{\varphi }}_i\right|dt}}$ (5.1)

$A={\displaystyle {\int }_0^T{\displaystyle \sum _{i=1}^n\left|Q_i{\dot{q}}_i\right|dt}}$. (5.2)

Здесь T – время движения (T = t_k), $M_i$ – моменты, развиваемые приводами аппарата, k – количество независимых приводов в шарнирах (k = 6), $Q_i,{\dot{q}}_i$ – обобщенные силы и обобщенные скорости, n – количество участков переменной длины (n = 3). Знак модуля используется вследствие отсутствия механизма рекуперации энергии. Так как затраты энергии подсчитываются как интегралы, т.е. суммы бесконечно малых работ, то, если бы не было знака модуля, отрицательные работы при торможении звеньев вычитались бы так, если бы приводы работали в режиме генераторов энергии, что не предусмотрено конструкцией данной модели экзоскелета.

В результате применения формул (5.1)–(5.2) к расчету энергетических затрат приводов, осуществляющих повороты звеньев при найденных управляющих моментах и силах, получаем: А = 260.1 Дж, являющиеся суммарными энергозатратами механизма.

6. Анимационная визуализация движения модели опорной ноги экзоскелета. Для наглядного представления о движении модели полезно строить анимацию движения механизма. Для предложенной модели приведем несколько кадров (рис. 8) из анимационного ролика, в котором механизм движется в соответствии с тем, как заданы изменения углов и длин звеньев формулами (4.1)–(4.9).

Рис. 8. Кадры анимации движения механизма.

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

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

Работа выполнена при финансовой поддержке за счет гранта Российского научного фонда № 22-21-00491, https://rscf.ru/project/22-21-00491/.