数值分析算法 MATLAB 实践 常微分方程求解

简介: 数值分析算法 MATLAB 实践 常微分方程求解

数值分析算法 MATLAB 实践 常微分方程求解

Euler 法及改进算法

function [x,y] = euler(fun,a,b,h,y0)
%一阶常微分方程的一般表达式的右端函数:fun
% 显示欧拉格式
% f是带求函数的一阶导形式
% a,b分别是自变量取值上下限
% y0 是初始条件y(0)
% h是步长
    s = (b - a) / h; % 求步数
    X = zeros(1, s+1);
    Y = zeros(1, s+1);
    X = a:h:b;
    Y(1) = y0;
    for k = 1:s
        Y(k+1) = Y(k) + h * fun(X(k), Y(k))
    end
    x = X';
    y = Y';
end
%% euler求解微分方程
% dfun1 = y^2-y^3;
 [x, y] = euler(@dfun1, 0,5,0.01,0.1);
% dfun1 = y;
% [x, y] = euler(@dfun1, 0,1,0.1,1);
figure
plot(x, y);
title('显示欧拉格式');

%% 微分方程
function dfun1 = dfun1(t,y)
    dfun1 = y^2-y^3;
    %dfun1 = y;
end
function[x,y]=imp_euler(func,a_start,b_end,h_step,y0)
%一阶常微分方程的一般表达式的右端函数:fun
% 显示欧拉格式
% func是带求函数的一阶导形式
% a_start,b_end分别是自变量取值上下限
% y0是初始条件y(0)
% h_step是步长
x = a_start : h_step : b_end;
N = length(x);
y = zeros(1, N);
y(1) = y0;
for i = 2:N
    % 显式 Euler 作为初始值迭代计算
    yi_0 =  y(i-1) + h_step * func(x(i-1), y(i-1));
    yi_1 = y(i - 1) + h_step * func(x(i), yi_0);
    while abs(yi_1 - yi_0) > 1e-6
        yi_0 = yi_1;
        yi_1 = y(i - 1) + h_step * func(x(i), yi_0);
    end
    y(i) = yi_1;
end
function[x,y]=improve_euler(func,a_start,b_end,h_step,y0)
%一阶常微分方程的一般表达式的右端函数:fun
% 显示欧拉格式
% func是带求函数的一阶导形式
% a_start,b_end分别是自变量取值上下限
% y0是初始条件y(0)
% h_step是步长
x = a_start : h_step : b_end;
N = length(x);
y = zeros(1, N);
y(1) = y0;
for i = 2:N
    yp = y(i-1) + h_step * func(x(i-1), y(i-1));
    yq = y(i-1) + h_step * func(x(i), yp);
    y(i) = 0.5 * (yp + yq);
end

Runge-Kutta 算法

4阶-单变量龙格库塔公式

rk45
rk45

4阶-多变量龙格库塔公式

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

% 单变量龙格库塔Runge_kutta 经典法
%一阶常微分方程的一般表达式的右端函数:func
% func是带求函数的一阶导形式
% a_start,b_end分别是自变量取值上下限
% y0是初始条件y(0)
% h_step是步长
function[x,y]=Runge_kutta(func,a_start,b_end,h_step,y0)
x = a_start : h_step : b_end;
N = length(x);
y = zeros(1, N);
y(1) = y0;
for i = 2:N
    k1 = func(x(i-1), y(i-1));
    k2 = func(x(i-1) + h_step/2, y(i-1) + h_step/2*k1);
    k3 = func(x(i-1) + h_step/2, y(i-1) + h_step/2*k2);
    k4 = func(x(i-1) + h_step, y(i-1) + h_step*k3);
    y(i) = y(i-1) + h_step/6*(k1 + 2*k2 + 2*k3 + k4);
end
function[x,y]=Runge_kutta45(dyfunc,xspan,y0,h_step)
% 多变量龙格库塔Runge_kutta45
% h_step是步长常选取为0.01;
% ufunc是函数名;
% x0是初始时间值;
% y0是初始化值; 
% n 是迭代步数;
      x = xspan(1):h_step:xspan(2);
      y = zeros(length(y0),length(x));
      y(:,1) = y0(:);
      %循环迭代数值求解部分
     for n=1 : (length(x)-1)
          k1=feval(dyfunc, x(n),y(:,n));
          k2=feval(dyfunc, x(n)+h_step/2,y(:,n)+h_step/2*k1);
          k3=feval(dyfunc, x(n)+h_step/2,y(:,n)+h_step/2*k2);
          k4=feval(dyfunc, x(n+1),y(:,n)+h_step*k3);
          y(:,n+1)=y(:,n)+h_step*(k1+2*k2+2*k3+k4)/6; 
          %按照4阶多变量龙格库塔方法进行数值求解
     end
