UVa 199 Partial Differential Equations
题目分析
本题要求通过有限差分法（Finite Difference Method\texttt{Finite Difference Method}Finite Difference Method）求解二维泊松方程（Poisson Equation\texttt{Poisson Equation}Poisson Equation）的边界值问题。给定一个[0,1]×[0,1][0, 1] \times [0, 1][0,1]×[0,1]的正方形区域，将其离散化为(n+1)×(n+1)(n+1) \times (n+1)(n+1)×(n+1)个网格点，步长h=1nh = \frac{1}{n}h=n1。边界上的函数值由b1,b2,b3,b4b_1, b_2, b_3, b_4b1,b2,b3,b4给出，内部点的值未知。
核心公式是将泊松方程∂2u∂x2+∂2u∂y2=f(x,y)\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} = f(x, y)∂x2∂2u+∂y2∂2u=f(x,y)离散化为：
1h2(ui−1,j+ui,j−1+ui+1,j+ui,j+1−4ui,j)=fi,j \frac{1}{h^2}(u_{i-1,j} + u_{i,j-1} + u_{i+1,j} + u_{i,j+1} - 4u_{i,j}) = f_{i,j}h21(ui−1,j+ui,j−1+ui+1,j+ui,j+1−4ui,j)=fi,j
更一般地，题目引入两个3×33 \times 33×3矩阵vvv和ggg，将方程扩展为：
1h2∑p=−11∑q=−11vp+1,q+1⋅ui+p,j+q=∑p=−11∑q=−11gp+1,q+1⋅fi+p,j+q \frac{1}{h^2} \sum_{p=-1}^{1} \sum_{q=-1}^{1} v_{p+1,q+1} \cdot u_{i+p, j+q} = \sum_{p=-1}^{1} \sum_{q=-1}^{1} g_{p+1,q+1} \cdot f_{i+p, j+q}h21p=−1∑1q=−1∑1vp+1,q+1⋅ui+p,j+q=p=−1∑1q=−1∑1gp+1,q+1⋅fi+p,j+q
对于边界附近的点，当索引超出[0,n][0, n][0,n]范围时，使用给定的边界条件替换。
解题思路
离散化：内部点共有(n−1)2(n-1)^2(n−1)2个，按行优先顺序编号（从左到右、从上到下）。
矩阵构建：对于每个内部点(i,j)(i, j)(i,j)（其中1≤i,j≤n−11 \le i, j \le n-11≤i,j≤n−1对应离散化后的索引，但注意代码中使用i,ji, ji,j从222到nnn表示内部点），遍历3×33 \times 33×3邻域中的999个点。如果邻域点在边界上，则使用边界条件值并移到方程右侧；如果邻域点在内部，则将其系数加入系数矩阵的相应列。
系数处理：由于1h2=n2\frac{1}{h^2} = n^2h21=n2是整数，所有计算均可使用整数运算。vvv矩阵的系数乘以n2n^2n2后作为系数矩阵的元素。
右侧向量：bbb向量由两部分组成：
- 边界条件贡献：−左值×n2-\text{左值} \times n^2−左值×n2
- fff函数贡献：gp,q×fi+p,j+qg_{p,q} \times f_{i+p, j+q}gp,q×fi+p,j+q
注意事项
- 输入中函数fff的矩阵表示：每行对应固定的yyy值（从y=0y=0y=0到y=1y=1y=1），每列对应固定的xxx值（从x=0x=0x=0到x=1x=1x=1）。
- 边界条件b3b_3b3和b4b_4b4的索引需注意对应关系：b3(y)b_3(y)b3(y)对应左边界x=0x=0x=0，b4(y)b_4(y)b4(y)对应右边界x=1x=1x=1。
- 输出格式要求：矩阵aaa按行输出(n−1)2(n-1)^2(n−1)2行，每行(n−1)2(n-1)^2(n−1)2个整数；向量bbb输出一行(n−1)2(n-1)^2(n−1)2个整数。
代码实现
// Partial Differential Equations// UVa ID: 199// Verdict: Accepted// Submission Date: 2016-03-20// UVa Run Time: 0.000s//// 版权所有（C）2016，邱秋。metaphysis # yeah dot net#include<bits/stdc++.h>usingnamespacestd;inta[35][35],b[1000],v[3][3],g[3][3],boundary[4][35],f[35][35],temp[1000],n;intmain(){intm;cin>>m;while(m--){cin>>n;// 读取 3x3 矩阵 vfor(inti=0;i<3;i++)for(intj=0;j<3;j++)cin>>v[i][j];// 读取 3x3 矩阵 gfor(inti=0;i<3;i++)for(intj=0;j<3;j++)cin>>g[i][j];// 读取四个边界条件：b1(上边界), b2(下边界), b3(左边界), b4(右边界)for(inti=0;i<4;i++)for(intj=0;j<=n;j++)cin>>boundary[i][j];// 读取函数 f 在网格点上的值for(inti=0;i<=n;i++)for(intj=0;j<=n;j++)cin>>f[i][j];intcounter=0;// 遍历所有内部点 (i, j)，i 和 j 从 2 到 n 对应内部点的实际索引// 内部点在离散化区域中对应 (1..n-1) 的索引，但这里 i,j 从 2 开始是为了方便处理 3x3 邻域for(inti=2;i<=n;i++)for(intj=2;j<=n;j++){// 临时数组存储系数矩阵的一行，长度为 (n-1)^2fill(temp,temp+1000,0);intright=0;// 遍历当前点 (i, j) 周围的 3x3 邻域for(intx=i-2;x<=i;x++)for(inty=j-2;y<=j;y++){intleft=0;// 根据邻域点的位置判断处理方式if(x==0)// 上边界（y 方向，对应 b2？需注意索引）{left=v[x+2-i][y+2-j]*boundary[1][y];}elseif(x==n)// 下边界（对应 b1）{left=v[x+2-i][y+2-j]*boundary[0][y];}elseif(y==0)// 左边界（对应 b3）{left=v[x+2-i][y+2-j]*boundary[2][n-x];}elseif(y==n)// 右边界（对应 b4）{left=v[x+2-i][y+2-j]*boundary[3][n-x];}else// 内部点，记录系数{// 内部点的线性索引：(x-1)*(n-1) + (y-1)temp[(x-1)*(n-1)+y-1]=v[x+2-i][y+2-j]*n*n;}// 右侧项：减去边界贡献 * n^2，加上 g * f 的贡献right+=-left*n*n+g[x+2-i][y+2-j]*f[x][y];}b[counter++]=right;// 输出系数矩阵的当前行for(intk=0;k<(n-1)*(n-1);k++)cout<<(k>0?" ":"")<<temp[k];cout<<endl;}// 输出右侧向量 bfor(inti=0;i<(n-1)*(n-1);i++)cout<<(i>0?" ":"")<<b[i];cout<<endl;}return0;}复杂度分析
- 时间复杂度：O(n4)O(n^4)O(n4)，因为对于每个内部点（共(n−1)2(n-1)^2(n−1)2个），需要遍历999个邻域点并更新系数矩阵的一行。由于n≤30n \le 30n≤30，完全可接受。
- 空间复杂度：O(n2)O(n^2)O(n2)，用于存储边界条件、函数fff和临时数组。