二维椭圆边值问题的差分格式

二维椭圆边值问题的差分格式一．问题介绍考虑Poisson方程：-Au二f(x,y),(x,y)gG,G是xy平面上一有界区域，其边界r为分段光滑曲线。在r上u满足下列边值条件之一:u|「=a(x,y)(第一边值条件)，Qu=0(x,y)(第二边值条件)，onrQu+ku=y(x,y)(第三边值条件)，on rf(x,y),a(x,y),卩(x,y),y(x,y)及k(x,y)都是连续函数，k>0。本节讨论逼近方程(1)及相应边值条件的差分格式。．区域剖分丄取定沿X轴和y轴方向的步长h和h，h=(h2+h2)2。作两族与坐标轴平行的直线:1212x=ih,i=0,±1,...1y=jh,j=0,±1,.2两族直线的交点(ih,jh)称为网点或节点，记为(x,y)或(i,j)。说两个节点(x,y)和TOC\o"1-5"\h\z1 2 ij ij(x,y)是相邻的，如果•f -fijT j-=1或”一1T j-=1或”一1|+|j—j|=10h2―i i—+h1以G={x,y)ggI表示所有属于G内部的节点集合，并称如此的节点为内点。以r表h ij h示网线x=x或y=y与r的交点集合，并称如此的点为界点。令G=GUr，则Gi j hhh h就是代替域G=GUr的网点集合。若内点(x,y)的四个相邻点都属于G，就称为正则ij h内点；否则称为非正则内点。三．离散格式五点差分格式假定(x,y)为正则内点。沿x，y方向分别用二阶中心差商代替u和u，则有差分ij xxyy

