标签:常微分 向量场 mu dz matlab dx dy gdl
过去有画过常微分方程的向量场,通过向量场能够很形象的看出方程解的状态。
最近过节在家刷视频刷到了3Blue1Brown介绍微分方程的视频。
视频中对钟摆建立的微分方程组通过向量场的形式也很形象的表达了系统状态。
这里用matlab也实现一下,同时对三维情况也做了一个实现。
绘制的方法就是计算方程在二维或三维某个点的方向,然后把方向归一化,画出归一化的向量即可。
二维微分方程组如下:
三维微分方程组如下:
matlab代码如下:
二维情况:
clear all;close all;clc; mu = 0.1; gdl = 1; x = -2:0.4:16; y = -5:0.4:5; [x,y] = meshgrid(x,y); dx = y; dy= -mu*y-gdl*sin(x); d = sqrt(dx.^2+dy.^2); dx = dx./d; dy = dy./d; quiver(x,y,dx,dy); hold on; [t,h] = ode45(@test,[0 100],[0.01 3]); plot(h(:,1),h(:,2),'r') [t,h] = ode45(@test,[0 100],[0.01 2]); plot(h(:,1),h(:,2),'g') grid on; axis equal;
test.m:
function dy=test(t,x) mu = 0.1; gdl = 1; dy=[x(2); -mu*x(2)-gdl*sin(x(1))];
结果:
红色和绿色的线为方程组的两个特解。
三维情况:
clear all;close all;clc; a = 16; b = 4; c = 45; x = -40:6:40; y =-40:6:40; z = 0:6:80; [x,y,z] = meshgrid(x,y,z); dx = a*(y-x); dy = c*x - x.*z-y; dz = x.*y-b*z; d = sqrt(dx.^2+dy.^2+dz.^2); dx = dx./d; dy = dy./d; dz = dz./d; quiver3(x,y,z,dx,dy,dz); hold on; [t,h] = ode45(@test2,[0 30],[12 4 0]); plot3(h(:,1),h(:,2),h(:,3)); grid on; axis equal;
test2.m:
function dy=test2(t,y) a = 16; b = 4; c = 45; dy=[a*(y(2)-y(1)); c*y(1)-y(1)*y(3)-y(2); y(1)*y(2)-b*y(3)];
结果:
标签:常微分,向量场,mu,dz,matlab,dx,dy,gdl 来源: https://www.cnblogs.com/tiandsp/p/14402159.html
本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享; 2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关; 3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关; 4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除; 5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。