有限差分求解扩散的数值解

1、有限差分法求解扩散问题的 数值解 王智文 2013.5.15提 纲1.有限差分原理（一维坐标系）2.边界条件处理（二维柱坐标系）3.实验室应用举例（三维直角坐标系）1.有限差分的基本原理（一维坐标系）微分方程转化为代数方程：菲克第二定律的一维方程：1.有限差分的基本原理 （变量离散）X轴t轴m=L/k，n=T/h,平面内有mn个网格点LCi+1,jCi,j-1Ci,jCi-1,jCi,j+1tj=j*h在xt平面将距离和时间离散化：时间步长t取h，距离步长x取k，步数分别为n和m。Xi=i*k设法在xt平面网格点上逼近浓度C(x,t)值1.有限差分的基本原理（ 将微分方程转化为代数方程） 令，

2、得代数迭代格式由上式可知，某一时刻的浓度值，可以由前一个时刻的浓度值计算得出，即可从初始零时刻逐步迭代计算任一时刻任一离散点的浓度值注意：要使迭代式收敛稳定，须使0s 0.5Ci,jCi,j+1Ci+1,jCi-1,j高浓度端 10-5低浓度端 01.有限差分的基本原理（ 一维扩散例子） 设定参数 扩散系数D： 10-7cm2/s 透氧距离L： 0.01cm 高浓度C1： 10-5 mol/ml 低浓度C2： 0 距离步长k： 0.0005cm 距离等分数： 20 时间步长h： 1s 初始浓度： 0迭代代数式，Matlab编程计算1.有限差分的基本原理（ 一维扩散例子） 计算结果 稳态下扩散出

3、口绝缘 环扩散进口圆柱上底面和侧面为高浓度面（红色）圆柱下底面同心小圆为低浓度面（白色）圆柱下底面圆环为绝缘面 （蓝色）2.边界条件处理（ 二维柱坐标扩散例子） rZ0圆心对称2.边界条件处理（ 二维柱坐标扩散例子） 扩散微分方程： 柱坐标下取体积微元由质量衡算推导扩散偏微分方程rr+ rZ体积微元质量变化量移项化简得即代入可得，将2.边界条件处理（ 二维柱坐标扩散例子） Cki,jCki,j-1Cki,j+1Cki+1,jCki-1,j0变量离散，微分方程化为代数方程：Cki-1,jCki,j-1Cki,j+1Cki,jCki+1,j2.边界一Ck1,1Ck1,n=C0Ck1,nCkm,n=

4、C0圆柱上表面和侧面边界扩散出口边界条件扩散出口边界Ckm,1Ckm,n1=02.边界二下表面绝缘边界条件：Ckm-1,jCkm,j-1Ckm,j+1Ckm,j虚拟点：Ckm+1,j=Ckm,j2.边界三Cki-1,jCki,j-1Cki,j+1Cki,jCki+1,j圆柱轴中心线上格点Cki-1,1Cki+1,1Cki,1Cki,2无法用2.边界四圆柱轴中心线上格点微元体积质量衡算Cki-1,1Cki,1Cki+1,1Cki,22.边界四体积微元质量变化量移项化简得得到轴中心线上网格点迭代格式：2.Matlab计算结果50秒浓度及其等高线2.Matlab计算结果200秒浓度及其等高线2.Ma

5、tlab计算结果500秒浓度及其等高线2.Matlab计算结果1200秒浓度及其等高线3.扩散模拟表征多孔介质物质输运性能 实验室制备多孔陶瓷膜A：样品实物图B：SRCT三维重构图C：孔结构图ABCDCBA固相孔A：与样品上下表面皆连通的孔 B：与样品上表面连通的孔C：封闭在样品内的孤立孔 D：与样品下表面连通的孔 3.扩散模拟表征多孔介质物质输运性能 孔连通性DCBA固相孔通过判断稳态时各网格点浓度值区分四类孔初始浓度设为非0非1的数，以区分孤立孔死孔包含B、C和D上表面浓度：1下表面浓度：00Cs1气体扩散测定曲折因子C1C2C1C2无孔结构多孔介质C1C2J0Jeff曲折因子其中为孔隙率

6、SRCT构建孔结构数组灰度图二值图三维结构重构图A黑色：孔白色：实体孔隙率：0.353.扩散模拟表征多孔介质物质输运性能（相转换流延制备） 模拟计算原理原理：三维扩散方程计算方法：有限差分令x,y,z方向步长值相等，简化上式得Ci,j-1,kCi,j,k+1Ci-1,j,kCi,j+1,kCi,j,k-1Ci+1,j,kzyx3.扩散模拟表征多孔介质物质输运性能 Ci,j-1,kCi,j,k+1Ci-1,j,kCi,j+1,kCi,j,k-1Ci+1,j,k含孔结构计算处理技巧计算前，先判断该点三维方向六个点Ai 1,j1,k1的结构，如果周围某个点为实体，则令其浓度等于计算点Ai,j,k的浓

7、度Ci,j,k例如： Ai +1,j,k ， Ai -1,j,k ， Ai,j+1,k， Ai,j-1,k四个点为实体，则令其Ci +1,j,k = Ci -1,j,k = Ci,j+1,k= Ci,j-1,k = Ci,j,k则迭代时简化为为一维扩散方程离散代数表达式xyzCi,j-1,kCi,j,k+1Ci-1,j,kCi,j+1,kCi,j,k-1Ci+1,j,k3.扩散模拟表征多孔介质物质输运性能（相转换流延制备） 根据以上迭代式逐个求每个时间点网格点的浓度，直到进口流量和出口流量相等，即到达稳态。浓度不变的网格点构成的孔为死孔稳态时浓度分布图内部孔结构3.扩散模拟表征多孔介质物质输运性能 纵向剖面浓度分布纵向剖面结构图3.扩散模拟表征多孔介质物质输运性能（相转换流延制备） 3.扩散模拟表征多孔介质物质输运性能（孔随机分布） 孔随机分布三维结构 黑色：孔白色：实体孔隙率：0.353.扩散模拟表征多孔介质物质输运性能（孔随机分布） 稳态时浓度分布图纵向剖面浓度分布纵向剖面结构图3.扩散模拟表征多孔介质物