matlab如何绘制三维瀑布图

Question

如图四图五的程序所示,求解微分方程并得出时域响应图,见图一和图二,现在想要改变参数(齿侧间隙)b的值,b的范围为1e-5到2e-4,每隔1e-5增加一次,得到总共二十个值并作为y轴,x轴为时间t,z轴为振幅响应(y_1(:,2)*rb2-y_1(:,3)*rb3-b),目标是得到类似图三的瀑布图,请问有什么建议吗?非常感谢您的回答

```c++
clear;clc;close all;
%初始参数，该参数变化，myfun程序要对应调整
syms t;
% 起步加速
wm_func = @(t)((9+2.67*t*3.6)*9550/(3600*0.4064/(0.93*14)))*pi*2/60;
zi=4;%内转子极对数
zo=14;%外转子极对数k
z1=32;%小齿轮齿数
z2=99;%大齿轮齿数
ms=4.62e-3;
r2=ms*z1/2;%小齿轮分度圆
r3=ms*z2/2;%大齿轮分度圆
arf0=20/180*pi;
rb2=r2*cos(arf0);%小齿轮基圆
rb3=r3*cos(arf0);%大齿轮基圆
I=z2/z1;%机械齿轮速比
G=(zo+zi)/zi;%磁齿轮速比
b=1e-4;
%b=1e-5:1e-5:2e-4;
%一、计算过程
tmax=3.76;
bc=1e-4;
y0=[1400;1400/G;1400/(G*I);wm_func(1);wm_func(1)/G;wm_func(1)/(G*I)];
[t,y]=odeRK4sys('myfun3',tmax,bc,y0);
start_index = find(t >= 0, 1, 'first');
t_1 = t(start_index:end);
y_1 = y(start_index:end, :);

figure;plot(t_1,y_1(:,2)/I-y_1(:,3));
signal1 = y_1(:,2)/I - y_1(:,3);
set(gca,'FontSize',16);
title('时域响应图','FontSize',16);
xlabel('时间 (s)','FontSize',16);
ylabel('扭转振幅(rad)','FontSize',16);%磁齿轮的相对弹性角度
figure;plot(t_1,y_1(:,2)*rb2-y_1(:,3)*rb3-b);
signal2=y_1(:,2)*rb2-y_1(:,3)*rb3-b;
set(gca,'FontSize',16);
title('时域响应图','FontSize',16);
xlabel('时间 (s)','FontSize',16);
ylabel('线位移振幅(m)','FontSize',16);%机械齿轮的相对弹性转角;

function dydt=myfun3(t,y,wm_func)
wm_func = @(t)((9+2.67*t*3.6)*9550/(3600*0.4064/(0.93*14)))*pi*2/60;
wm=wm_func(t);
zi=4;%内转子极对数
zo=14;%外转子极对数
z1=32;%小齿轮齿数
z2=99;%大齿轮齿数
arf0=20/180*pi;%压力角
ms=4.62e-3;
r2=ms*z1/2;%小齿轮分度圆
r3=ms*z2/2;%大齿轮分度圆
rb2=r2*cos(arf0);%小齿轮基圆
rb3=r3*cos(arf0);%大齿轮基圆
I=z2/z1;%机械齿轮速比
G=(zo+zi)/zi;%磁齿轮速比
Tm =(1460*9.8*0.012+1.04*1460*2.67+0.32*2.25*(9+2.67*t*3.6)^2/21.15)*0.4064/(0.93*14);%起步加速
TL=Tm*I*G;%变速器负载
%3 转动惯量
IM=0.081197;%电机转子转动惯量kgm2
I0=0.181197;%磁齿轮低速级转动惯量kgm2
I1=0.013663;%磁齿轮高速级转动惯量kgm2
I2=0.0230596;%机械小齿轮转动惯量
I3=1.9825692;%机械大齿轮转动惯量
IL=21;%整车等效转动惯量

%4 刚度与阻尼
%齿轮的阻尼和刚度
w1=wm*G;%机械小齿轮输入转速
kp=10.51e8;
kb=0.25e8;%两个齿轮的刚度曲线为正弦曲线，
km=kp+kb*sin(w1*t*z1);
sg=0.16;%齿轮阻尼比
cm=2*sg*sqrt(kp*I2*I3/(I2*rb3^2+I3*rb2^2));
 b=1e-4;
 xx=rb2*y(2)-rb3*y(3);
 dxx=rb2*y(5)-rb3*y(6);

  if xx-b>0
      fx=xx-b;
  elseif xx+b<0
      fx=xx+b;
  else
      fx=0;
  end
F23=km*fx+cm*dxx;
T0=1.7*Tm*G;
Tc1=0.2*T0*sin(wm*t);
Tc=0.03*(zi*y(4)-(zi+zo)*y(5));
T21=T0*sin(zi*y(1)-(zi+zo)*y(2))+Tc+Tc1;
y(7)=(Tm-T21/G)/(I0+IM);
y(8)=(T21-r2*F23)/(I2+I1);
y(9)=(-TL+r3*F23)/(I3+IL);
dydt=[y(4)
      y(5)
      y(6)
      y(7)
      y(8)
      y(9)];
end

Answer 1

以下内容由AIGC及阿里嘎多学长共同生成、有用望采纳：

---

要使用 MATLAB 绘制三维瀑布图，你可以参考以下步骤和代码示例：

1. 准备数据：首先，你需要为不同的参数 b 值准备数据。这包括计算每个 b 值对应的时域响应。

2. 使用循环：通过循环遍历 b 的所有值，并为每个值计算响应。

3. 计算响应：在每次循环中，调用你的微分方程求解器 odeRK4sys 并计算响应。

4. 存储数据：将每个 b 值的时域响应存储在矩阵中。

5. 绘制三维图：使用 surf 或 mesh 函数来绘制三维图。

6. 自定义视图：可以使用 view 函数来调整视角，使瀑布图效果更加明显。

以下是根据你的要求修改后的 MATLAB 代码示例：

clear; clc; close all;

% 初始参数设置
% ...（省略了你提供的初始参数设置代码）

% 参数 b 的范围
b_values = 1e-5:1e-5:2e-4;
num_bs = length(b_values);

% 初始化存储结果的数组
t_values = zeros(1,1000); % 假设每个 b 值的时间向量长度为 1000
y_values = zeros(num_bs, 1000); % 存储所有 b 值的响应

% 循环遍历 b 的值
for idx = 1:num_bs
    b = b_values(idx);
    
    % 计算过程
    % ...（省略了你提供的计算过程代码）
    
    % 计算振幅响应
    signal = y_1(:,2)*rb2 - y_1(:,3)*rb3 - b;
    
    % 存储结果
    t_values = t_1; % 更新时间向量
    y_values(idx, :) = signal; % 存储当前 b 值的响应
end

% 绘制三维瀑布图
figure;
[X, Y] = meshgrid(t_values, b_values);
surf(X, Y, y_values);
xlabel('时间 (s)');
ylabel('参数 b');
zlabel('振幅响应');
title('三维瀑布图');

% 自定义视图以获得瀑布图效果
view(0, 90); % 视角从上方看下来
colorbar; % 显示颜色条

请注意，这段代码假设你已经有了 odeRK4sys 函数和 myfun3 函数的实现。此外，你可能需要根据实际的数据调整数组的大小和其他参数。这段代码将生成一个三维图，显示不同参数 b 值下的振幅响应随时间的变化。