Matlab/Simulink Bouc-Wen模型搭建的磁流变阻尼器

【类型】：磁流变阻尼器
【软件】：Matlab/Simulink
【状态】：模型搭建完成，
【描述】：根据论文搭建的Bouc-Wen模型


基于Bouc-Wen模型的磁流变阻尼器（MR Damper）通常包含力学模型和磁滞非线性两部分。Simulink模型通常由积分器、乘法器和加法器组成，但在代码层面，我们主要关注描述其动态特性的微分方程。

以下是基于Bouc-Wen模型的磁流变阻尼器的MATLAB Function代码（可直接放入Simulink的MATLAB Function Block中）以及参数说明。

核心微分方程

Bouc-Wen模型的标准形式通常如下：

ddot{x} + 2zetaomega_ndot{x} + omega_n^2 x = F_{MR}(t)

其中，磁流变阻尼力 F_{MR} 由以下非线性微分方程描述：

F_{MR} = c_0(dot{z} - dot{x}) + k_0(z - x) + alpha z

dot{z} = Adot{x} - beta|dot{x}| |z|^n text{sgn}(z) - gamma dot{x} |z|^n

MATLAB Function 代码 (Simulink可用)

你可以将以下代码放入 MATLAB Function Block 中，输入为位移 x 和速度 dot{x}，输出为阻尼力 F_{MR}。

function Fmr = fcn_BoucWen(x, dx, i) % x:位移, dx:速度, i:电流
% 说明：Bouc-Wen磁流变阻尼器模型
% 输入: x (位移), dx (速度), i (控制电流)
% 输出: Fmr (阻尼力)

% ==========================================
% 1. 模型参数 (需要根据论文或实验数据标定)
% ==========================================
% 力学参数
c0 = 1500; % 油孔阻尼系数 (Ns/m)
k0 = 1000; % 油孔刚度系数 (N/m)
alpha = 1e5; % 磁滞刚度增益 (N/m)

% Bouc-Wen 磁滞参数
A = 1.0; % 幅值参数
beta = 0.5; % 硬化/软化参数 (正为硬化)
gamma = 0.5; % 磁滞形状参数
n = 2; % 光滑度指数 (通常取2)

% 电流相关增益 (假设力与电流呈线性或非线性关系)
% 有些模型会将 alpha 或 A 乘以电流函数
k_i = 1 + 0.5 * i; % 示例：力随电流增大而增大

% ==========================================
% 2. 积分状态 (需要在函数中保持状态)
% ==========================================
% z 是内部状态变量，需要记忆上一时刻的值
persistent z dz
if isempty(z)
z = 0; % 初始状态
end

% ==========================================
% 3. Bouc-Wen 非线性滞回计算
% ==========================================
% 计算 dz/dt
% dot{z} = Adot{x} - beta|dot{x}| |z|^n text{sgn}(z) - gamma dot{x} |z|^n
dz = A * dx …
beta * abs(dx) * (abs(z) ^ n) * sign(z) …
gamma * dx * (abs(z) ^ n);

% 使用欧拉法积分得到 z (在Simulink中通常用Integrator模块)
% 这里为了代码演示，简单积分，实际建议用Simulink Integrator
% z = z + dz * dt; (此处省略dt，Simulink会自动处理)

% ==========================================
% 4. 计算阻尼力
% ==========================================
% F = c0(dx - dz) + k0(z - x) + alpha * z
Fmr = c0 * (dx - dz) + k0 * (z - x) + alpha * k_i * z;

% ==========================================
% 5. 更新状态 (用于下一次调用)
% ==========================================
% 注意：在Simulink中，通常将 dz 作为输出连接到 Integrator 模块，
% 而不是在这里手动积分。
end

Simulink 模型搭建建议

在Simulink中搭建时，通常不建议在MATLAB Function中做积分，而是使用标准模块：

输入：连接 Position (x) 和 Velocity (dot{x})。
计算 z 的导数 (dot{z})：
使用 Math Function 模块计算 |dot{x}|。
使用 Abs 和 Sign 模块处理非线性项。
使用 Product 和 Sum 模块组合公式 dot{z} = …。
积分：将计算出的 dot{z} 输入到 Integrator 模块，得到状态 z。
计算力 F_{MR}：将 x, dot{x}, z, dot{z} 代入力方程计算最终输出力。

参数标定参考

如果你的论文中有具体的参数表，可以直接替换上述代码中的 c0, k0, alpha, A, beta, gamma。通常，这类模型的拟合效果如下图所示（滞回曲线）：

  F
  ^
  |     /\
  |    /  \
  |   /    \

------|/------____> x

这个模型的逻辑是：
输入一个激励信号（波形发生器）。
信号经过一个微分模块（frac{Delta u}{Delta t}），这通常代表加速度（如果是位移输入）或直接代表控制电压的变化率。
位移信号和处理后的信号（加速度或电压）输入到一个封装的阻尼器模块中。
模块输出阻尼力。

