ODE45求解微分方程 出现问题 积分公差要求无法满足。

桢心

function dx=spring3(t,x)
M0=13.5e-3;M1=0.3645;M2=1.5;
Pn=2.1e4;F0=5659.5;oumiga=418.67;
R=16.6e-3;R1=5.4e-3;R2=3.4e-3;d=8.6e-3;h=1.5e-4;
f=0.1;k=20;lanmada=100;
miu=8.7e-3;
cp=2.3e3;rou=890;beidae=1e9;alafa=7.8e-4;

%Pn 为柱塞底部液压油压力
%Fs为弹簧力
Fs=5486.5+k*(x(1)+x(3))+lanmada*(x(1)+x(3))^3;

%Fgn为九个柱塞组件的惯性力
%fai为柱塞随时间的变化角度
%fai=oumiga*t;
%Fgn=M0*R*oumiga*0.3249*(1.118*sin(fai+1.107))=0.0039;

%Nn为斜盘对滑靴的支撑力
%静压力n1=pi/2*(R1^2-R2^2)/ln(R1/R2)*Pn
n1=1254.37;
%动压力n2=miu*3.14(R1^2-R2^2)*oumiga*R1/(h-x(1))
%n2=miu*3.14*(R1^2-R2^2)*oumiga*R1/(h-x(1));
n2=1.0870/(h-x(1));
%挤压力n3=3*pi*miu*x(2)/(2*h^3)*((R1^4-R2^4)-((R1^2-R2^2)^2)/ln(R1/R2))
n3=3*3.14*miu*x(2)/(2*h^3)*((R1^4-R2^4)-((R1^2-R2^2)^2)/4.43);
%热锲力n4=16*beida*alafa*miu*oumiga*R2*(2/3-0.78+1/3*0.78^3+0.35^2*R1^3*0.78)/(cp*rou*(h-x(1))^2);
%n4=16*beidae*alafa*miu*oumiga*R1*(2/3-0.78+1/3*0.78^3+0.35^2*R2^3*0.78)/(cp*rou*(h-x(1))^2);
N1=5*(n1+n2+n3);
N2=4*(n2+n3);
Nn=-N1+N2;

%Ff为柱塞与缸体间的摩擦力
%Ff=9*f*(P1+P2)
s1=3.93;
s2=1.02*x(3);
Ff=-1*0.1*0.29*N1*s1-1*0.1*R*oumiga^2*s2*0.95/(h^2);

dx=zeros(4,1);
%dx=[x(2)
    %(1/M1)*(9*3.14*d^2*Pn/4+Fs+Ff-0.95*Nn)
    %x(4)
    %(1/M2)*(F0-Fs)];
dx(1)=x(2);
    dx(2)=(1/M1)*(9*3.14*d^2*Pn/4+Fs+Ff-0.95*Nn);
dx(3)=x(4);
    x(4)= (1/M2)*(F0-Fs);



clc
clear
options=odeset('reltol',0.00001,'abstol',5.344224e-34);
[t,x]=ode23tb(@GuiYi,[0:0.0001:100],[1e-5;0;1e-5;0],options);
x1=x(:,1);
x2=x(:,2);
x3=x(:,3);
x4=x(:,4);
plot(x1,x3);
figure
subplot(1,4,1);
plot(x2,x1);
subplot(1,4,2);
plot(x4,x3);
subplot(1,4,3);
plot(t,x1);
subplot(1,4,4);
plot(t,x3);
% x1=x(:1);
% x2=x(1001:10000,:);
% x3=x(1001:10000,:);
% x4=x(1001:10000,:);
% t=t(1001:10000,:);
% subplot(1,4,1);
% plot(x2,x1);
% subplot(1,4,2);
% plot(x4,x3);
% subplot(1,4,3);
% plot(t,x1);
% subplot(1,4,4);
% plot(t,x3);


哪位大神帮我看看 求助求助
Q863224163