end
clc;
clear all;
y0=[0,2,9];%初值
xspan = [0,200];%求解区间
h_step = 0.001;%ode45是变步长的算法
[x,y] = Runge_kutta45(@lorenz_diff,xspan,y0,h_step);
figure(1);
plot3(y(1,:),y(2,:),y(3,:),'.');title("x-y-z");
figure(2);
plot3(y(1,:),y(3,:),y(2,:),'.');title("x-z-y");
figure(3);
plot3(y(2,:),y(1,:),y(3,:),'.');title("y-x-z");
function dydt = lorenz_diff(t,y)
%{
    x-->y(1),y-->y(2),z-->y(3)
%}
dydt = [        10*(y(2)-y(1));
             -y(1)*y(3)+30*y(1)-y(2)
              y(1)*y(2)-8/3*y(3)];  
end

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

Matlab 函数库求解

[t, Xt] = ode45(odefun, tspan, X0)
odefun是函数句柄,可以是函数文件名,匿名函数句柄或内联函数名
tspan是区间 [t0 tfinal] 或者一系列散点[t0,t1,…,tf]
X0是初始值向量
t返回列向量的时间点
Xt返回对应T的求解列向量

Lorenz系统
在这里插入图片描述

function dydt = lorenz_diff(t,y)
%{
    x-->y(1),y-->y(2),z-->y(3)
%}
dydt = [        10*(y(2)-y(1));
             -y(1)*y(3)+30*y(1)-y(2)
              y(1)*y(2)-8/3*y(3)];  
end
clc;
y0 = [0,2,9];
[t,y] = ode45('lorenz_diff',[0,200],y0); 
%% 调用ode45绘制Lorenz系统 2D
figure(1);
plot(y(:,1),y(:,3),'.');
xlabel('x');ylabel('z');title("x-z");
figure(2);
plot(y(:,1),y(:,2),'.');
xlabel('x');ylabel('y');title("x-y");
figure(3);
plot(y(:,2),y(:,3),'.');
ylabel('y');zlabel('z');title("y-z");

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

%% 火焰传播数学模型求解
clc;
clear all;
delta=0.01;
f=@(t,y)y^2-y^3;
opts=odeset('Reltol',1.e-4);
[t1,y1]=ode45(f,[0  2/delta], delta, opts);
figure(1)
plot(t1,y1,'-','Marker','.');
title('数值解曲线');
ylabel('y'); xlabel('t');

在这里插入图片描述