Simulink 模型配置说明

在搭建 Simulink 模型时，你需要以下模块：

输入端：Sine Wave 或 From Workspace。
微分：使用 Derivative 模块（注意：Simulink 的 Derivative 模块在仿真中可能不稳定，建议用 Transfer Fcn s/(s+1e-3) 代替）。
核心计算：使用 MATLAB Function 模块（代码如下）。
积分器：为了得到“位移”，你需要在输入信号后串联一个 Integrator 模块。
微分器：为了得到“速度”，你需要将位移信号通过一个 Derivative 模块。

MATLAB Function 代码 (Bouc-Wen 模型)

这段代码实现了经典的 Bouc-Wen 模型，输入为位移和电流，输出为阻尼力。它会自动计算内部状态变量 z。

function Fmr = fcn(x, dx, i)
% fcn - 磁流变阻尼器 Bouc-Wen 模型
% 输入:
% x : 位移 (m)
% dx : 速度 (m/s)
% i : 控制电流 (A)
% 输出:
% Fmr: 阻尼力 (N)

% ==========================================
% 1. 参数初始化 (请根据你的论文或实验数据修改)
% ==========================================
% 力学参数 (线性部分)
c0 = 1500; % 油孔阻尼系数 (Ns/m)
k0 = 1000; % 油孔刚度系数 (N/m)

% Bouc-Wen 磁滞参数
A = 1.0; % 幅值参数 (无量纲)
beta = 0.5; % 磁滞形状参数
gamma = 0.5; % 磁滞形状参数
n = 2; % 光滑度参数 (通常取 1 或 2)
alpha = 1e5; % 磁滞刚度增益 (N/m)

% 电流相关增益 (假设力与电流线性相关，也可用非线性函数)
% 这里假设电流从 0 到 2A，力从 0 到 2000N
k_i = 1000 * i;

% ==========================================
% 2. 状态变量 z 的微分方程 (Bouc-Wen 核心)
% dz/dt = dx - beta|dx|z|^nsgn(z) - gammdx|z|^n
% ==========================================
% 注意：这里使用 persistent 变量来存储上一时刻的 z 值
% 在 Simulink 中，更推荐使用 “State” 定义或在 Simulink 层面使用 Integrator

% 由于 MATLAB Function 无法直接积分，这里演示计算 dzdt
% 实际使用时，你需要在 Simulink 中用 Integrator 积分 dzdt 得到 z
% 这里为了完整性，假设我们已经有了 z (在复杂模型中，z 是状态变量)

% 假设 z 是上一时刻的值（简化演示，实际需用 Integrator）
% 在真实模型中，dzdt 应该连接到 Integrator，输出 z 再反馈回来
persistent z_val
if isempty(z_val)
z_val = 0; % 初始状态
end

% 计算 dzdt
abs_dx = abs(dx);
abs_z = abs(z_val);
sgn_z = sign(z_val);

dzdt = A * dx - beta * abs_dx * (abs_z^n) * sgn_z - gamma * dx * (abs_z^n);

% 更新 z (欧拉法近似，实际用 Simulink Integrator 更好)
% z_val = z_val + dzdt * dt; % (Simulink 会自动处理积分)

% 这里我们假设 z_val 已经由外部积分器提供，或者我们简单更新
% 为了演示，我们直接计算 z (仅用于测试)
% 实际部署时，请移除此行，改为输入端口输入 z
z_val = z_val + dzdt * 0.001;

% ==========================================
% 3. 计算阻尼力 Fmr
% ==========================================
% Fmr = cdx + k0x + alpha*z + k_i * (非线性项)
% 简化版本：考虑电流对刚度的影响
Fmr = c0 * dx + k0 * x + alpha * z_val + k_i * z_val;

end

搭建建议

积分器处理：Bouc-Wen 模型中的 z 是一个状态变量（dot{z} = f(x, dot{x})）。在 Simulink 中，不要在 MATLAB Function 里用 persistent 变量做积分，那样精度低且不稳定。
正确做法：在 Simulink 层面，使用一个 Integrator 模块。输入是上面代码中的 dzdt，输出是 z。然后将 z 反馈回 MATLAB Function 或用于力的计算。
微分模块：图中的 frac{Delta u}{Delta t} 如果是位移信号，则代表加速度。在 Simulink 中，直接用 Derivative 模块即可，但要注意噪声问题。
输入信号：图中的波形发生器通常是正弦波（Sine Wave），用于模拟路面激励或振动台输入。

完整的 Simulink 连接逻辑

Sine Wave (幅值 0.01m, 频率 1Hz) -> Integrator (得到位移 x)。
Integrator 输出 -> Derivative (得到速度 dx)。
Constant (电流 0~2A) -> 输入到阻尼器模块。
位移 x 和 速度 dx 以及 电流 -> 输入到 MATLAB Function (即上面的代码)。
MATLAB Function 输出 力 Fmr -> Scope。