旋转双心室辅助装置的无传感器生理控制、防吸和流量平衡算法（Matlab代码实现）

最新推荐文章于 2025-04-09 10:30:59 发布

yyds_2201

最新推荐文章于 2025-04-09 10:30:59 发布

阅读量611

点赞数 8

文章标签：算法 matlab 数据结构

本文链接：https://blog.csdn.net/2201_75535220/article/details/143769352

版权

💥💥💞💞欢迎来到本博客❤️❤️💥💥

🏆博主优势：🌞🌞🌞博客内容尽量做到思维缜密，逻辑清晰，为了方便读者。

⛳️座右铭：行百里者，半于九十。

📋📋📋本文目录如下：🎁🎁🎁

目录

⛳️赠与读者

💥1 概述

📚2 运行结果

🎉3 参考文献

🌈4 Matlab代码实现

⛳️赠与读者

👨‍💻做科研，涉及到一个深在的思想系统，需要科研者逻辑缜密，踏实认真，但是不能只是努力，很多时候借力比努力更重要，然后还要有仰望星空的创新点和启发点。当哲学课上老师问你什么是科学，什么是电的时候，不要觉得这些问题搞笑。哲学是科学之母，哲学就是追究终极问题，寻找那些不言自明只有小孩子会问的但是你却回答不出来的问题。建议读者按目录次序逐一浏览，免得骤然跌入幽暗的迷宫找不到来时的路，它不足为你揭示全部问题的答案，但若能让人胸中升起一朵朵疑云，也未尝不会酿成晚霞斑斓的别一番景致，万一它居然给你带来了一场精神世界的苦雨，那就借机洗刷一下原来存放在那儿的“躺平”上的尘埃吧。

或许，雨过云收，神驰的天地更清朗.......🔎🔎🔎

💥1 概述

旋转双心室辅助装置的无传感器生理控制、防吸和流量平衡算法研究

目的：
旋转双心室辅助装置（BiVAD）是一种专为双心室衰竭患者设计的机械泵，通过植入左心室和右心室来泵送血液并提供机械循环支持。本研究旨在开发和测试一种创新的无传感器控制算法，该算法旨在实现以下目标：提供生理控制（即BiVAD流量满足心脏需求）、防止心室抽吸以及提供平衡的左右（全身和肺部）流量，同时避免在非线性、时变和不连续的循环系统中使用植入式流量或压力传感器。

方法：
控制算法的核心是两个增益调度比例积分控制器，分别用于左心室和右心室辅助装置。该算法仅依赖于固有的泵参数（如速度和功率）来维持差动泵速度（ΔRPM_L和ΔRPM_R）高于用户定义的阈值，以防止心室抽吸。同时，通过计算平均参考压头（ΔP_L和ΔP_R）来确保生理灌注和左右侧流量的平衡。为估算ΔP_L和ΔP_R，本研究采用了基于模型的方法，并结合了扩展卡尔曼滤波器和格雷-萨维茨基滤波器。在模拟休息和运动测试条件下，通过计算机仿真评估了该算法的有效性和鲁棒性，测试场景包括：1）过大的ΔP_L和/或ΔP_R设定点；2）肺血管或腔静脉阻力迅速增加三倍；3）从运动到休息的过渡；以及4）心室颤动。

结果与结论：
所提出的无传感器BiVAD算法在所有测试条件下均表现出色，成功防止了心室抽吸，恢复了生理灌注，并保持了左右两侧的流量平衡。这一算法不仅满足了生理控制的需求，还显著提高了系统的安全性和可靠性。

意义：
由于该算法无需对设备进行任何修改，因此可以轻松地集成到当前的临床BiVAD系统中。这一创新不仅为双心室衰竭患者提供了更有效的治疗方案，还为未来的心脏辅助装置研发开辟了新的道路。

目的：旋转双心室辅助装置（BiVAD）是一种机械泵，植入双心室衰竭患者的左心室和右心室，用于泵血并提供机械循环支持。本文的目的是开发和测试一种新颖的无传感器控制算法，该算法同时满足提供生理控制（BiVAD流量满足心脏需求）、防止心室抽吸和提供平衡的左右（全身和肺部）流量的目标，而无需在非线性、时变和不连续的循环系统中使用植入式流量或压力传感器。方法：控制算法由用于左心室和右心室辅助装置的两个增益调度比例积分控制器组成，只需要固有的泵参数（速度和功率）来维持差动泵速度（ Δ RPM L 和 Δ RPM R ）高于用户定义的阈值，以防止心室抽吸，以及平均参考压头（ ΔPL ， ΔPR ），以提供生理灌注和平衡左右侧流量。采用基于模型的方法，结合扩展卡尔曼滤波器和格雷-萨维茨基滤波器来估算ΔP L和ΔP R。在模拟休息和运动测试条件下，在计算机上评估了该算法的有效性和鲁棒性，用于：1）过大的ΔPL和/或ΔPR设定点；2）肺血管或腔静脉阻力迅速增加三倍；3)从运动到休息的过渡；以及4)心室颤动。结果与结论：所提出的无传感器BiVAD算法成功防止了抽吸，恢复了生理灌注，并在所有测试条件下固有地保持了左右两侧的平衡。意义：所提出的算法不需要任何设备修改，可以集成到当前的临床BiVAD中。

