-
如何使用MATLAB求解微分方程组_非线性微分方程组求解
2020-07-09 08:31:36TMU_BME_2013 Topic: 如何使用 MATLAB 求 解常微分方程组 a.What ? 微分方程 指描述未知函数的导数与自变 量之间的关系的方程未知函数是一元函 数的微分方程称作 常微分方程 未知函数 是多元函数的微分方程称作 偏... -
偏微分方程 python_有关一个偏微分方程组的求解?
2020-12-24 14:04:14写了个MATLAB的小程序,用特征线法求解 偏微分方程组。[XX,TT]=meshgrid(0:0.4:4,0:0.1:1);N=size(XX,2);T=size(XX,1);u1=zeros(T,N);u2=zeros(T,N);dx=0.4;dt=0.1;ds=0.1;u10=1:N;u20=(N:-1:1); % initial valueu...题目写的看不懂。我也瞎答。
写了个MATLAB的小程序,用特征线法求解 偏微分方程组。
[XX,TT]=meshgrid(0:0.4:4,0:0.1:1);
N=size(XX,2);
T=size(XX,1);
u1=zeros(T,N);
u2=zeros(T,N);
dx=0.4;
dt=0.1;
ds=0.1;
u10=1:N;u20=(N:-1:1); % initial value
u10t=1:T;u20t=1:T; % value at a given node point
u1(1,:)=u10;
u2(1,:)=u20;
u1(:,1)=u10t;
u2(:,1)=u20t;
for tt=2:T % 特征线法求解线性偏微分方程组
for ii=2:N
u1(tt,ii)=(-400*u1(tt-1,ii-1)+400*u2(tt-1,ii-1))*ds;
u2(tt,ii)=(0.35*u1(tt-1,ii)-0.35*u2(tt-1,ii))*dt;
end
end
quiver(XX,TT,u1,u2)
-
用python怎样解偏微分方程组_用Python数值求解偏微分方程
2020-12-29 02:16:39然而,只有很少的微分方程可以解析求解,尤其对于偏微分方程,能解析求解的种类更是寥寥可数。更多的微分方程可以采用数值法进行求解,只要精度足够高,就可以满足科学和工程上的需求。数值求解微分方程的基本思路是...1 引言
微分方程是描述一个系统的状态随时间和空间演化的最基本的数学工具之一,其在物理、经济、工程、社会等各方面都有及其重要的应用。然而,只有很少的微分方程可以解析求解,尤其对于偏微分方程,能解析求解的种类更是寥寥可数。更多的微分方程可以采用数值法进行求解,只要精度足够高,就可以满足科学和工程上的需求。
数值求解微分方程的基本思路是先把时间和空间离散化,然后将微分化为差分,建立递推关系,然后利用计算机强大的重复计算能力,快速得到任意格点处的值。Python的Numpy、Scipy工具包可以很好地实现此功能,Matplotlib工具包则可以将求解结果画为非常直观的图形。
接下来,我们先以常微分方程的数值求解为例,引入差分的思想,再将其推广到偏微分方程中。
2 常微分方程的差分求解
一般地,一阶常微分方程可以写为
首先,将连续的变量
和
离散化,连续的函数
和
化为离散的序列
和
,则上述微分方程可以化为差分方程【1】
从而我们得到递推关系
有了递推关系和初始条件以后,就可以利用Python的强大计算功能,得到任意的
的值了,下面我们通过一个具体的例子来说明。
2.1 RC回路放电问题
对于一个
回路,我们有
其中,
分别为电流,电阻,电量和电容,利用
,并定义
,我们得到一个含初始条件的一阶常微分方程
这个方程当然可以解析求解,得到
。另一方面按照差分法,可以得到递推关系
下面我们用Python进行数值求解,并和解析结果进行比较。
import numpy as np
import matplotlib.pyplot as plt
rc = 2.0 #设置常数
dt = 0.5 #设置步长
n = 1000 #设置分割段数
t = 0.0 #设置初始时间
q = 1.0 #设置初始电量
#先定义三个空列表
qt=[] #用来盛放差分得到的q值
qt0=[] #用来盛放解析得到的q值
time = [] #用来盛放时间值
for i in range(n):
t = t + dt
q1 = q - q*dt/rc #qn+1的近似值
q = q - 0.5*(q1*dt/rc + q*dt/rc) #差分递推关系
q0 = np.exp(-t/rc) #解析关系
qt.append(q) #差分得到的q值列表
qt0.append(q0) #解析得到的q值列表
time.append(t) #时间列表
plt.plot(time,qt,'o',label='Euler-Modify') #差分得到的电量随时间的变化
plt.plot(time,qt0,'r-',label='Analytical') #解析得到的电量随时间的变化
plt.xlabel('time')
plt.ylabel('charge')
plt.xlim(0,20)
plt.ylim(-0.2,1.0)
plt.legend(loc='upper right')
plt.show()
图1:用改进的欧拉法差分得到的结果和解析结果的比较
在图1中给出了差分法得到的结果与解析法得到结果的比较,发现两者符合得很好,这说明对于这个问题,改进的欧拉法已经可以给出足够精确的结果。需要注意的是,这个微分方程本身比较简单,可以解析求解,而对于复杂得多的微分方程,没法解析求解,但是上述数值求解差分方法仍然是适用的。
3 偏微分方程的差分求解
有了差分代替微分的思想,接下来我们将其推广到偏微分方程的求解中。以一般二阶抛物型偏微分方程为例,一般的可以写为
仍然是将时间和空间离散化,将微分化为差分,即
其中
和
分别为空间步长和时间步长,
和
分别标记空间指标和时间指标,则我们得到差分方程
由此得到递推关系
下面我们考察一个具体的例子,一维热传导方程的求解。
3.1 一维热传导方程的求解
一维热传导方程是一个典型的抛物型二阶偏微分方程。设
表示在时间
,空间
处的温度,则根据傅里叶定律(单位时间内流经单位面积的热量和该处温度的负梯度成正比),可以导出热传导方程【2】
其中
称为热扩散率,
分别为热导率,比热和质量密度,都是由系统本身确定的常量。
为了具体,设
,设边界条件为
设步长为:
,从而
,所以递推关系为
图2:递推关系示意图
图2直观地给出了差分法求解偏微分方程的过程。先把时空坐标都离散化,根据递推关系,由下一行的三个蓝点的值可以给出上一行的一个红点的值,由于边界条件和初始条件(即最下方和两边的绿线)已知,所以按这个递推关系可以得到网格中的所有值。下面我们用Python代码来实现。
import numpy as np
import matplotlib.pyplot as plt
h = 0.1#空间步长
N =30#空间步数
dt = 0.0001#时间步长
M = 10000#时间的步数
A = dt/(h**2) #lambda*tau/h^2
U = zeros([N+1,M+1])#建立二维空数组
Space = arange(0,(N+1)*h,h)#建立空间等差数列,从0到3,公差是h
#边界条件
for k in arange(0,M+1):
U[0,k] = 0.0
U[N,k] = 0.0
#初始条件
for i in arange(0,N):
U[i,0]=4*i*h*(3-i*h)
#递推关系
for k in arange(0,M):
for i in arange(1,N):
U[i,k+1]=A*U[i+1,k]+(1-2*A)*U[i,k]+A*U[i-1,k]
上述代码中,我们首先把位于0-3中的空间等分为30份,位于0-1的时间等分为10000份,然后通过设置初始条件、边界条件和递推关系并借助for循环就得到了1个30*10000的二维数组,里面放着每个离散的时空点的温度值。
为了直观地展现温度随时空的变化关系,接下来开始画图,首先画出不同时刻温度随空间坐标的变化
#不同时刻的温度随空间坐标的变化
plt.plot(Space,U[:,0], 'g-', label='t=0',linewidth=1.0)
plt.plot(Space,U[:,3000], 'b-', label='t=3/10',linewidth=1.0)
plt.plot(Space,U[:,6000], 'k-', label='t=6/10',linewidth=1.0)
plt.plot(Space,U[:,9000], 'r-', label='t=9/10',linewidth=1.0)
plt.plot(Space,U[:,10000], 'y-', label='t=1',linewidth=1.0)
plt.ylabel('u(x,t)', fontsize=20)
plt.xlabel('x', fontsize=20)
plt.xlim(0,3)
plt.ylim(-2,10)
plt.legend(loc='upper right')
plt.show()
图3:不同时刻的温度随空间坐标的变化
从图3可以看到,温度关于
呈现轴对称分布,这是由初始条件造成的。另外,对每一点的空间坐标,随着时间的推移,温度越来越低。
接下来,我们来画出温度等高线来描述温度随任意时空点的变化
#温度等高线随时空坐标的变化,温度越高,颜色越偏红
extent = [0,1,0,3]#时间和空间的取值范围
levels = arange(0,10,0.1)#温度等高线的变化范围0-10,变化间隔为0.1
plt.contourf(U,levels,origin='lower',extent=extent,cmap=plt.cm.jet)
plt.ylabel('x', fontsize=20)
plt.xlabel('t', fontsize=20)
plt.show()
图4:温度随时空坐标的变化,温度越高,颜色越红
在图4中,利用颜色的深浅来标记温度,温度越高,颜色越红。从中同样可以看到,温度随空间的分布关于
轴对称,而且随着时间的推移,温度越来越低。
4 总结
在本文中,我们利用Python数值求解了常微分方程和偏微分方程,基本思想是先将连续的坐标离散化,然后将微分化为差分,由差分方程得到递推关系,然后利用计算机强大的重复计算能力得到任意格点处的函数值。虽然上面只算了两个例子,但是这种方法完全可以推广到任意偏微分方程的求解中。
在量子色动力学(QCD)中,由于强相互作用具有渐进自由的特性,所以在低能情况下没办法像QED那样使用微扰论计算,这时就要采用格点QCD的方法计算。其基本思想也是将时空离散化,然后从第一性原理的路径积分出发去计算。由于时空被离散化了,相当于人为地引入了一个最小时空距离
,在傅里叶变换到动量空间后相当于引入了一个最大的动量截断
,所以计算结果不会出现紫外发散,从而可以算到很高的精度,在一些情况下,格点的计算结果甚至比实验更精确。所以,将连续参数离散化,把微分化为差分的思想,是极其重要的。
【1】不同的算法对于方程右边具体取什么形式并不一样,从而精度也不一样。例如,欧拉法右边取得是
;改进的欧拉法右边取的是
;二阶Runge-Kutta法右边取的是
;四阶Runge-Kutta法右边取的是
,其中
,
,
,
。这些算法的差别在于计算精度不同,并不改变差分的本质思想。为了具体,我们这里采用改进的欧拉法。
【2】傅里叶定律告诉我们单位时间通过单位面积的热量和该处的温度负梯度成正比,即
,其中
是热流,即单位时间通过单位面积的热量,
为热导率,
为温度。能量守恒定律告诉我们单位时间流出某闭合曲面的热量等于其内部减少的能量,即
,其中
表示单位体积的热容,而
和
分别为质量密度和单位质量的热容(比热),联合散度定理
,我们得到
,再将傅里叶定律带入,就得到了热传导方程
,其中
称为热扩散率。在一维的情况下,热传导方程就退化到了正文中的形式,即
。
-
Matlab求解微分方程(组)及偏微分方程(组)
2020-05-01 06:30:26Matlab求解微分方程(组)及偏微分方程(组) -
Matlab求解微分方程(组)及偏微分方程(组).pdf
2020-03-07 17:32:05(1)PDEtool(GUI)求解偏微分方程的一般步骤 在Matlab命令窗口输入pdetool,回车,PDE工具箱的图形用户界面(GUI)系统就启动了.从定义一个偏微分方程问题到完成解偏微分方程的定解,整个过程大致可以分为六个阶段 -
Matlab偏微分方程求解方法
2019-05-15 10:19:14非稳态的偏微分方程组是一个比较难解决的问题,也是在热质交换等方面的常常遇到的问题,因此需要一套程序来解决非稳态偏微分方程组的数值解。 -
DuFort-Frankel格式求解椭圆-抛物型偏微分方程组
2020-05-13 08:00:57DuFort-Frankel格式求解椭圆-抛物型偏微分方程组,matlab程序,其中椭圆用积分公式,抛物用DuFort-Frankel格式,多多指教 -
用python怎样解偏微分方程组_Python能解偏微分方程吗
2020-12-18 13:43:55由菲克第二定律可以得到动态扩散的偏微分方程。求解可以得到浓度分布和流出曲线。不确定这个问题有没有解析解,不过数值求解是一种较为通用的解决方法。fipy是目前难得的还活着的PDE求解python包,作者根据官方示例...fipy
菲克定律是指在不依靠宏观的混合作用发生的传质现象时,描述分子扩散过程中传质通量与浓度梯度之间关系的定律。菲克定律是阿道夫·菲克(Adolf Fick)于1855年提出。
由菲克第二定律可以得到动态扩散的偏微分方程。求解可以得到浓度分布和流出曲线。
不确定这个问题有没有解析解,不过数值求解是一种较为通用的解决方法。
fipy是目前难得的还活着的PDE求解python包,作者根据官方示例改写本程序。
相关推荐:《python视频教程》
问题
一个二维平板,顶端1摄氏度(100也可以,只是一个系数),另外三个边缘0摄氏度,初始时刻整个板子都是0摄氏度,随之时间的推进,热量在板子上传递,最后达到平衡态,我们不仅希望知道平衡态的温度分布,也希望知道温度随时间是如何变化的。热量的传递由微分方程给出,可以简单地理解为热量按照温度降低最快的方向进行传递。
公式右边是温度的梯度,左边是温度随时间的变化
最后整个板子的温度分布大致呈现怎样
只有一个包需要导入import fipy as fp
确定求解区域,一个20*20的格点
#求解区域nx = 20ny = 20dx = 1.dy = dxL = dx * nxmesh = fp.Grid2D(dx=dx, dy=dy, nx=nx, ny=ny)
phi = fp.CellVariable(name = "solution variable", mesh = mesh, value = 0.)
创建微分方程
设立边界条件
创建画图
求解
建议在命令行里面运行,命令行里面可以获得动图,ipython里面只有最后一张图
spyder的ipython里面只有最后的一张图片
-
微分方程_微分方程 | 线性微分方程组的求解(上)
2020-12-31 08:27:37附:微分方程的分类:常微分方程和偏微分方程。1、常微分方程(ODE)是指微分方程的自变量只有一个的方程。最简单的常微分方程,未知数是一个实数或是复数的函数,但未知数也可能是一个向量函数或是矩阵函数,后者可...附:
微分方程的分类:常微分方程和偏微分方程。
1、常微分方程(ODE)是指微分方程的自变量只有一个的方程。最简单的常微分方程,未知数是一个实数或是复数的函数,但未知数也可能是一个向量函数或是矩阵函数,后者可对应一个由常微分方程组成的系统。
2、偏微分方程(PDE)是指微分方程的自变量有两个或以上,且方程式中有未知数对自变量的偏微分。偏微分方程的阶数定义类似常微分方程,但更细分为椭圆型、双曲线型及抛物线型的偏微分方程,尤其在二阶偏微分方程中上述的分类更是重要。
有些偏微分方程在整个自变量的值域中无法归类在上述任何一种型式中,这种偏微分方程则称为混合型。
常微分方程及偏微分方程都可以分为线性微分方程及非线性微分方程二类。
一个常微分方程是不是有特解呢?如果有,又有几个呢?这是微分方程论中一个基本的问题,数学家把它归纳成基本定理,叫做存在和唯一性定理。因为如果没有解,而我们要去求解,那是没有意义的;如果有解而又不是唯一的,那又不好确定。因此,存在和唯一性定理对于微分方程的求解是十分重要的。 -
matlab偏微分方程工具箱求解
2020-12-14 18:13:50偏微分方程组的matlab求解语句 该命令用以求解以下的PDEPDEPDE方程式: c(x,t,u,∂u∂x)∂u∂t=x−m∂(xmf(x,t,u,∂u∂x))∂x+s(x,t,u,∂u∂x) c(x,t,u,\frac{\partial u }{\partial x})\frac{\partial u}{\... -
运用偏微分方程(PDE)方法进行图像处理的matlab程序,包括图像滤波、图像分割、插值、图像增强、恢复及一些...
2019-08-19 15:10:52运用偏微分方程(PDE)方法进行图像处理的matlab程序,包括图像滤波、图像分割、插值、图像增强、恢复及一些方程组求解等在偏微分方法处理图像处理 -
python解隐式方程_用隐式方法求解偏微分方程组
2021-01-30 02:15:19我有一个偏微分方程组(PDE),特别是应用于传热和对流的扩散-平流-反应方程,我用有限差分法求解。在我的模拟环境中,我有很多不同的部分,比如管道,储能器,热交换器等等。。。根据零件的不同,每个零件的PDE可以是... -
matlab常微分方程数值求解
2020-09-01 21:08:13首先,凡含有参数,未知函数和未知函数导数 (或微分) 的方程,称为微分方程,有时简称为方程,未知函数是一元函数的微分方程称作常微分方程,未知函数是多元函数的微分方程称作偏微分方程。微分方程中出现的未知函数... -
Matlab中的PDEPE求解"瞬态型"或"发展型"非线性偏微分方程组
2016-02-11 19:03:11背景求解真实问题中建模得到的非线性偏微分方程组, 尽可能少手写代码,如果用matlab, 可能的选项很有限。pdetoolbox有限元方法能够求解一些,但是要求对微分方程的类型以及pdetoolbox特有的书写方式非常熟悉,并能够把... -
二阶偏微分方程组 龙格库塔法_常微分法数值计算方法
2021-01-09 13:48:38假设一阶常微分方程组有下式给出 其中矢量 为状态变量 构成的向量,即 ,常称为系统的状态向量,n称为系统的阶次,而 为任意函数数,t为时间变量,这样就可以采用数值方法求解常微分方程组了。另外任意高阶微分方程都... -
matlab 三弯矩方程_数值计算(三)matlab求解一般的偏微分方程组
2020-12-31 08:28:461 pdepe()函数的一般调用格式 [1]sol=pdepe(m,@pdefun,@pdeinit,@pdebound,x,t),其中pdefun是偏微分方的描述函数,标准形式为: 初始条件满足: 边界条件满足: 说明:a,b表示上、下边界;ua,ub是上、下边界的近似... -
基于微分变换方法的具有初值和边界条件的偏微分方程组新算法
2020-06-03 00:18:48该方法可以快速有效地求解具有初始边界值(IBVP)的线性和非线性偏微分方程。 根据边界条件,将初始条件扩展为傅立叶级数。 之后,将IBVP转换为K域中的迭代关系。 可以得到级数解或精确解。 通过比较FDTM获得的结果... -
偏微分方程数值解法python_利用MATLABpdetool解偏微分方程——以数学建模2018A题为例...
2020-12-31 08:41:052.1 建立偏微分方程组 2.2 利用PDE Modeler求解 2.3 误差分析 附录 MATLAB代码1 PDE Modeler使用方法介绍物理学中的偏微分方程(PDE)无处不在,如热传导方程、扩散方程、电磁场方程,甚至量子力学中也能大量... -
一阶微分方程的物理意义_一阶线性偏微分方程特征线解法.pdf
2020-12-19 22:11:16• 偏微分方程理论 主要研究具有实际背景的偏微分方程或偏微分方程组。是数学的基础学科之一。第 1页二、《数学物理方程》课程的特点:1、数学理论、解题方法与物理实际有机结合。可以学到:如何根据物理现象建立偏... -
基于R语言的微分方程求解
2015-11-23 16:04:59采用R语言实现微分方程,偏微分方程以及差分方程及方程组的求解方法 -
求解偏微分方程开源有限元软件deal.II学习--Step 8
2016-10-18 17:55:01求解偏微分方程开源有限元软件deal.II学习-...真实生活中,大多数的偏微分方程都是一组方程,相应地,解也通常是矢量场。跟之前单方程的求解以及解是标量场相比,这种问题只是在组装矩阵和右端项时有些复杂,其他都一样 -
偏微分方程-特征线法(转)
2012-06-03 10:55:57数学中的特征线法是求解偏微分方程的一种方法,适用于准线性偏微分方程的求解。只要初始值不是沿着特征线给定,即可通过特征线法获得偏微分方程的精确解。 其基本思想是通过把双曲线型的准线性偏微分方程转化为... -
运用偏微分方程(PDE)方法进行图像处理的matlab程序
2020-04-28 10:19:38运用偏微分方程(PDE)方法进行图像处理的matlab程序,包括图像滤波、图像分割、插值、图像增强、恢复及一些方程组求解等在偏微分方法处理图像处理领域常用且重要的处理程序 -
方程组的直接解法和迭代法 python_数值偏微分方程-差分法(Python)
2020-12-06 20:12:41类似的方法可以拓展到解偏微分方程的问题,这里整理有限差分法的相关笔记。计算机无法处理连续介质的无限多点,因此只能通过网格划分的方式,将带求解区域内连续分布的函数进行离散化,这样只用求解有限多的格点即可... -
二阶偏微分方程组 龙格库塔法_数值积分公式的推广:龙格库塔法(一)
2021-01-11 23:56:111. 写在前面的话龙格库塔法是求解常微分方程的数值计算方法。它始于1895年,历史悠久,后人以龙格(Runge)和库塔(Kutta)两位计算数学家的名字为它命名。关于它的研究文献有很多,至今仍有新的论文不断涌现,卷帙之... -
求解偏微分方程开源有限元软件deal.II学习--Step 11
2016-10-18 17:58:03求解偏微分方程开源有限元软件deal.II学习--Step 11 Posted on 2016-09-16 | In computational material science | 暂无评论 引子 这个例子的亮点:介绍了最简单的Laplace方程的组装矩阵和右端项的一键... -
偏微分方程数值解的Matlab 实现
2013-05-31 20:22:36工程中有许多问题可以归结为偏微分方程问题,如弹塑性力学中研究对象(结构、边坡等)内部的应力应变问题、地下水渗流问题等。这些由偏微分方程及边界条件、初始条件等组合...Matlab采用有限元法求解偏微分方程的数值解 -
Matlab中求解双曲椭圆一维初边值偏微分方程(组)的pdepe
2014-03-09 09:51:32http://www.mathworks.cn/help/matlab/ref/pdepe.html;jsessionid=143095f6a8df9e5c4b9dfd0f8be0 有很多默认的模式,还不如自己手动解; ...Solve initial-boundary value problems for parabolic-elliptic PD
-
获取USB摄像头名字和device ID等信息
-
同学Linux,同成长
-
Glasterfs 分布式网络文件系统
-
在线支付公司Stripe的服务发现架构设计过程分享
-
莫队 学习笔记
-
FTP 文件传输服务
-
SecureCRT 连接 GNS3/Linux 的安全精密工具
-
基于java的图形化的小说阅读器
-
access应用的3个开发实例
-
Galera 高可用 MySQL 集群(PXC v5.7+Hapro)
-
C++MFC开发远程控制软件教程(VS2013)
-
牛客19935[CQOI2014]通配符匹配
-
nmap大全
-
基于双向分布反射函数的红外偏振特性分析
-
jlink-v8 固件.zip
-
激光等离子体X射线偏振度探测
-
Windows系统管理
-
向狂神学习的第3天
-
Windows系统关机故障
-
LeetCode867. 转置矩阵