目录
相关文章
|
4天前
|
算法 数据安全/隐私保护 计算机视觉
基于FPGA的图像双线性插值算法verilog实现,包括tb测试文件和MATLAB辅助验证
本项目展示了256×256图像通过双线性插值放大至512×512的效果,无水印展示。使用Matlab 2022a和Vivado 2019.2开发,提供完整代码及详细中文注释、操作视频。核心程序实现图像缩放,并在Matlab中验证效果。双线性插值算法通过FPGA高效实现图像缩放,确保质量。
|
3天前
|
算法
基于SOA海鸥优化算法的三维曲面最高点搜索matlab仿真
本程序基于海鸥优化算法(SOA)进行三维曲面最高点搜索的MATLAB仿真,输出收敛曲线和搜索结果。使用MATLAB2022A版本运行,核心代码实现种群初始化、适应度计算、交叉变异等操作。SOA模拟海鸥觅食行为,通过搜索飞行、跟随飞行和掠食飞行三种策略高效探索解空间,找到全局最优解。
|
1天前
|
传感器 算法
基于GA遗传算法的多机无源定位系统GDOP优化matlab仿真
本项目基于遗传算法(GA)优化多机无源定位系统的GDOP,使用MATLAB2022A进行仿真。通过遗传算法的选择、交叉和变异操作,迭代优化传感器配置,最小化GDOP值,提高定位精度。仿真输出包括GDOP优化结果、遗传算法收敛曲线及三维空间坐标点分布图。核心程序实现了染色体编码、适应度评估、遗传操作等关键步骤,最终展示优化后的传感器布局及其性能。
|
3天前
|
算法 数据可视化 数据安全/隐私保护
一级倒立摆平衡控制系统MATLAB仿真,可显示倒立摆平衡动画,对比极点配置,线性二次型,PID,PI及PD五种算法
本课题基于MATLAB对一级倒立摆控制系统进行升级仿真,增加了PI、PD控制器,并对比了极点配置、线性二次型、PID、PI及PD五种算法的控制效果。通过GUI界面显示倒立摆动画和控制输出曲线,展示了不同控制器在偏转角和小车位移变化上的性能差异。理论部分介绍了倒立摆系统的力学模型,包括小车和杆的动力学方程。核心程序实现了不同控制算法的选择与仿真结果的可视化。
31 15
|
5天前
|
传感器 算法 物联网
基于粒子群算法的网络最优节点部署优化matlab仿真
本项目基于粒子群优化(PSO)算法,实现WSN网络节点的最优部署,以最大化节点覆盖范围。使用MATLAB2022A进行开发与测试,展示了优化后的节点分布及其覆盖范围。核心代码通过定义目标函数和约束条件,利用PSO算法迭代搜索最佳节点位置,并绘制优化结果图。PSO算法灵感源于鸟群觅食行为,适用于连续和离散空间的优化问题,在通信网络、物联网等领域有广泛应用。该算法通过模拟粒子群体智慧,高效逼近最优解,提升网络性能。
|
2天前
|
机器学习/深度学习 算法 安全
基于深度学习的路面裂缝检测算法matlab仿真
本项目基于YOLOv2算法实现高效的路面裂缝检测,使用Matlab 2022a开发。完整程序运行效果无水印,核心代码配有详细中文注释及操作视频。通过深度学习技术,将目标检测转化为回归问题,直接预测裂缝位置和类别,大幅提升检测效率与准确性。适用于实时检测任务,确保道路安全维护。 简介涵盖了算法理论、数据集准备、网络训练及检测过程,采用Darknet-19卷积神经网络结构,结合随机梯度下降算法进行训练。
|
5天前
|
机器学习/深度学习 数据采集 算法
基于GWO灰狼优化的CNN-GRU-SAM网络时间序列回归预测算法matlab仿真
本项目基于MATLAB2022a,展示了时间序列预测算法的运行效果(无水印)。核心程序包含详细中文注释和操作视频。算法采用CNN-GRU-SAM网络,结合灰狼优化(GWO),通过卷积层提取局部特征、GRU处理长期依赖、自注意力机制捕捉全局特征,最终实现复杂非线性时间序列的高效预测。
|
6月前
|
安全
【2023高教社杯】D题 圈养湖羊的空间利用率 问题分析、数学模型及MATLAB代码
本文介绍了2023年高教社杯数学建模竞赛D题的圈养湖羊空间利用率问题,包括问题分析、数学模型建立和MATLAB代码实现,旨在优化养殖场的生产计划和空间利用效率。
272 6
【2023高教社杯】D题 圈养湖羊的空间利用率 问题分析、数学模型及MATLAB代码
|
6月前
|
存储 算法 搜索推荐
【2022年华为杯数学建模】B题 方形件组批优化问题 方案及MATLAB代码实现
本文提供了2022年华为杯数学建模竞赛B题的详细方案和MATLAB代码实现,包括方形件组批优化问题和排样优化问题,以及相关数学模型的建立和求解方法。
162 3
【2022年华为杯数学建模】B题 方形件组批优化问题 方案及MATLAB代码实现
|
6月前
|
数据采集 存储 移动开发
【2023五一杯数学建模】 B题 快递需求分析问题 建模方案及MATLAB实现代码
本文介绍了2023年五一杯数学建模竞赛B题的解题方法,详细阐述了如何通过数学建模和MATLAB编程来分析快递需求、预测运输数量、优化运输成本,并估计固定和非固定需求,提供了完整的建模方案和代码实现。
138 0
【2023五一杯数学建模】 B题 快递需求分析问题 建模方案及MATLAB实现代码

热门文章

最新文章