偏微分方程数值解法

一、问题ffinoninuuu2tu0uxu0二、问题分析第一步利用Green公式，求出方程的变分形式变分形式为：求uL20,T;H1，使得0u,vu,vu2,vf,v vH1 （）t 0和 uxu.0第二步对空间进行离散，得出半离散格式对区域

,

,,

，那么（*）的Galerkin逼近为：

h 1 2 NGt0,T，求ut,xVh

H1，使得0d,vdt h h

t,vut2,vf,vh h h

Vh h

（**）和u0u

，u0,

为初始条件u0

在V中的逼近，设uh

为u0,h

在V.h则t0，有u

tNGt，u

=NG

，代人（**）即可取得一常微分方程组.h i1

0,

i1

0i i第三步进一步对时刻进行离散，取得全离散的逼近格式对du等分为n个小区间,tdt i1

，其长度tti

ti1

T,tn

T. 如此把求t时刻的近似记为ui

的近似.那个地址对（**）采纳向后的欧拉格式，即i1

h 0 ui1ui,v

ui1,v ui1 ,

f

,v V

(***)t h h h

h h h h

ih h hi=0,1,2…,n-1. u0=h

0,h由于向后欧拉格式为隐式格式且含有非线性项，故相邻两时刻步之间采纳牛顿迭代，即：u un1 k h ,

,v

2u

u,

2,

f,v h k h

kk h

kh n hu其顶用un1作为牛顿迭代的初值.hu三、计算流程取u100t1xyx1y1作为方程的准确解，算得初值函数：u 100xyx1y10右端函数：f100xyx1y1200t1xx1200t1yy1+10000t2x2y2x2y2.取y;0x1,0y对x,y方向对[0,1]N等分并将每一个小正方形沿对角线分为两个三角形.对节点进行整体编码（那个地址按先沿x方向再沿y方向进行顺序编码..阵...n个时刻步利用牛顿迭代法算得第n+1个时刻步的数值解，取u0h

u .0,h计算各个时刻步的有限元数值解和L2误差. L2误差utn

unh0

的阶数为th2

其中t为时刻步长，h为空间步长那个地址取N=25,那么元素总数LEE=2*N*N=1250,节点总数NG=(N+1)(N+1)=676,取时刻步长TSTP=,时刻步数TN=100,即在t=1时，算得结果为：即当t=,h=时，L2误差为阶数为th2附C程序#include""#include""#include""#defineN25 =1;b[i].y=(i/(N+1)-1)*N);}else{}}

b[i].x=(i%(N+1)-1)*N);b[i].y=(i/(N+1))*N);voidboundnote(intbd[],structxyb[])==0||b[i].x==||b[i].y==0||b[i].y=={bd[j]=i;j++;}}}void dealwithbd(double **uk,int bd[])+b[a[i][2]].x+b[a[i][3]].x)/,(b[a[i][1]].y+b[a[i][2]].y+b[a[i][3]].y)/,(n+1)*TSTP);}voidUFE(doubleufe[],double**fe,int**a),b[i].y);double*aryn1;for(n=0;n<TN;n++){aryn1=(double*)calloc(NG+1,sizeof(double));NEWTONITERATIVE(aryn1,aryn,n,a,b,bd);for(i=0;i<NG+1;i++)aryn[i]=aryn1[i];free(aryn1);}for(i=1;i<NG+1;i++)fesolution[i]=aryn[i];printf("\n第%d步的结果为：\n\n",n);printf("节点整体编码 有限元值 精准值 误差 for(j=N/2;j<NG+1;j+=2*(N+1))printf("%8d%15f%15f%15f\n",j,fesolution[j],u(b[j].x,b[j].y,TSTP*TN),fesolution[j]-u(b[j].x,b[j].y,TSTP*TN