t=[0:0.04:1.48];
y=[7.86 7.84 7.82 7.77 7.72 7.68 7.61 7.51 7.42 7.33 7.21 7.07 6.94 6.79 6.64 6.48 6.29 6.11 5.92 5.72 5.50 5.27 5.03 4.78 4.53 4.25 3.98 3.69 3.40 3.10 2.78 2.43 2.09 1.77 1.42 1.09 0.68 0.30];
求dy/dt并绘图....谢谢啦....
MATLAB中没有直接提供求数值导数的函数,只能通过差分估算。
DX=diff(X) 计算向量X的向前差分,DX(i)=X(i+1)-X(i),0<i<n。
t=[0:0.04:1.48];
y=[7.86 7.84 7.82 7.77 7.72 7.68 7.61 7.51 7.42 7.33 7.21 7.07 6.94 6.79 6.64 6.48 6.29 6.11
5.92 5.72 5.50 5.27 5.03 4.78 4.53 4.25 3.98 3.69 3.40 3.10 2.78 2.43 2.09 1.77 1.42 1.09 0.68
0.30];
Dy=diff(y)/.04;
Dt=[0.04:0.04:1.48];
p=polyfit(t,y,2);
Y=polyval(p,t);
DY=diff(Y)/.04;
plot(t,y,'b*',t,Y,'r-',Dt,Dy,'bs',Dt,DY,'g-')
xlabel('t')
ylabel('y&y''')
legend('y值','y拟合2阶多项式','y''差分值','y''拟合值')
扩展资料
clc;clear all
h=0.01;
%x属于【a,b】
a=-5;b=5
x=a:h:b;
n=length(x);
%定义y
y=sin(0.3*x).*cos(3*x);
hold on
grid on
yx=zeros(1,n);
yxx=zeros(1,n);
for i=2:n-1
yx(i-1)=(y(i+1)-y(i-1))/(2*h);
yxx(i-1)=(y(i+1)+y(i-1)-2*y(i))/h^2;
end
plot(x,y,'r','linewidth',2)
plot(x(2:n-1),yx(1:n-2),'g','linewidth',2);
plot(x(2:n-1),yxx(1:n-2),'b','linewidth',2);
legend('原函数','差分一阶导数','差分二阶导数')
xlabel('x','Interpreter','latex','color','r','fontsize',28);
ylabel('y','Interpreter','latex','color','r','fontsize',28);
参考资料:百度百科 差分法
分析如下:
MATLAB中没有直接提供求数值导数的函数,只能通过差分估算。
DX=diff(X) 计算向量X的向前差分,DX(i)=X(i+1)-X(i),0<i<n。
t=[0:0.04:1.48];
y=[7.86 7.84 7.82 7.77 7.72 7.68 7.61 7.51 7.42 7.33 7.21 7.07 6.94 6.79 6.64 6.48 6.29 6.11 5.92 5.72 5.50 5.27 5.03 4.78 4.53 4.25 3.98 3.69 3.40 3.10 2.78 2.43 2.09 1.77 1.42 1.09 0.68 0.30];
Dy=diff(y)/.04;
Dt=[0.04:0.04:1.48];
p=polyfit(t,y,2);
Y=polyval(p,t);
DY=diff(Y)/.04;
plot(t,y,'b*',t,Y,'r-',Dt,Dy,'bs',Dt,DY,'g-')
xlabel('t')
ylabel('y&y''')
legend('y值','y拟合2阶多项式','y''差分值','y''拟合值')
资料拓展:
1、MATLAB 是美国MathWorks公司出品的商业数学软件,用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境,主要包括MATLAB和Simulink两大部分。
2、MATLAB是matrix&laboratory两个词的组合,意为矩阵工厂(矩阵实验室)。是由美国mathworks公司发布的主要面对科学计算、可视化以及交互式程序设计的高科技计算环境。它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中,为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案。
资料来源:百度词条matlab
本回答被网友采纳能不能解释一下
Dy=diff(y)/.04;
Dt=[0.04:0.04:1.52];
Dt=[0.04:0.04:1.48];
还有画出来的dt-dy图不光滑,如何把dt-dy图拟合成一条光滑曲线?
谢谢您的回答...
Dt=[0.04:0.04:1.52]; 这句不要,直接复制的历史命令,忘删了。。。。
数值导数就是在每点t(i)处△y(i)/△t(i)阿,我这里diff只是计算△y,△t是已知的不用算了,当然你要写成/diff(t)那也是一样的。。。
图不光滑是自然的,原始数据有偏差阿。虽然看图y是光滑的,但你用
bar(diff(y))看一下会比较明显,等距△t下△y序列不光滑。再用bar(diff(y,2))看一下二阶差分,一直在正负波动,plot绘的连接折线自然是非凸的了。
至于y看起来光滑,y'看起来会很不光滑,这也正体现了导数对误差的敏感性,因为△t很小,△y的小偏差会急剧放大。
所以如果你想要光滑的y'曲线,就不能直接用数值差分估计导数,应该先拟合原始数据,再用拟合表达式来直估导数。这样导数估计会准确些。
看了下你这个数据用二阶多项式回归拟合已经很好了。重写代码如下:
=========================================
t=[0:0.04:1.48];
y=[7.86 7.84 7.82 7.77 7.72 7.68 7.61 7.51 7.42 7.33 7.21 7.07 6.94 6.79 6.64 6.48 6.29 6.11 5.92 5.72 5.50 5.27 5.03 4.78 4.53 4.25 3.98 3.69 3.40 3.10 2.78 2.43 2.09 1.77 1.42 1.09 0.68 0.30];
Dy=diff(y)/.04;
Dt=[0.04:0.04:1.48];
p=polyfit(t,y,2);
Y=polyval(p,t);
DY=diff(Y)/.04;
plot(t,y,'b*',t,Y,'r-',Dt,Dy,'bs',Dt,DY,'g-')
xlabel('t')
ylabel('y&y''')
legend('y值','y拟合2阶多项式','y''差分值','y''拟合值')
=========================================
很感谢您的详细回答
将p=polyfit(t,y,2);改为p=polyfit(t,y,6);是一条光滑的曲线,
如果另有y1=3t+5
如何将dy1/dt的图像和以上Dt-DY图画在同一个坐标系中?谢谢回答........
。。。。。很无语。。。。。。
你弄个6阶多项式有什么必要????2阶不论是y还是y'都已经拟合得很光滑了阿
建模应该尽可能减少不必要的复杂度。
要添加图
t=[0:0.04:1.48];
y=[7.86 7.84 7.82 7.77 7.72 7.68 7.61 7.51 7.42 7.33 7.21 7.07 6.94 6.79 6.64 6.48 6.29 6.11 5.92 5.72 5.50 5.27 5.03 4.78 4.53 4.25 3.98 3.69 3.40 3.10 2.78 2.43 2.09 1.77 1.42 1.09 0.68 0.30];
Dy=diff(y)/.04;
Dt=[0.04:0.04:1.48];
p=polyfit(t,y,2);
Y=polyval(p,t);
DY=diff(Y)/.04;
y1=3*t+5;
plot(t,y,'b*',t,Y,'r-',Dt,Dy,'bs',Dt,DY,'g-',t,y1,'b--')
xlabel('t')
ylabel('y&y''&y1')
hleg=legend('y值','y拟合2阶多项式','y''差分值','y''拟合值','y1=3t+5');
set(hleg,'Location','SouthWest')
我不要添加t-y1图,而是要添加dy1/dt图(即t-v1)图
因为对于y1=3*t+5(这里y1是我随便举的例子,我需要的比这个更复杂);
要求t-v1,需要用到ezplot函数,所以我就不会了。我没有学过matlab,希望您能帮我耐心解答。谢谢.....
事儿真多阿。。。。这个问题到此为止了。
=========================================
t=[0:0.04:1.48];
y=[7.86 7.84 7.82 7.77 7.72 7.68 7.61 7.51 7.42 7.33 7.21 7.07 6.94 6.79 6.64 6.48 6.29 6.11 5.92 5.72 5.50 5.27 5.03 4.78 4.53 4.25 3.98 3.69 3.40 3.10 2.78 2.43 2.09 1.77 1.42 1.09 0.68 0.30];
Dy=diff(y)/.04;
Dt=[0.04:0.04:1.48];
p=polyfit(t,y,2);
Y=polyval(p,t);
DY=diff(Y)/.04;
plot(t,y,'b*',t,Y,'r-',Dt,Dy,'bs',Dt,DY,'g-')
hold on
syms x
v1=diff(3*x+5);
ezplot(v1)
axis([0,1.5,-12,12])
xlabel('t,x')
ylabel('y,y'',y1''')
hleg=legend('y值','y拟合2阶多项式','y''差分值','y''拟合值','y1''符号求导');
set(hleg,'Location','SouthWest')