📚2 运行结果

全部运行结果：
链接: 百度网盘请输入提取码

提取码: b57t
--来自百度网盘超级会员v6的分享

部分代码：

C2=1/.290; % Compliance of Pulmonary Artery
C6=1/1.03; % Compliance of Aorta
v1=144.9; % Initial condition of Right Heart Volume (RHV)
v2=50.8; % Initial condition of pulmonary artery volume
v3=297.0; % Initial condition of pulmonary arterial volume
v4=408.2; % Initial condition of pulmonary venous volume
v5=146.1; % Initial condition of Left Heart Volume (LHV)
v6=75.5; % Initial condition of Aorta volume
v7=3221.7; % Initial condition of systemic circulation volume
v8=504.1; % Initial condition of Vena Cava volume
v9=22.5; % Initial condition of coronary volume
v10=2.5; % Initial condition of subclavian volume
v15=v2./C2; % Initial condition of Pulmonary Artery Pressure (PAP)
v16=v6./C6; % Initial condition of Aoric Pressure (AoP)
y0 = [v1 v2 v3 v4 v5 v6 v7 v8 v9 v10 0 0 0 0 v15 v16]'; % Initial conditions

%********** pump parameters **********%
Kb=.003;
j=.916e-6;
b=.66e-6;
a0=.738e-12;
a1=.198e-10;
b0=-.296;
b1=-.027;
b2=.0000933;
B=b2/b1;
E=-1/b1;
F=3*Kb/(2*j);
J=-b/j;
K=-a0/j;
L=-a1/j;
A1=b0/b1;
%********** pump parameters **********%

dpR=20; % Right deltaP setpoint
dpsR=665; % Right deltaRPM setpoint

dpL=110; % Left deltaP setpoint
dpsL=1535; % Left deltaRPM setpoint

Fs=100;
Fc=5;
[BB,AA]=butter(2,Fc/(Fs/2)); % lowpass Butterworth filter

Ti=5; % coefficient in the control law
cofac=3; % coefficient for
Rcorin=100; % coronary inlet resistance

c11=1/.016; % initial/final value of compliance of right heart
c22=1/.335; % initial/final value of compliance of left heart
C3=1/.049; % Compliance of Pulmonary arterial
C8=1/.0071; % Compliance of Vena Cava
C4=1/.033; % Compliance of Pulmonary venous

bt_time = 0.75; % beat time (cardiac cycle)
Rfac = 1.0; % resistance tuning
cfacr = 0.953; % compliance tuning for right heart
cfacL = 5.293; % compliance tuning for left heart

%resetting the seed for the pseudo random number generator to compare performance bet. 2 runs.
rand('state',0);
% gg=input('number of beats: '); % number of heart beats
g = 80; % number of seconds
dt=.01; % step size
n=g/dt; % number of points
t0=0;
tf=t0+dt;
y(1,:)=y0';
P1(1)= v1./c11/cfacr; % initial right heart pressure
P2(1)= v2./C2; % initial pulmonary artery pressure
P5(1)= v5./c22/cfacL; % initial left heart pressure
P6(1)= v6./C6; % initial aortic pressure

measR(1,1)=y(1,12); % initial Right rpm
measR(1,2)=y(1,15)-P1(1,1); % initial Right Delta P
delppR(1)=measR(1,2);
delpeR(1)=delppR(1);
delpR=delppR(1);
ye1R(1,1)=0;
ye1R(1,2)=0;
IRnew=0;
mesR(1,1)=measR(1,1)+2*sqrt(.0001)*(rand-.5)*measR(1,1);% initial measured Right rpm

measL(1,1)=y(1,14); % initial Left rpm
measL(1,2)=y(1,16)-P5(1,1); % initial Left Delta P
delppL(1)=measL(1,2);
delpeL(1)=delppL(1);
delpL=delppL(1);
ye1L(1,1)=0;
ye1L(1,2)=0;
ILnew=0;
mesL(1,1)=measL(1,1)+2*sqrt(.0001)*(rand-.5)*measL(1,1);% initial measured Left rpm

