MATLAB求分数阶微分的数值解，G

2023-09-15 19:30| 来源: 网络整理| 查看: 265

分数阶微积分学是整数阶微积分学的直接拓展，将一阶导数、二阶导数、一重积分、二重积分等整数阶微积分拓展到0.75阶导数、 $\sqrt{2}$ 阶导数等实数甚至是复数阶的导数或积分。这无疑拓展了微积分学的深度。

对于整数阶微积分，一般可以具有简洁明确的物理意义，比如位移、速度和加速度可以很好地解释一个信号与其整数阶导数之间的关系。然而分数阶微积分却没有那么简洁易懂的物理解释。目前对于分数阶微积分的定义，比较应用广泛的是G-L定义，R-L定义和Caputo定义。

Grünwald-Letnikov定义：

用MATLAB语言编写出Grünwald-Letnikov分数阶微积分数值计算的函数：

function dy = glfdiff(y,t,gam) if strcmp(class(y),'function_handel') y = y(t); end h = t(2)-t(1); w = 1; y = y(:); t = t(:); for j = 2:length(t) w(j) = w(j-1)*(1-(gam+1)/(j-1)); end for i = 1:length(t) dy(i) = w(1:i)*[y(i:-1:1)]/h^gam; end

该函数的调用格式为：y1=glfdiff(y,t,r)。其中t为等间距的时间向量，r为阶次，y可以为函数f(t)在t各点的上的函数值，也可以为f(t)的函数句柄，若r为负值，则计算原函数的-r阶积分。

Riemann-Liouville定义：

f(t)的分数阶积分定义：

f(t)的分数阶微分定义应基于积分定义给出：

下面给出一种R-L定义的分数阶微积分的数值解法的函数：

function [dy,t] = rlfdiff(y,t,r) h = t(2)-t(1); n = length(t); dy1 = zeros(1,n); y = y(:); t = t(:); if r>-1, m = ceil(r)+1; p = m-r; y3 = t.^(p-1); elseif r==-1 n = length(t); yy = 0.5*(y(1:n-1)+y(2:n)).*diff(t); dy(1) = 0; for i=2:n dy(i) = dy(i-1)+yy(i-1); end return else m=-r; y3 = t.^(m-1); end for i = 1:n dy1(i) = y(i:-1:1).'*(y3(1:i)); end if r>-1 dy = diff(dy1,m)/(h^(m-1))/gamma(p); t = t(1:n-m); else dy = dy1*h/gamma(m); end

该函数的调用格式为[y1,t]=rlfdiff(y,t,r)。其调用格式类似刚刚的glfdiff()函数，不同的是该函数返回了两个向量y1和t。因为该函数使用了MATLAB求差分的函数diff()，向量y1的实际长度可能小于y的长度。

Caputo定义：

事实上Caputo分数阶微积分的定义与R-L定义是相同的，对于r阶微积分，记m=[r]，则二值之间的关系为：

Caputo分数阶微积分的数值解函数：

function dy = caputo(y,t,alfa,y0,L) dy = glfdiff(y,t,alfa); if nargin0 q = ceil(alfa); if alfa

【本文地址】

公司简介

联系我们