偏微分方程数值解上机实验报告(matlab做的).doc

1、偏微分方程数值解法上机报告（一）一、实验题目：用ritz-galerkin方法求解边值问题的第n次近似，基函数.二、实验目的：通过本次上机实验，理解求解初值问题的变分问题的最重要的近似解法ritz-galerkin方法，以便为学习有限元法打好基础。此外，要熟悉用matlab解决数学问题的基本编程方法，提高运用计算机解决问题的能力。三、实验代码：n=5;syms x;for i=1:n p(i)=sin(i*pi*x); q(i)=-i2*pi2*sin(i*pi*x);endfor i=1:n b(i)=2*int(p(i),0,1); for j=1:n a(i,j)=int(-q(j)+p

2、(j)*p(i),0,1); endendt=inv(a)*b四、运行结果：t= 2251799813685248/3059521645650671/pi 0 281474976710656/9481460623939047/pi 0 281474976710656/43582901062631895/pi五、总结：通过本次上机，我了解了ritz-galerkin方程 ，明白了用ritz-galerkin方法解决边值问题的变分问题的基本原理，并接近一步提高自己的编程动手能力，受益匪浅。偏微分方程数值解法上机报告（二）一、 实验题目：用线性元求下列边值问题的数值解二、 实验目的：通过本次上机，熟

3、悉和掌握用galerkin法观点出发导出的求解处置问题数值解的线性有限元法。增强用matlab解决数学问题的能力。三、 实验代码：n=10; a=0;b=1;h=(b-a)/n; p=1;q=pi2/4; syms s;f=2*sin(pi/2*s); x=0:(b-a)/n:1; b=; for i=1:n b(i)=h*int(f*(x(i)+h*s)*s,a,b)+h*int(f*(x(i+1)+h*s)*(1-s),a,b);end a=; for i=1:n-1 for j=1:n if i-j=-1 a(i,j)=neiji(1,j,n); elseif i-j=0 a(i,j)=

4、neiji(2,j,n); elseif i-j=1 a(i,j)=neiji(3,j,n); end endenda(n,n-1)=neiji(3,n-1,n);a(n,n)=neiji(4,n,n); u=inv(a)*b; ufunction t=neiji(index,j,n)p=1;q=pi2/4;a=0;b=1;h=(b-a)/n;syms s;x=0:h:1;if index=1 t=int(-p*(x(j)+h*s)/h+h*q*(x(j)+h*s)*(1-s)*s,a,b);elseif index=2t=int(-p*(x(j)+h*s)/h+h*q*(x(j)+h*s)*

5、s*s,a,b)+int(-p*(x(j+1)+h*s)/h+h*q*(x(j+1)+h*s)*(1-s)*(1-s),a,b);elseif index=3 t=int(-p*(x(j+1)+h*s)/h+h*q*(x(j+1)+h*s)*(1-s)*s,a,b);elseif index=4 t=int(p*(x(10)+h*s)/h+h*q*(x(10)+h*s)*s*s,a,b);end四、 运行结果：ans = -0.0086 0.0029 -0.0097 0.0036 -0.0101 0.0038 -0.0101 0.0037 -0.0100 0.0034五、 总结：通过本次上机，

6、使我理解了线性有限元法的基本原理和方法。另外，我也懂得了按galerkin方法推导有限元方程的优点，它比ritz法更加方便直接。我也对虚功原理有了初步的认识。因为galerkin方法基于虚功原理，所以不但可用于保守场问题，也可使用于非保守场即非驻定问题。偏微分方程数值解法上机报告（三）实验题目：用线性元求下列问题的数值解（精确到小数点后第四位）实验目的：通过本次上机，掌握二阶椭圆方程的有限元法，进一步熟悉有限元计算的有关问题。实验步骤： 1.在matlab中输入pdetool 2.在弹出的pdetool工具箱中输入求解区域，在object dialog对话框中输入left为-1，bottom为

7、-1，width为2，height2，单击ok按钮。 3. 设置边界条件：左、右边界用neumann条件，左边界输入g为1，q为0，右边界输入g为0，q为0；上、下边界用dirichlet条件，输入h为1，r为0，作网格剖分。设置方程类型为椭圆形，键入c=-1，a=0，f=-2，d=0。 4. 网格剖分 单击工具，或者单击mesh菜单中initialize mesh选项，可进行初始网格剖分。 5. 解方程 单击工具，显示方程色彩解。如图：6单击mesh菜单中export mesh，选择默认值。7. 输出解的数值 单击solve菜单中export solution选项，在打开的export对话框

8、中输入u，单击ok按钮确定。部分节点如下： columns 23 through 33 0.8000 0.6000 0.4000 0.2000 0 -0.2000 -0.4000 -0.6000 -0.8000 -1.0000 -1.0000 -1.0000 -1.0000 -1.0000 -1.0000 -1.0000 -1.0000 -1.0000 -1.0000 -1.0000 -0.8000 -0.6000 columns 45 through 55 -0.8710 -0.4358 -0.4426 0.2619 0.4640 -0.6894 0.8228 0.6909 0.0468 0

9、.3862 -0.0469 0.8741 0.4026 -0.3048 -0.4881 0.2602 -0.8257 -0.6889 0.8274 0.5240 -0.0691 -0.3720部分数值解如下23：330 0 0 0 0 0 0 0 0 -0.0024 0.097645：550.0667 0.5645 0.6069 0.6724 0.8447 0.1668 0.4912 0.2952 0.6128 0.8902 0.7050总结：有限元计算的有关问题有：把初值问题化为变分形式，对求解域作网络分割，构造基函数（或单元形状函数），形成有限元方程。通过本次实验，我懂得了用有限元方法求解

10、初值问题的基本数学思想和方法，也增强了编程能力，提高了用计算机解决数学问题的兴趣。偏微分方程数值解法上机报告（四）一、 实验题目：设g是一个十字形区域，有五个相等的单位正方形组成，用五点差分格式求下列边值问题的数值解：二、 实验目的：通过本次上机，掌握椭圆型方程的有限差分法，熟悉其计算过程与基本的思想。三、 实验代码：h=0.125;a=zeros(6,14);for i=1:14a(1,i)=0;a(6,i)=0;endfor j=1:6a(j,1)=0;a(j,14)=0;enda(2,11)=0;a(2,12)=0;a(2,13)=0;a(3,12)=0;a(3,13)=0;a(4,13

11、)=0;n=0;for i=2:5for j=2:9+i-1 a(i,j)=h2/4+(a(i-1,j)+a(i,j-1)+a(i+1,j)+a(i,j+1)/4; n=n+1;endendan四、运行结果：a = columns 1 through 8 0 0 0 0 0 0 0 0 0 0.0039 0.0049 0.0051 0.0052 0.0052 0.0052 0.0052 0 0.0049 0.0063 0.0068 0.0069 0.0069 0.0069 0.0069 0 0.0051 0.0068 0.0073 0.0075 0.0075 0.0075 0.0075 0 0.0052 0.0069 0.0075 0.0076 0.0077 0.0077 0.0077 0 0 0 0 0 0 0 0 columns 9 through 14 0 0 0 0 0 00.0052 0.0052 0 0 0 00.0069 0.0069 0.0056 0 0 00.