tic
for count=1:n
[t,yn]=ode15s('model',[t0 (t0+tf)/2 tf],y0,[],ILnew,IRnew,bt_time,cofac,Rcorin,C2,C3,C4,C6,C8,Rfac,cfacr,cfacL,c11,c22); % solver
tdh0R=ye1R(count,:)';
tdh0L=ye1L(count,:)';
% non-linear prediction
[tomR,dickR]=ode15s('wfp',[t0 (t0+tf)/2 tf],tdh0R,[],IRnew,delpR);
[tomL,dickL]=ode15s('wfp',[t0 (t0+tf)/2 tf],tdh0L,[],ILnew,delpL);

% getting compliances
tt=(count)*dt;

c=.065; %compliace tuning for subclavian

Crvp1=[c11 1/.03]*cfacr; % part 1
Crvp2=[1/.03 1/.046 1/.070 1/.089 1/.106 1/.120 1/.124 1/.116 1/.05]*cfacr; % part 2
Crvp3=[1/.05 c11]*cfacr; % part 3

Clvp1=[c22 1/.75]*cfacL; % part 1
Clvp2=[1/.75 1/1.18 1/1.59 1/1.92 1/2.17 1/2.32 1/2.35 1/1.99 1/1.1]*cfacL; % part 2;
Clvp3=[1/1.1 c22]*cfacL; % part 3

% t1,t2,t3= time points for curve fitting for the 3 parts
t1=bt_time*[0 0.04]; % the first part is valid only for the first 0.04 seconds of the beat
t2=bt_time*[0.04 0.08 0.12 0.16 0.2 0.24 0.28 0.32 0.36]; % the second from .04 to .36secs
t3=bt_time*[0.36 0.4]; % and the third from .36 to .4 seconds.
warning off
%curve fitting for the 3 parts
[PCrv1,S1]=polyfit(t1,Crvp1,2);
[PClv1,S2]=polyfit(t1,Clvp1,2);
[PCrv2,S1]=polyfit(t2,Crvp2,20);
[PClv2,S2]=polyfit(t2,Clvp2,20);
[PCrv3,S1]=polyfit(t3,Crvp3,2);
[PClv3,S2]=polyfit(t3,Clvp3,2);

%interpolation between points for non-linear compliances for active blocks
if rem(tt,1*bt_time)<=.04*bt_time
Crv=polyval(PCrv1,rem(tt,1*bt_time));
Clv=polyval(PClv1,rem(tt,1*bt_time));
Rcorin=100;
C11=1/(1.4*cofac);
end
if and(rem(tt,1*bt_time)>.04*bt_time,rem(tt,1*bt_time)<=.36*bt_time)
Crv=polyval(PCrv2,rem(tt,1*bt_time));
Clv=polyval(PClv2,rem(tt,1*bt_time));
Rcorin=100;
C11=1/(.8*cofac);
end
if and(rem(tt,1*bt_time)>.36*bt_time,rem(tt,1*bt_time)<=.4*bt_time)
Crv=polyval(PCrv3,rem(tt,1*bt_time));
Clv=polyval(PClv3,rem(tt,1*bt_time));
Rcorin=100;
C11=1/(.8*cofac);
end
% If from 0.4 secs to 1 second the compliances are constant for the beat.
if rem(tt,1*bt_time)>.4*bt_time
Crv=c11*cfacr;
Clv=c22*cfacL;
Rcorin=37;
C11=1/(.7*cofac);
end

C10=c*1/2; % subclavian compliance
Rscin=.1; % sibclavian inlet resistance

t0=tf;
tf=t0+dt;
y(count+1,:)=yn(end,:);
y0=yn(end,:);
dickR=dickR(end,:);
dickL=dickL(end,:);
P1(count+1,1)=y(count+1,1)/Crv; % RH Pres.
P2(count+1,1)=y(count+1,2)/C2; % PA Pres.
P5(count+1,1)=y(count+1,5)/Clv; % LH Pres.
P6(count+1,1)=y(count+1,6)/C6; % AoP Pres.

if (-P6(count+1,1)+P5(count+1,1)) > 0
R6i(count+1,1)=Rfac*80; % aorta resistance for aortic valve is open
else
R6i(count+1,1)=0; % aorta resistance for aortic valve is closed
end
AoF(count+1,1)=(-P6(count+1,1)+P5(count+1,1))*R6i(count+1,1)*60/1000; % aortic flow

if (-P2(count+1,1)+P1(count+1,1)) > 0
R2i(count+1,1)=Rfac*220; % pulmonary artery resistance for pulmonary valve is open
else
R2i(count+1,1)=0; % pulmonary artery resistance pulmonary valve is closed
end
PAF(count+1,1)=(-P2(count+1,1)+P1(count+1,1))*R2i(count+1,1)*60/1000; % pulmonary artery flow