CC BY
66
НАУЧНОЕ ИЗДАНИЕ МГТУ ИМ. Н. Э. БАУМАНА
НАУКА и ОБРАЗОВАНИЕ
Эл № ФС77 • 48211. Государственная регистрация №0421200025. ISSN 1994-0408
электронный научно-технический журнал
Численное моделирование операции стернотомии с помощью
хирургического робота-манипулятора при движении по
заданной траектории
# 05, май 2013
DOI: 10.7463/0513.0599870
Гуськов А. М., Кузнецова М. С.
УДК 531, 53.072, 617-089
Россия, МГТУ им. Н.Э. Баумана gouskov [email protected] [email protected]
Введение
Применение роботизированных хирургических систем имеет много преимуществ. Роботизированные системы предлагают хирургам трехмерное изображение операционного поля и ощущение глубины, а также предоставляют в их распоряжение инструменты, которые в отличие от человеческих рук позволяют проводить длительные операции [1]. Если говорить о послеоперационных обострениях, то риски осложнений, связанных с искусственной вентиляцией легких, инфицированием раны и операционным травматическим стрессом, определяются фактором времени выполнения каждого этапа оперативного вмешательства [2]. Так как по сравнению с человеческой рукой инструменты роботизированных систем имеют управляемые с высокой точностью позиционирования степени свободы и меньшие размеры, то достигается большая эффективность операций. Результаты операций менее зависят от квалификации проводящего их специалиста, тем самым снижается так называемый «человеческий фактор» при проведении длительных и требующих максимальной точности хирургических вмешательств.
В данной работе сделана попытка смоделировать процесс стернотомии (от др.-греч. отеруоу - грудина и тоц^ - разрез, рассечение) -хирургического маневра, заключающегося в рассечении грудины [2]. Стернотомия выполняется для обеспечения доступа к органам и патологическим образованиям переднего средостения: сердцу, отходящим от него
крупным кровеносным сосудам и др. Задача численного моделирования стернотомии решена в программном пакете MATLAB® Robotics Toolbox + Simulink®.
Нужно отметить, что при проведении робот-ассистированных операций управление роботом осуществляется исходя из усилий, действующих на инструмент. Данные усилия, возникающие при разрезании ткани, могут вызывать перемещения инструмента в направлениях, не совпадающих с направлением приложения рабочей силы. Это явление принято называть координатной связью [3].
Если инструмент в процессе операции отклоняется от прогнозируемой траектории, то в оперируемой ткани могут возникнуть серьезные нарушения. Поэтому оценка и учет координатной связи является одной из важнейших задач численного моделирования и проектирования робот-ассистированных операций [4]. В работе рассмотрена рабочая конфигурация робота, для которой произведена оценка перемещения конечной точки
робота-манипулятора под действием нагрузки и вычислена матрица податливости.
®
Решение получено методом конечных элементов в программном пакете ABAQUS CAE .
1. Постановка задачи моделирования операции стернотомии 1.1. Общие сведения
Как было сказано выше, стернотомия представляет собой хирургическую операцию, заключающуюся в рассечении грудины. При проведении полной продольной стернотомии разрез производится по пунктирной линии (рис. 1.1, а). Операция выполняется при помощи специального прибора - электростернотома (рис. 1.1, б)
Рисунок 1.1 - а) Грудина; б) разрез с помощью электростернотома
Моделируемый в данной работе робот-ассистент должен провести операцию стернотомии в программном режиме с учетом координатной связи.
® ♦
1.2. Численное моделирование в MATLAB Robotics Toolbox
За простейшую основу моделируемого робота был взят прототип робота-манипулятора Puma 560 (рис. 1.2), рассмотренный в [3].
После добавления третьего связующего звена к приведенному выше рисунку создана следующая кинематическая схема (рис. 1.3, а). Схема представляет собой систему из семи вращательных шарниров и трех связующих звеньев (Ц, Ц2, Ц3). При моделировании численные значения длин связующих звеньев равны: Ц = 0.2 м, Ь2 = 0.6 м, Ь3 = 0.4 м. Согласно примерам, рассмотренным в [3], был смоделирован робот-манипулятор (рис. 1.3, б). Листинг программы приведен в приложении 1.
а)
б)
Рисунок 1.3 - а) Кинематическая схема робота; б) результат моделирования в MATLAB
Robotics Toolbox
®
Далее рассматривается поведение робота при стернотомии. При этом жесткий Г-образный инструмент нагружается продольной силой F и моментом M (рис. 1.4)
Рисунок 1.4 - Упрощенная модель стернотома
Так как толщина кости по длине грудины меняется незначительно ([5, 6]), то усилие, действующее на стернотом, принимается постоянным и равным 2 Н.
Траектория движения задается следующим образом. В ходе операции стернотом проходит расстояние из точки А с координатами (0.2, 0, 0.55) до точки В с координатами
(0.3, 0, 0.55) за 2 секунды. Т.е. инструмент совершает перемещение на расстояние 0.1 м вдоль оси X. После моделирования робота и задания движения можно вывести на экран траекторию движения стернотома в процессе операции (рис. 1.5).
Рисунок 1.5 - Траектория движения инструмента в проекции на плоскость ^х, г)
2. Численное моделирование методом конечных элементов
Для оценки поведения точки приложения сил под действием нагрузки производится расчет перемещения инструмента и определение матрицы податливости для выбранной конфигурации робота.
На рис. 2.1 изображен общий вид рассматриваемого робота-ассистента.
1 - Зажим инструмента с установленным инструментом; 2 - блок приводов ориентации; 3 - предплечье; 4 - плечо; 5 -привод плеча; 6 - опорно-поворотное устройство
Рисунок 2.1 - Общий вид робота
Как было отмечено в пункте 1.2, робот имеет семь степеней свободы и состоит из шести звеньев. Первое звено соответствует установленному инструменту и смоделировано в данной работе упрощенно. Исследуемая конфигурация робота приведена на рис. 2.2.
Рисунок 2.2 - Рассматриваемая конфигурация робота
2.1. Определение поля перемещений
При определении перемещений необходимо учитывать, что к каждому звену приложена нагрузка в виде сосредоточенной силы в центре масс, соответствующей весу звена. Силы тяжести звеньев приведены в таблице 1.
Таблица 1 - Вес звеньев
Номер звена Вес звена, Н
1 10
2 14,5
3 16
4 82
5 56
6 62
При определении перемещений робота под действием только сил тяжести, силы, действующие на инструмент, не рассматривались.
В качестве типа конечного элемента выбран четырехузловой линейный элемент С3Б4. Разбиение на конечные элементы представлено на рис. 2.3.
Рисунок 2.3 - Конечноэлементная модель робота. Количество элементов при данной разбивке равно 489110, число степеней свободы 2742846
Поле перемещений робота под действием сил тяжести всех звеньев представлено на рис. 2.4.
Рисунок 2.4 - Поле перемещений точек робота под действием сил тяжести
Компоненты вектора перемещений точки приложения рабочих сил инструмента по осям глобальной системы координат {X, У, Z} представлены в таблице 2.
Таблица 2 -Перемещения на конце робота для рассматриваемой конфигурации
Направление перемещения Перемещение, м
по оси Х 2.55816-10-5
по оси У -4.48448 • 10-4
по оси Ъ -2.32425 • 10-5
2.2. Расчет матрицы податливости
Определяются компоненты матрицы податливости робота в точке приложения рабочих сил - точке А (рис. 2.3), - для рассматриваемой конфигурации робота. При получении матрицы податливости используется следующее представление:
и = АЖ (1)
где и - вектор перемещений, А - матрица податливости, Ж - вектор обобщенных сил.
3x1 3x3 3x1
Матрица А симметричная и определенно положительная. Матричные вычисления проведены в глобальном декартовом базисе {X, У, Z }.
Вычисление матрицы податливости производится по следующей схеме:
1. К выбранной конфигурации робота прикладываются по очереди три единичные силы в точке А в направлении соответствующих осей (рис. 2.5).
Рисунок 2.5 - Положение системы координат конечной точки инструмента
Точка А - точка приложения силы при расчете матрицы податливости, {х, у, z}-оси локальной системы координат, расположенной в точке А.
2. Вычисляются три компоненты вектора перемещений иф для каждой силы.
3. Заполняется соответствующая колонка матрицы податливости
А = и®
(2)
Правильность вычислений доказывается симметричностью полученной матрицы податливости.
В соответствии с проведенными расчетами, ниже посчитана матрица податливости для рассматриваемой конфигураций робота.
" 2.29 -1.43 2.01"
А = -1.43 11.97 1.61 •10-6
2.01 1.61 8.06
Как видно, матрица податливости является симметричной (в пределах допустимой погрешности вычислений), что доказывает правильность расчета.
Пусть единичный вектор сосредоточенной силы, приложенной к инструменту
равен
0.4959 0.5639 -0.5877
(3)
Тогда согласно уравнению (1), единичный вектор смещения точки приложения силы оказывается равным
еи =•!
-0.14464 0.86484 -0.48079
(4)
Относительное расположение векторов в пространстве представлено на рис. 2.6.
Рисунок 2.6 - Неколлинеарность векторов силы и перемещения
Из рис. 2.6 видно, что в данном случае "несовпадение" направлений силы и перемещения точки приложения силы оказывается существенным. Эту особенность нужно учитывать при выборе рабочей конфигурации.
2.3. Моделирование движения робота в SolidWorks®
Для моделирования операции стернотомии рассмотрена произвольная рабочая конфигурация робота, представленная на рис. 2.7, а. Траектория движения инструмента показана на рис. 2.7, б. В соответствии с параграфом 1.2, длина разреза выбрана равной 0,1 м.
а) б)
Рисунок 2.7 - а) рабочая конфигурация робота; б) траектория движения инструмента
Для определения матрицы податливости в процессе движения робота траектория движения инструмента в оперируемом объекте была разбита на 11 последовательных положений робота, включая начальную и конечную точки разреза. Каждое последующее положение взято на расстоянии 0,01 м от предыдущего. Данные положения представлены на рис. 2.8, где точка А является начальной точкой разреза, а точка В - конечной.
Алгоритм определения матриц податливости для каждого положения описан в параграфе 2.2. Матричные вычисления проведены в локальной системе координат {х, у, z} (см. рис. 2.5).
В соответствии с проведенными расчетами, ниже посчитаны матрицы податливости для 11 указанных положений робота. Индекс указывает на номер положения.
А 0 =
52.635 -4.2272 7.7556 -4.2280 87.489 -4.9174 7.7559 -4.9168 24.222
10-7 А1 =
" 49.788 -4.5026 4.4588 "
А 2 = -4.5033 82.137 -4.8637
4.4591 -4.8632 23.068
" 46.632 -4.7686 1.4897 "
А4 = -4.7692 77.445 -4.7800
1.4899 -4.7797 22.863
" 43.040 -4.9985 -1.0462"
А 6 = -4.9990 72.745 -4.6389
-1.0460 -4.6385 22.900
" 39.157 -5.2366 -3.1972"
А 8 = -5.2369 68.134 -4.4558
-3.1970 -4.4554 23.228
" 35.066
А10 = -5.4611
-4.8285
10-7 А з =
10-7 А5 =
10-
А 7 =
10-
А 9 =
51.229 -4.3706 6.08 "
- 4.37 84.7 - 4.9 •10-
6.08 -4.9 23.5
48.228 -4.6421 2.9420 "
-4.6426 79.703 -4.8357 • 10-7
2.9423 -4.8351 22.822
44.893 -4.8709 0.1601 "
-4.8714 75.005 -4.6991 • 10-7
0.1603 -4.6987 22.862
41.134 -5.1232 -2.1923"
-5.1237 70.418 -4.5533 • 10-7
-2.1921 -4.5529 23.052
37.130 -5.3296 -4.0750"
-5.3300 65.987 -4.3492 • 10-7
-4.0748 -4.3490 23.547
-4.8287"
-4.2394 •Ю-7
23.944
(5)
В пределах допустимой погрешности вычислений матрицы податливости являются симметричными.
Для указанных положений робота построен график отклонения реального перемещения (с учетом координатной связи) от заданного. Отклонение от заданной траектории представлено на рис. 2.9: красная линия и красные точки соответствуют программному положению инструмента (длина пути, равная 0.1 м, отображена в единицу), синие отрезки показывают отклонение конца инструмента от программного в одиннадцати равномерно расположенных узлах (за единицу принято отклонение в конечном положении, равное 4.4 мкм).
N
О ■0.2 0.4 0.6 ■0.8
-1
.2
.4
0
0.5
1.5 0
У
х
Рисунок 2.9 - Отклонение траектории инструмента от заданного движения для 11
положений робота
Заключение
Полученные результаты говорят о том, что предложенная модель робота позволяет выполнить операцию стернотомии. Исходя из данного исследования, можно заключить, что при приложении усилий к инструменту робота смещение точки приложения усилий неколлинеарно вектору приложенной силы на всем пути движения робота. Эту особенность нужно учитывать при выборе рабочей конфигурации.
Благодарность
Работа выполнена при поддержке Министерства образования и науки РФ в рамках федеральной целевой научно-технической программы «Исследования и разработки по приоритетным направлениям развития научно-технологического комплекса России» на 2007 - 2012 годы, номер контракта 16.523.11.3011.
Список литературы
1. Богданова Ю.В., Гуськов А.М., Жидких И.В., Нарайкин О.С., Шурхай В.А.
Проблемы разработки стереотаксических систем наведения хирургического
инструмента на примере биопсии головного мозга // Наука и образование. МГТУ
им. Н.Э. Баумана. Электрон. журн. 2012. № 12. БОТ: 10.7463/0113.0486770
2. Информация о позвоночнике. Режим доступа: http://www.pozvonochnik.info/text/2/content/3/ (дата обращения 20.02.2013).
3. Corke P. Robotics, Vision and Control: Fundamental Algorithms in MATLAB. Springer, 2011. 558 c. DOI: 10.1007/978-3-642-20144-8
4. Банин Е.П., Барышева О.П., Батанов А.Ф., Богданова Ю.В., Гуськов А.М., Козубняк С. А., Нарайкин О.С., Пуценко А.А. Координатная связь медицинского робота-манипулятора // Наука и образование. МГТУ им. Н.Э. Баумана. Электрон. журн. 2012. № 12. DOI: 10.7463/0113.0520630
5. Tillier Y., Paccini A., Durand-Reville M., Chenot J. Finite element modeling for soft tissue surgery based on linear and nonlinear elasticity behavior // Computational Aided Surgery. 2006. Vol. 11, no. 2. P. 63-68.
6. Shields T.W., Ponn R.B., Rusch V.W. General Thoracic Surgery. Lippincott Williams & Wilkins, 2004. 3040 p.
7. Bathe K.J. Finite Element Procedures. Prentice Hall, Englewood Cliffs, New Jersey, 1996. 1037 p.
Приложение 1
Программа численного моделирования условий работы хирургических инструментов с биотканями. Программный комплекс MATLAB Toolbox.
Содержание программы:
% Numerical Simulation of the three-link robotic arm for sternotomy incision
startup rvc;
function from matlab sourse
L1 = .2;
L2 = .6;
L3 = .4;
length of link 1, [m] length of link 2, [m] length of link 3, [m]
tool len = .2;
tool length, [m]
z = .55;
x_startOut = .2; x finalOut = .3;
i; depth of laparotomy incision
i start coordinate for laparotomy incision, [m] end coordinate for laparotomy incision, [m]
% length of the incision = 5 sm
y = 0; % laparotomy incision will be in (XZ) plane
t stepOut = .05; % time step for outside motion, [s]
t_finalOut = 2;
s = 'Rz(q1) Tz(L1) Ry(q2) Rz(q3) Tz(L2) Ry(q4) Rz(q5) Tz(L3) Rz(q6) Ry(q7)';
%10 items
dh = DHFactor(s);
dh.display;
cmd = dh.command('three-link robot arm');
robot = eval(cmd);
q0 = [0 0 0 0 0 0 0];
robot.tool = transl(tool len,0,0);
figure('Name','three-link robot arm');
robot.plot(q0);
robot.teach();
% initial conditions for 7 items % tool creating *rpy2tr(0,pi/2,0); % robot plot
T1 out = transl(x startOut,y,z); % robot trajectory
T2 out = transl(x finalOut,y,z);
q1 = robot.ikine(T1 out);
q2 = robot.ikine(T2 out);
t out = (0:t stepOut:t finalOut);
Ts = ctraj(T1_out,T2_out,length(t_out));
qi = robot.ikine(Ts);
T = robot.fkine(qi);
p = transl(T);
figure('Name','Displacement of the cutting tool in ZX plane'); holdon; grid on; box on;
xlabel('X'); ylabelCZ');
plot(p(:,1), p(:,3),'LineWidth',2); plot(p(1,1), p(1,3),'ko',p(end,1), p(1,3) set(gca,'ylim',[0.54 0.56]);
'*k','LineWidth',2); %to avoid unstability
figure('Name','Joint Coordinates'); holdon; grid on; box on; xlabel('time'); ylabel('q');
plot(t out,qi,'LineWidth',2);
legend!'q1','q2','q3','q4','q5','q6','q7');
x start = x finalOut; % start/final coordinates for the
incision x final = 0.35;
dist = abs(x final-x start); % distance for tool motion
x step = .01; % step along trajectory
step num = dist/x step; % number of steps
% start/final/step time for inside motion
t start = t finalOut;
t_final = t_finalOut*2;
t step = (t final-t start)/step num;
% forces
Fx = 2;
Fy = 0;
Mx = 0;
My = Fx*tool_len; Mz = 0;
T1 = transl(x start,y,z); T2 = transl(x final,y,z); q1 = robot.ikine(T1); q2 = robot.ikine(T2);
t = (t start:t step:t final)'; % time step, [s]
Ts = ctraj(T1,T2,length(t));
qqi = robot.ikine(Ts);
TT = robot.fkine(qqi);
pp = transl(TT);
% assign the load
Q = [Fx,Fy,Fz,Mx,My,Mz];
Torques = robot.jacob0(q2)'*Q';
% conversion qqi (cell-type) to matrix for animation qi = [qi; qqi]; % animation
figure('Name','incision of 5 sm length'); robot.plot(qi);
SCIENTIFIC PERIODICAL OF THE BAUMAN MSTU
SCIENCE and EDUCATION
EL № FS77 - 48211. №0421200025. ISSN 1994-0408
electronic scientific and technical journal
Numerical simulation of robot-assistant operation of sternotomy
with predetermined path
# 05, May 2013
DOI: 10.7463/0513.0599870
Gus'kov A.M., Kuznecova M.S.
Bauman Moscow State Technical University, 105005, Moscow, Russian Federation
gouskov [email protected] [email protected]
This article presents a simulation of kinetostatics of 7-DOF surgical robotic manipulator by the example of sternotomy. Motion of the robot was performed along the predetermined path in the global coordinate system. Numerical simulation was conducted in MATLAB®Robotics Toolbox + Simulink program package. The same motion was investigated in SolidWorks® for calculating the flexibility matrices with a use of the finite-element method. Flexibility matrices for the point of load application were found for 11 successive position of the robot. Compliance matrices were obtained using ABAQUS CAE® software.
Publications with keywords: FEM, surgical robotic manipulator, trajectory tracking, MATLAB® Robotics Toolbox, sternotomy, flexibility matrices, ABAQUS CAE®, , SolidWorks®
Publications with words: FEM, surgical robotic manipulator, trajectory tracking, MATLAB® Robotics Toolbox, sternotomy, flexibility matrices, ABAQUS CAE®, , SolidWorks®
References
1. Bogdanova Yu.V., Gus'kov A.M., Zhidkikh I.V., Naraykin O.S., Shurkhay V.A. Problemy razrabotki stereotaksicheskikh sistem navedeniya khirurgicheskogo instrumenta na primere biopsii golovnogo mozga [Problems of development of stereotaxic guidance systems designed to control surgical instruments, illustrated by brain biopsy]. Nauka i obrazovanie MGTUim. N.E. Baumana [Science and Education of the Bauman MSTU], 2012, no. 12. DOI: 10.7463/0113.0486770
2. Informatsiya opozvonochnike [Information about spine]. Available at: http://www.pozvonochnik.info/text/2/content/3/ , accessed 2013 February 20.
3. Corke P. Robotics, Vision and Control: Fundamental Algorithms in MATLAB. Springer, 2011. 558 p. DOI: 10.1007/978-3-642-20144-8
4. Banin E.P., Barysheva O.P., Batanov A.F., Bogdanova Yu.V., Gus'kov A.M., Kozubnyak
5.A., Naraikin O.S., Pucenko A.A., 2012. Koordinatnaya svyaz' meditsinskogo robota-manipulyatora [Coordinate relationship of the medical robot-manipulator]. Nauka i obrazovanie MGTUim. N.E. Baumana [Science and Education of the Bauman MSTU], 2012, no. 12. DOI: 10.7463/0113.0520630
5. Tillier Y., Paccini A., Durand-Reville M., Chenot J. Finite element modeling for soft tissue surgery based on linear and nonlinear elasticity behavior. Computational Aided Surgery, 2006, vol. 11, no. 2, pp. 63-68.
6. Shields T.W., Ponn R.B., Rusch V.W. General Thoracic Surgery. Lippincott Williams & Wilkins, 2004. 3040 p.
7. Bathe K.J. Finite Element Procedures. Prentice Hall, Englewood Cliffs, New Jersey, 1996. 1037 p.