方程：—Auhiju—2u+uu—2u+uh21方程：—Auhiju—2u+uu—2u+uh21h22]=fijij式中寫表示节点GJ上的函数值。九点差分格式u—2u+u利用Taylor展式，将一^3 匸h21化简就得到逼近Poisson方程的九点差分格式:—2u+uj h22j在u处展开，然后相加ij—Au u—2(u+u+u+uhij12ij i—1,j i+j—1 i+1,j i,j+1 i—1,j—1f"(x,y)+h2f〃(x,y)]ij(h2+h2)/h2h21212二f+ 2”一.ij12 1xxij2yy+ui+1,j—1+u+u]i+1,j+1 i—1,j+1四．格式稳定性五点差分格式的收敛阶为O(h2)九点差分格式的收敛阶为O(h4)五．数值例子例1令u(x,y)=sin(兀x)sin(兀y),区间［0,1］。程序结果如下：输入划分区间的点数n(输入0结束程序)：5输入划分区间的点数m(输入0结束程序)：5xiyj准确值u(x,y)近似值u[i][j]误差err[i]0.1666670.1666670.2500000.2557910.0057910.1666670.3333330.4330130.4430430.0100300.1666670.5000000.5000000.5115810.0115810.1666670.6666670.4330130.4430430.0100300.1666670.8333330.2500000.2557910.0057910.3333330.1666670.4330130.4430430.0100300.3333330.3333330.7500000.7673720.0173720.3333330.5000000.8660250.8860850.0200600.3333330.6666670.7500000.7673720.0173720.3333330.8333330.4330130.4430430.0100300.5000000.1666670.5000000.5115810.0115810.5000000.3333330.8660250.8860850.0200600.5000000.5000001.0000001.0231630.023163

0.5000000.6666670.8660250.8860850.0200600.5000000.8333330.5000000.5115810.0115810.6666670.1666670.4330130.4430430.0100300.6666670.3333330.7500000.7673720.0173720.6666670.5000000.8660250.8860850.0200600.6666670.6666670.7500000.7673720.0173720.6666670.8333330.4330130.4430430.0100300.8333330.1666670.2500000.2557910.0057910.8333330.3333330.4330130.4430430.0100300.8333330.5000000.5000000.5115810.0115810.8333330.6666670.4330130.4430430.0100300.8333330.8333330.2500000.2557910.005791误差与步长的2-范数e[0]：0.000805n为5和m为5时其最大误差:0.023163输入划分区间的点数n（输入0结束程序）：10输入划分区间的点数m（输入0结束程序）：10xiyj准确值u（x,y）近似值u[i][j]误差err[i]0.0909090.0909090.0793730.0799150.0005420.0909090.1818180.1523160.1533560.0010400.0909090.2727270.2129190.2143720.0014530.0909090.3636360.2562730.2580220.0017490.0909090.4545450.2788650.2807680.0019030.0909090.5454550.2788650.2807680.0019030.0909090.6363640.2562730.2580220.0017490.0909090.7272730.2129190.2143720.0014530.0909090.8181820.1523160.1533560.0010400.0909090.9090910.0793730.0799150.0005420.1818180.0909090.1523160.1533560.0010400.1818180.1818180.2922920.2942870.0019950.1818180.2727270.4085890.4113780.0027890.1818180.3636360.4917840.4951410.0033560.1818180.4545450.5351380.5387900.0036520.1818180.5454550.5351380.5387900.0036520.1818180.6363640.4917840.4951410.0033560.1818180.7272730.4085890.4113780.0027890.1818180.8181820.2922930.2942870.0019950.1818180.9090910.1523160.1533560.0010400.2727270.0909090.2129190.2143720.0014530.2727270.1818180.4085890.4113780.0027890.2727270.2727270.5711570.5750560.0038980.2727270.3636360.6874540.6921460.004692

0.2727270.4545450.7480570.7531630.0051060.2727270.5454550.7480570.7531630.0051060.2727270.6363640.6874540.6921460.0046920.2727270.7272730.5711570.5750560.0038980.2727270.8181820.4085890.4113780.0027890.2727270.9090910.2129190.2143720.0014530.3636360.0909090.2562730.2580220.0017490.3636360.1818180.4917840.4951410.0033560.3636360.2727270.6874540.6921460.0046920.3636360.3636360.8274300.8330780.0056470.3636360.4545450.9003730.9065180.0061450.3636360.5454550.9003730.9065180.0061450.3636360.6363640.8274300.8330780.0056470.3636360.7272730.6874540.6921460.0046920.3636360.8181820.4917840.4951410.0033560.3636360.9090910.2562730.2580220.0017490.4545450.0909090.2788650.2807680.0019030.4545450.1818180.5351380.5387900.0036520.4545450.2727270.7480570.7531630.0051060.4545450.3636360.9003730.9065180.0061450.4545450.4545450.9797460.9864330.0066870.4545450.5454550.9797460.9864330.0066870.4545450.6363640.9003730.9065180.0061450.4545450.7272730.7480570.7531630.0051060.4545450.8181820.5351380.5387900.0036520.4545450.9090910.2788650.2807680.0019030.5454550.0909090.2788650.2807680.0019030.5454550.1818180.5351380.5387900.0036520.5454550.2727270.7480570.7531630.0051060.5454550.3636360.9003730.9065180.0061450.5454550.4545450.9797460.9864330.0066870.5454550.5454550.9797460.9864330.0066870.5454550.6363640.9003730.9065180.0061450.5454550.7272730.7480570.7531630.0051060.5454550.8181820.5351380.5387900.0036520.5454550.9090910.2788650.2807680.0019030.6363640.0909090.2562730.2580220.0017490.6363640.1818180.4917840.4951410.0033560.6363640.2727270.6874540.6921460.0046920.6363640.3636360.8274300.8330780.0056470.6363640.4545450.9003730.9065180.0061450.6363640.5454550.9003730.9065180.0061450.6363640.6363640.8274300.8330780.0056470.6363640.7272730.6874540.6921460.004692

0.6363640.8181820.4917840.4951410.0033560.6363640.9090910.2562730.2580220.0017490.7272730.0909090.2129190.2143720.0014530.7272730.1818180.4085890.4113780.0027890.7272730.2727270.5711570.5750560.0038980.7272730.3636360.6874540.6921460.0046920.7272730.4545450.7480570.7531630.0051060.7272730.5454550.7480570.7531630.0051060.7272730.6363640.6874540.6921460.0046920.7272730.7272730.5711570.5750560.0038980.7272730.8181820.4085890.4113780.0027890.7272730.9090910.2129190.2143720.0014530.8181820.0909090.1523160.1533560.0010400.8181820.1818180.2922930.2942870.0019950.8181820.2727270.4085890.4113780.0027890.8181820.3636360.4917840.4951410.0033560.8181820.4545450.5351380.5387900.0036520.8181820.5454550.5351380.5387900.0036520.8181820.6363640.4917840.4951410.0033560.8181820.7272730.4085890.4113780.0027890.8181820.8181820.2922930.2942870.0019950.8181820.9090910.1523160.1533560.0010400.9090910.0909090.0793730.0799150.0005420.9090910.1818180.1523160.1533560.0010400.9090910.2727270.2129190.2143720.0014530.9090910.3636360.2562730.2580220.0017490.9090910.4545450.2788650.2807680.0019030.9090910.5454550.2788650.2807680.0019030.9090910.6363640.2562730.2580220.0017490.9090910.7272730.2129190.2143720.0014530.9090910.8181820.1523160.1533560.0010400.9090910.9090910.0793730.0799150.000542误差与步长的2-范数e[1]：0.000128那么log(e[l]/e[0])的值是：1.837767n为10和m为10时其最大误差：0.006687以上数值结果是通过五点差分格式用C++工具求得的;下面采用九点差分格式用C++工具求数值结果。输入划分区间的点数n(输入0结束程序)：5输入划分区间的点数m(输入0结束程序)：5xiyj准确值u(x,y)近似值xiyj准确值u(x,y)近似值u[i][j]误差err[i]0.1666670.1666670.2500000.2496780.0003220.1666670.3333330.4330130.4324550.0005580.1666670.5000000.5000000.4993560.0006440.1666670.6666670.4330130.4324550.0005580.1666670.8333330.2500000.2496780.0003220.3333330.1666670.4330130.4324550.0005580.3333330.3333330.7500000.7490340.0009660.3333330.5000000.8660250.8649100.0011150.3333330.6666670.7500000.7490340.0009660.3333330.8333330.4330130.4324550.0005580.5000000.1666670.5000000.4993560.0006440.5000000.3333330.8660250.8649100.0011150.5000000.5000001.0000000.9987120.0012880.5000000.6666670.8660250.8649100.0011150.5000000.8333330.5000000.4993560.0006440.6666670.1666670.4330130.4324550.0005580.6666670.3333330.7500000.7490340.0009660.6666670.5000000.8660250.8649100.0011150.6666670.6666670.7500000.7490340.0009660.6666670.8333330.4330130.4324550.0005580.8333330.1666670.2500000.2496780.0003220.8333330.3333330.4330130.4324550.0005580.8333330.5000000.5000000.4993560.0006440.8333330.6666670.4330130.4324550.0005580.8333330.8333330.2500000.2496780.000322误差与步长的2-范数e[0]：0.000002n为5和m为5时其最大误差:0.001288输入划分区间的点数n（输入0结束程序）：10输入划分区间的点数m（输入0结束程序）：10xiyj准确值u（x,y）近似值u[i][j]误差err[i]0.0909090.0909090.0793730.0793640.0000090.0909090.1818180.1523160.1522990.0000170.0909090.2727270.2129190.2128950.0000240.0909090.3636360.2562730.2562440.0000290.0909090.4545450.2788650.2788340.0000310.0909090.5454550.2788650.2788340.0000310.0909090.6363640.2562730.2562440.0000290.0909090.7272730.2129190.2128950.0000240.0909090.8181820.1523160.1522990.0000170.0909090.9090910.0793730.0793640.0000090.1818180.0909090.1523160.1522990.000017

0.1818180.1818180.2922920.2922600.0000330.1818180.2727270.4085890.4085430.0000460.1818180.3636360.4917840.4917290.0000550.1818180.4545450.5351380.5350780.0000600.1818180.5454550.5351380.5350780.0000600.1818180.6363640.4917840.4917290.0000550.1818180.7272730.4085890.4085430.0000460.1818180.8181820.2922930.2922600.0000330.1818180.9090910.1523160.1522990.0000170.2727270.0909090.2129190.2128950.0000240.2727270.1818180.4085890.4085430.0000460.2727270.2727270.5711570.5710940.0000640.2727270.3636360.6874540.6873770.0000770.2727270.4545450.7480570.7479730.0000840.2727270.5454550.7480570.7479730.0000840.2727270.6363640.6874540.6873770.0000770.2727270.7272730.5711570.5710940.0000640.2727270.8181820.4085890.4085430.0000460.2727270.9090910.2129190.2128950.0000240.3636360.0909090.2562730.2562440.0000290.3636360.1818180.4917840.4917290.0000550.3636360.2727270.6874540.6873770.0000770.3636360.3636360.8274300.8273380.0000930.3636360.4545450.9003730.9002730.0001010.3636360.5454550.9003730.9002730.0001010.3636360.6363640.8274300.8273380.0000930.3636360.7272730.6874540.6873770.0000770.3636360.8181820.4917840.4917290.0000550.3636360.9090910.2562730.2562440.0000290.4545450.0909090.2788650.2788340.0000310.4545450.1818180.5351380.5350780.0000600.4545450.2727270.7480570.7479730.0000840.4545450.3636360.9003730.9002730.0001010.4545450.4545450.9797460.9796370.0001100.4545450.5454550.9797460.9796370.0001100.4545450.6363640.9003730.9002730.0001010.4545450.7272730.7480570.7479730.0000840.4545450.8181820.5351380.5350780.0000600.4545450.9090910.2788650.2788340.0000310.5454550.0909090.2788650.2788340.0000310.5454550.1818180.5351380.5350780.0000600.5454550.2727270.7480570.7479730.0000840.5454550.3636360.9003730.9002730.0001010.5454550.4545450.9797460.9796370.000110

0.5454550.5454550.9797460.9796370.0001100.5454550.6363640.9003730.9002730.0001010.5454550.7272730.7480570.7479730.0000840.5454550.8181820.5351380.5350780.0000600.5454550.9090910.2788650.2788340.0000310.6363640.0909090.2562730.2562440.0000290.6363640.1818180.4917840.4917290.0000550.6363640.2727270.6874540.6873770.0000770.6363640.3636360.8274300.8273380.0000930.6363640.4545450.9003730.9002730.0001010.6363640.5454550.9003730.9002730.0001010.6363640.6363640.8274300.8273380.0000930.6363640.7272730.6874540.6873770.0000770.6363640.8181820.4917840.4917290.0000550.6363640.9090910.2562730.2562440.0000290.7272730.0909090.2129190.2128950.0000240.7272730.1818180.4085890.4085430.0000460.7272730.2727270.5711570.5710940.0000640.7272730.3636360.6874540.6873770.0000770.7272730.4545450.7480570.7479730.0000840.7272730.5454550.7480570.7479730.0000840.7272730.6363640.6874540.6873770.0000770.7272730.7272730.5711570.5710940.0000640.7272730.8181820.4085890.4085430.0000460.7272730.9090910.2129190.2128950.0000240.8181820.0909090.1523160.1522990.0000170.8181820.1818180.2922930.2922600.0000330.8181820.2727270.4085890.4085430.0000460.8181820.3636360.4917840.4917290.0000550.8181820.4545450.5351380.5350780.0000600.8181820.5454550.5351380.5350780.0000600.8181820.6363640.4917840.4917290.0000550.8181820.7272730.4085890.4085430.0000460.8181820.8181820.2922930.2922600.0000330.8181820.9090910.1523160.1522990.0000170.9090910.0909090.0793730.0793640.0000090.9090910.1818180.1523160.1522990.0000170.9090910.2727270.2129190.2128950.0000240.9090910.3636360.2562730.2562440.0000290.9090910.4545450.2788650.2788340.0000310.9090910.5454550.2788650.2788340.0000310.9090910.6363640.2562730.2562440.0000290.9090910.7272730.2129190.2128950.0000240.9090910.8181820.1523160.1522990.0000170.9090910.9090910.0793730.0793640.000009误差与步长的2-范数e[1]：0.000000那么log(e[l]/e[0])的值是：4.281081n为10和m为10时其最大误差：0.000110附录1#include<stdio.h>#include<stdlib.h>#include<math.h>#include<cmath>#defineMax_N3000doublea[Max_N][Max_N],f[Max_N],det[Max_N][Max_N],v[Max_N][Max_N],err[Max_N],u[Max_N],e[Max_N];intn,m; 〃将区间[A,B]分成n等分，区间[C,D]分成m等分voidGussMethod()//高斯消元法求解五点差分格式；u(x,y)=sin(pi*x)sin(pi*y);{intk,i,j;doubletemp;intt=n*m;for(k=1;k<=t-1;k++){if(a[k][k]==0)exit(0);else{for(i=k+1;i<=t;i++){v[i][k]=a[i][k]/a[k][k];f[i]=f[i]-v[i][k]*f[k];for(j=k;j<=t;j++)a[i][j]=a[i][j]-v[i][k]*a[k][j];}}}u[t]=f[t]/a[t][t];for(i=t-1;i>=1;i--){temp=0;for(j=i+1;j<=t;j++)temp=temp+a[i][j]*u[j];u[i]=(f[i]-temp)/a[i][i];}}intmain(){doubleh,H;//步长doubleA,C;doubleB,D;inti,j;intk=0;printf("输入划分区间的点数n(输入0结束程序)：\n");if(scanf(”％d",&n))printf("输入区间左端点的值A： \n");if(scanf("%lf",&A))printf(”输入区间右端点的值B(输入0结束程序)：\n");if(scanf("%lf",&B))printf("输入划分区间的点数m(输入0结束程序)：\n");if(scanf("%d",&m))printf('输入区间左端点的值C：\n");if(scanf("%lf",&C))printf(”输入区间右端点的值D(输入0结束程序)：\n");while(scanf("%lf",&D)&&n&&D&&B&&m){doublepi=3.1415926;h=(B-A)/(n+1);H=(D-C)/(m+1);for(i=l,j=l;iv=n*m,jv=n*m;i++,j++)〃给系数矩阵赋值{a[i][j]=2/(h*h)+2/(H*H);}for(i=l,j=2;i<=n*m-l,j<=n*m;i++,j++){if(j%m==l)a[i][j]=0;elsea[i][j]=-l.0/(H*H);}for(i=l,j=m+l;i<=(n-l)*m,j<=n*m;i++,j++){a[i][j]=-l/(h*h);}for(i=m+l,j=l;i<=n*m,j<=n*(m-l);i++,j++){a[i][j]=-l/(h*h);}for(i=2,j=l;i<=n*m,j<=n*m-l;i++,j++){if(j%m==0)a[i][j]=0;elsea[i][j]=-l/(H*H);}for(i=1;iv=n;i++)〃给f[i][j]赋值{for(j=1;j<=m;j++){if(i==1&&j==1){f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))+1/(h*h)*sin(pi*(A+(i-1)*h))*sin(pi*(C+j*H))+1/(H*H)*sin(pi*(A+i*h))*sin(pi*(C+(j-1)*H));}elseif(i==1&&j==m){f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))+1/(h*h)*sin(pi*(A+(i-1)*h))*sin(pi*(C+j*H))+1/(H*H)*sin(pi*(A+i*h))*sin(pi*(C+(j+1)*H));}elseif(i==n&&j==1){f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))+1/(h*h)*sin(pi*(A+(i+1)*h))*sin(pi*(C+j*H))+1/(H*H)*sin(pi*(A+i*h))*sin(pi*(C+(j-1)*H));}elseif(i==n&&j==m){f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))+1/(h*h)*sin(pi*(A+(i+1)*h))*sin(pi*(C+j*H))+1/(H*H)*sin(pi*(A+i*h))*sin(pi*(C+(j+1)*H));}elseif(i==1&&j!=1&&j!=m){f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))+1/(h*h)*sin(pi*(A+(i-1)*h))*sin(pi*(C+j*H));}elseif(i==n&&j!=1&&j!=m){f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))+1/(h*h)*sin(pi*(A+(i+1)*h))*sin(pi*(C+j*H));}elseif(i!=1&&i!=n&&j==1){f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))+1/(H*H)*sin(pi*(A+i*h))*sin(pi*(C+(j-1)*H));}elseif(i!=1&&i!=n&&j==m){f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))+1/(H*H)*sin(pi*(A+i*h))*sin(pi*(C+(j+1)*H));}else{f[(i-1)*m+j]=2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H));}}}GussMethod();printf("xi yj准确值u(x,y)近似值u[i][j]误差err[i]\n");doubleT=0;doubletemp=0;doubletem=0;for(i=1;i<=n*m;i++){if(i%m!=0){temp=sin(pi*(A+(i/m+1)*h))*sin(pi*(C+(i%m)*H));//求真实值}else{temp=sin(pi*(A+(i/m)*h))*sin(pi*(C+m*H));}err[i]=u[i]>temp?u[i]-temp:temp-u[i];if(i%m!=0)printf("%lf%lf%lf%lf%lf\n",A+(i/m+1)*h,C+(i%m)*H,temp,u[i],err[i]);elseprintf("%lf%lf%lf%lf%lf\n",A+(i/m)*h,C+m*H,temp,u[i],err[i]);if(err[i]>T)T=err[i];〃求最大误差if(h==H)tem=tem+err[i]*err[i]*h;〃求范数}if(h==H){printf("\n");e[k]=tem;printf("误差与步长的2-范数e[%d]：%lf\n",k,e[k]);if(k>=l)printf("那么log(e[%d]/e[%d])的值是:％lf\n",k,k-l,log(e[k-l]/e[k]));}printf("n为％d和m为％d时其最大误差：%lf\n",n,m,T);for(i=l;i<=n*m;i++)//给系数矩阵清零{for(j=l;j<=n*m;j++){a[i][j]=0;}k++;printf("\n");printf("输入划分区间的点数n(输入0结束程序)：\n");if(scanf(”％d",&n))printf('输入区间左端点的值A： \n");if(scanf(”％lf",&A))printf(”输入区间右端点的值B(输入0结束程序)：\n");if(scanf("%lf",&B))printf(”输入划分区间的点数m(输入0结束程序)：\n");if(scanf("%d",&m))printf('输入区间左端点的值C：\n");if(scanf("%lf",&C))printf(”输入区间右端点的值D(输入0结束程序)：\n");}system("PAUSE");return0;}附录2#include<stdio.h>#include<stdlib.h>#include<math.h>#include<cmath>#defineMax_N3000doublea[Max_N][Max_N],f[Max_N],det[Max_N][Max_N],v[Max_N][Max_N],err[Max_N],u[Max_N],e[Max_N];intn,m; 〃将区间[A,B]分成n等分，区间[C,D]分成m等分voidGussMethod()//高斯消元法求解九点差分格式；u(x,y)=sin(pi*x)sin(pi*y);{intk,i,j;doubletemp;intt=n*m;for(k=1;k<=t-1;k++){if(a[k][k]==0)exit(0);else{for(i=k+1;i<=t;i++){v[i][k]=a[i][k]/a[k][k];f[i]=f[i]-v[i][k]*f[k];for(j=k;j<=t;j++)a[i][j]=a[i][j]-v[i][k]*a[k][j];}}}u[t]=f[t]/a[t][t];for(i=t-1;i>=1;i--)temp=0;for(j=i+1;j<=t;j++)temp=temp+a[i][j]*u[j];u[i]=(f[i]-temp)/a[i][i];}}intmain(){doubleh,H;//步长doubleA,C;doubleB,D;inti,j;intk=0;printf("输入划分区间的点数n(输入0结束程序)：\n");if(scanf(”％d",&n))printf("输入区间左端点的值A： \n");if(scanf("%lf",&A))printf(”输入区间右端点的值B(输入0结束程序)：\n");if(scanf("%lf",&B))printf("输入划分区间的点数m(输入0结束程序)：\n");if(scanf("%d",&m))printf('输入区间左端点的值C： \n");if(scanf("%lf",&C))printf(”输入区间右端点的值D(输入0结束程序)：\n");while(scanf("%lf",&D)&&n&&D&&B&&m){doublepi=3.1415926;h=(B-A)/(n+1);H=(D-C)/(m+1);for(i=l,j=l;iv=n*m,jv=n*m;i++,j++)〃给系数矩阵赋值{a[i][j]=40;}for(i=l,j=2;i<=n*m-l,j<=n*m;i++,j++){if(j%m==l)a[i][j]=0;elsea[i][j]=4-24*h*h/(h*h+H*H);}for(i=l,j=m+l;i<=(n-l)*m,j<=n*m;i++,j++){a[i][j]=4-24*h*h/(h*h+H*H);}for(i=m+l,j=l;i<=n*m,j<=n*(m-l);i++,j++){a[i][j]=4-24*H*H/(h*h+H*H);}for(i=2,j=l;i<=n*m,j<=n*m-l;i++,j++){if(j%m==0)a[i][j]=0;elsea[i][j]=4-24*H*H/(h*h+H*H);}for(i=1,j=m+2;i<=n*m-1,j<=n*m;i++,j++){if(j%m==1)a[i][j]=0;elsea[i][j]=-2;}for(i=2,j=m+1;i<=n*m,j<=n*m-1;i++,j++){if(j%m==0)a[i][j]=0;elsea[i][j]=-2;}for(i=n+1,j=2;i<=n*m-1,j<=n*m;i++,j++){if(j%m==1)a[i][j]=0;elsea[i][j]=-2;}for(i=n+2,j=1;i<=n*m,j<=n*m-1;i++,j++){if(j%m==0)a[i][j]=0;elsea[i][j]=-2;}for(i=l;iv=n;i++)〃给f[i][j]赋值{for(j=l;j<=m;j++){if(i==l&&j==l){f[(i-l)*m+j]=24.0*h*h*H*H/(h*h+H*H)*2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-2.0*h*h*H*H*2.0*pi*pi*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-4.0*sin(pi*(A+(i-l)*h))*sin(pi*(C+j*H))-4.0*sin(pi*(A+i*h))*sin(pi*(C+(j-l)*H))+2.0*sin(pi*(A+(i-l)*h))*sin(pi*(C+(j-l)*H))+2.0*sin(pi*(A+(i+l)*h))*sin(pi*(C+(j-l)*H))+2.0*sin(pi*(A+(i-l)*h))*sin(pi*(C+(j+l)*H));}elseif(i==l&&j==m){f[(i-l)*m+j]=24.0*h*h*H*H/(h*h+H*H)*2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-2.0*h*h*H*H*2.0*pi*pi*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-4.0*sin(pi*(A+(i-1)*h))*sin(pi*(C+j*H))+4.0*sin(pi*(A+i*h))*sin(pi*(C+(j+1)*H))+2.0*sin(pi*(A+(i-1)*h))*sin(pi*(C+(j-1)*H))+2.0*sin(pi*(A+(i+1)*h))*sin(pi*(C+(j+1)*H))+2.0*sin(pi*(A+(i-1)*h))*sin(pi*(C+(j+1)*H));}elseif(i==n&&j==1){f[(i-1)*m+j]=24.0*h*h*H*H/(h*h+H*H)*2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-2.0*h*h*H*H*2.0*pi*pi*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-4.0*sin(pi*(A+(i+1)*h))*sin(pi*(C+j*H))-4.0*sin(pi*(A+i*h))*sin(pi*(C+(j-1)*H))+2.0*sin(pi*(A+(i-1)*h))*sin(pi*(C+(j-1)*H))+2.0*sin(pi*(A+(i+1)*h))*sin(pi*(C+(j-1)*H))+2.0*sin(pi*(A+(i+1)*h))*sin(pi*(C+(j+1)*H));}elseif(i==n&&j==m){f[(i-1)*m+j]=24.0*h*h*H*H/(h*h+H*H)*2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-2.0*h*h*H*H*2.0*pi*pi*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-4.0*sin(pi*(A+(i+1)*h))*sin(pi*(C+j*H))-4.0*sin(pi*(A+i*h))*sin(pi*(C+(j+1)*H))+2.0*sin(pi*(A+(i-1)*h))*sin(pi*(C+(j+1)*H))+2.0*sin(pi*(A+(i+1)*h))*sin(pi*(C+(j-1)*H))+2.0*sin(pi*(A+(i+1)*h))*sin(pi*(C+(j+1)*H));}elseif(i==1&&j!=1&&j!=m){f[(i-1)*m+j]=24.0*h*h*H*H/(h*h+H*H)*2.0*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-2.0*h*h*H*H*2.0*pi*pi*pi*pi*sin(pi*(A+i*h))*sin(pi*(C+j*H))-4.0*sin(pi*(A+(i-1)*h))*sin(pi*(C+j*H))+2.0*sin(pi*(A+(i-1)*h))*sin(pi*(C+(j-1)*H))+2.0*sin(pi*(A+(i-1)*h))*sin(pi*(