🏆本文收录于《CSDN问答解惑》专栏，主要记录项目实战过程中的Bug之前因后果及提供真实有效的解决方案，希望能够助你一臂之力，帮你早日登顶实现财富自由🚀；同时，欢迎大家关注&&收藏&&订阅！持续更新中，up！up！up！！

问题描述

    matlab数值溢出该怎么解决?核心是解一个三自由度微分方程组,解出θ关于时间t的变化图像.设定一个变转速的工况,图中g(t)是功率随时间t变化函数,wm是根据功率求得的角速度.图一是主程序,图二是参数和方程,图三是求解函数,系统提示我数值溢出了,想问一下我这个错在哪里有没有修改办法,万分感谢

clear;clc;close all;
%初始参数，该参数变化，myfun程序要对应调整
syms t;
g(t)=(1460*9.8+1.04*1460*2.67+0.32*2.25*(10+2.67*t)^2/21.15)*(10+2.67*t)/3600;
wm=(g(t)*9550/((1460*9.8+1.04*1460*2.67+0.32*2.25*(10+2.67*t)^2/21.15)*0.4064/(0.93*14)))*pi*2/60;zi=4;%内转子极对数
zo=14;%外转子极对数
z1=32;%小齿轮齿数
z2=99;%大齿轮齿数
ms=4.62e-3;
% r2=ms*z1/2;%小齿轮分度圆
% r3=ms*z2/2;%大齿轮分度圆
I=z2/z1;%机械齿轮速比
% G=zo/zi;
% ms1=4.12e-3;
G=(zo+zi)/zi;%磁齿轮速比
%一、计算过程
tmax=20;
% tmax=900*2*pi/wm;% 计算时长，按输入轴的周期数计算
bc=2e-5;
y0=[0;0;0;wm;wm/G;wm/G/I];
[t,y]=odeRK4sys('myfun3',tmax,bc,y0);
figure;plot(t,y(:,1)/G-y(:,2));
title('时域响应图');
xlabel('时间 (s)');
ylabel('扭转振幅(rad)');%磁齿轮的相对弹性角度
figure;plot(t,y(:,2)/I-y(:,3));
title('时域响应图');
xlabel('时间 (s)');
ylabel('扭转振幅(rad)');%机械齿轮的相对弹性转角;
figure;plot(t,y(:,4),t,y(:,5),t,y(:,6));
title('时域响应图');
xlabel('时间 (s)');
ylabel('转速(rad/s)');
function dydt=myfun3(t,y)
syms t;
g(t)=(1460*9.8+1.04*1460*2.67+0.32*2.25*(10+2.67*t)^2/21.15)*(10+2.67*t)/3600;
wm=(g(t)*9550/((1460*9.8+1.04*1460*2.67+0.32*2.25*(10+2.67*t)^2/21.15)*0.4064/(0.93*14)))*pi*2/60;
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+1.04*1460*2.67+0.32*2.25*(10+2.67*t)^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=116;%整车等效转动惯量%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>0fx=xx-b;elseif xx+b<0fx=xx+b;elsefx=0;end
F23=km*fx+cm*dxx;
T0=1.7*Tm*G;
Tc1=0.1*T0*sin(wm*t);
Tc=0.1*(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
function [t,y] = odeRK4sys(diffeq,tn,h,y0)
t = (0:h:tn)';         %  Column vector of elements with spacing h
nt = length(t);        %  number of steps
neq = length(y0);      %  number of equations that are simultaneously advanced
y = zeros(nt,neq);     %  Preallocate y for speed
y(1,:) = y0(1,:)';       %  Assign IC. y0(:) is column, y0(:)' is row vectorh2 = h/2;  h3 = h/3;  h6 = h/6;   %  Avoid repeated evaluation of constants    
k1 = zeros(neq,1);  k2 = k1;      %  Preallocate memory for the Runge-Kutta
k3 = k1;  k4 = k1;  ytemp = k1;   %  coefficients and a temporary vector%  Outer loop for all steps:  j = time step index;  n = equation number index 
for j=2:nttold = t(j-1);   yold = y(j-1,:)';   %  Temp variables, yold is column vectork1 = feval(diffeq,told,yold);       %  Slopes at the starting pointfor n=1:neqytemp(n) = yold(n) + h2*k1(n);     %  Estimate all y's at midpointendk2 = feval(diffeq,told+h2,ytemp);    %  1st estimate of slopes at midpointfor n=1:neqytemp(n) = yold(n) + h2*k2(n);     %  2nd estimate of all y's at midpointendk3 = feval(diffeq,told+h2,ytemp);    %  2nd estimate of slopes at midpointfor n=1:neqytemp(n) = yold(n) + h*k3(n);      %  Estimate y at end pointendk4 = feval(diffeq,told+h,ytemp);     %  Estimate of slopes at endpointfor n=1:neq                          %  Simultaneously advance all y'sy(j,n) = yold(n) + h6*(k1(n)+k4(n)) + h3*(k2(n)+k3(n));end

报错截图如下：

如上问题有来自我自身项目开发，有的收集网站，有的来自读者，如有侵权，立马删除。

解决方案

  如下是上述问题的解决方案，仅供参考：

  在MATLAB中，数值溢出通常发生在计算过程中的数值超出了MATLAB能够表示的数值范围。根据您提供的代码，以下是一些可能导致数值溢出的原因和相应的解决方案：

1. 过大的数值：在您的代码中，g(t) 和 wm 计算中包含了大数值和复杂的表达式，这可能导致在某些 t 的值下计算结果过大。

   解决方案：检查 g(t) 和 wm 的表达式，尝试简化或分解表达式，或者在计算过程中使用更高精度的数据类型。

2. 变量初始化问题：在 odeRK4sys 函数中，t 是通过 (0:h:tn)'; 初始化的，如果 tn 或 h 过大，可能导致 t 数组中的数值超出范围。

   解决方案：确保 tn 和 h 的值是合理的，并且不会在循环中产生过大的数值。

3. 循环中的数值累积：在 odeRK4sys 函数的循环中，如果每一步的增量过大，可能导致 y(j,:) 的值在迭代过程中迅速增长并溢出。

   解决方案：检查 feval 调用的结果，确保每一步的增量不会过大。可能需要调整步长 h 或者使用更稳健的数值积分方法。

4. 符号计算：您的代码中使用了符号计算（syms t;），这在循环中可能导致性能问题和数值不稳定。

   解决方案：如果可能，尝试将符号表达式转换为数值表达式，使用 double 或 eval 函数来避免在每次迭代中重复符号计算。

5. 检查 myfun3 函数：myfun3 函数中可能存在某些计算导致数值溢出，特别是在计算 km，F23，T21 等表达式时。

   解决方案：检查这些计算，确保它们在所有 t 的值下都能产生合理的数值。

6. 使用条件语句：在 myfun3 函数中，您使用了 if-elseif-else 语句来确定 fx 的值。如果 xx 的值在 b 附近快速变化，可能导致 fx 的值不稳定。

   解决方案：考虑使用更平滑的过渡函数来代替这个条件语句，例如 smoothstep 函数。

7. 调试和测试：逐步测试代码的不同部分，以确定导致溢出的确切位置。

8. 使用MATLAB的警告和错误信息：MATLAB在发生数值溢出时会提供警告和错误信息，仔细阅读这些信息可以帮助您定位问题。

  最后，如果问题依然存在，您可能需要逐步调试代码，检查每个变量在每个时间步的值，以确定何时何地发生了数值溢出。使用MATLAB的调试工具，如 dbstop 命令，可以帮助您在特定条件下暂停执行并检查变量的值。

  希望如上措施及解决方案能够帮到有需要的你。

  PS：如若遇到采纳如下方案还是未解决的同学，希望不要抱怨&&急躁，毕竟影响因素众多，我写出来也是希望能够尽最大努力帮助到同类似问题的小伙伴，即把你未解决或者产生新Bug黏贴在评论区，我们大家一起来努力，一起帮你看看，可以不咯。

  若有对当前Bug有与如下提供的方法不一致，有个不情之请，希望你能把你的新思路或新方法分享到评论区，一起学习，目的就是帮助更多所需要的同学，正所谓「赠人玫瑰，手留余香」。

☀️写在最后

  ok，以上就是我这期的Bug修复内容啦，如果还想查找更多解决方案，你可以看看我专门收集Bug及提供解决方案的专栏《CSDN问答解惑》，都是实战中碰到的Bug，希望对你有所帮助。到此，咱们下期拜拜。

码字不易，如果这篇文章对你有所帮助，帮忙给bug菌来个一键三连(关注、点赞、收藏) ，您的支持就是我坚持写作分享知识点传播技术的最大动力。

同时也推荐大家关注我的硬核公众号:「猿圈奇妙屋」 ；以第一手学习bug菌的首发干货，不仅能学习更多技术硬货，还可白嫖最新BAT大厂面试真题、4000G Pdf技术书籍、万份简历/PPT模板、技术文章Markdown文档等海量资料，你想要的我都有！

📣关于我

我是bug菌，CSDN | 掘金 | InfoQ | 51CTO | 华为云 | 阿里云 | 腾讯云 等社区博客专家，C站博客之星Top30，华为云2023年度十佳博主，掘金多年度人气作者Top40，掘金等各大社区签约作者，51CTO年度博主Top12，掘金/InfoQ/51CTO等社区优质创作者；全网粉丝合计 30w+；硬核微信公众号「猿圈奇妙屋」，欢迎你的加入！免费白嫖最新BAT互联网公司面试真题、4000G PDF电子书籍、简历模板等海量资料，你想要的我都有，关键是你不来拿哇。

---