(x,y,z)需要解非线性方程组和作二元分片插值。
-----------------------------------------------------------------------*/
void main()
{
double u_0[6]={0,0.4,0.8,1.2,1.6,2};
double t_0[6]={0,0.2,0.4,0.6,0.8,1.0};
function( t_1, u_1, x_0);
chazhi(t_1,u_1,z,t_0,u_0,z_0);
k=nihe(z,c,delt);
for(i=0;i<=10;i++)
for(j=0;j<=20;j++)
{
printf("\n");printf("\t");
printf("xi=%3.2f
double b=0; int i,flag=k; b=fabs(x[k]); for(i=k+1;i<=3;i++)
if(fabs(x[i])>b) {
flag=i; b=fabs(x[i]); } return flag; } /*----------------------------------------------------------------------fabsmax_flag_1 函数：比较数据绝对值的大小，并返回按模最大值所在的行号 -----------------------------------------------------------------------*/ int fabsmax_flag_1(double x[11],int k,int l) { double b=0; int i,flag=k; b=fabs(x[k]); for(i=k+1;i<=l;i++) if(fabs(x[i])>b) { flag=i; b=fabs(x[i]); } return flag; } /*----------------------------------------------------------------------gauss 函 数 ： 列 主 元 素 Gauss 法 解 线 性 方 程 组 -----------------------------------------------------------------------*/ void gauss(double a_0[4][4],double b[4],double deltx[4]) { double a[4][4],x[4],aa; int i,j,k,flag; for(i=0;i<=3;i++) for(j=0;j<=3;j++) a[i][j]=a_0[i][j];
计算时的 k 值需要程序自动确定，要求最小的 k 值使精度达到：

10
=
20

(
f
(xi ,
yi )

p(xi ,
yi ))2
107
。
i=0 j 0
二、计算结果
1 ． 二 元 二 次 分 片 插 值 得 到 数 表 ： xi ， y j ， f xi, y j ，
i 0,1,,10; j 0,1,, 20 。
n
fn (x(k))
其中 ej 是第 j 个 n 维基本单位向量。迭代公式
x(k1) x(k ) J (x(k ) , h(k ) )1 F (x(k ) ), k 0,1,…… （2）二元二次分片插值
根据数表作二元二次分片插值得到函数 z g t,u 。分片插值能保证收
敛性，且有较高的精度。根据给定的数表，可将整个插值区域分成 16 个小
u_0[6],double t_0[6],double z_0[6][6],double c[11][11]);
/*-----------------------------------------------------------------------
main 函数：在给定区域上作二元拟合得到满足精度要求的 P(X,Y)，为确定拟合节点


f1(x1, x2 ,…,xn ) 0
F
(
x)

0
，一般形式为

f2
(
x1
,
x2
,…,xn
)

0
的非线性方程组，
fn (x1, x2 ,…,xn ) 0
令 h(k ) (h1(k ) , h2(k ) ,, hn(k ) ), hj(k ) 0, j 1, 2,n
其中：
1
B

1 1
x0 x1 xn

x0k x1k
xnk

,G

1 1 1
y0 y1 yn

y0k y1k
ynk



,U


(
f
(xi ,
y j ))nm
令 A (BT B)1 BTU ， DT G(GT G)1 ， 作 变 换 得 BT BA BTU ， GTGD GT ，可通过列主元素 Gauss 法解得矩阵 A 和 D，再算得 C ADT 。
yi=%3.2f
f(xi,yi)=%14.12e",0.08*i,(0.5+0.05*j),z[i][j]);
} printf("\n"); for(i=1;i<=5;i++) {
printf("\t"); printf("k=%d delt=%14.12e",i,delt[i-1]); printf("\n"); } printf("\n");printf("\t"); printf("k=5 delt=%14.12e",delt[4]); printf("\n");printf("\n"); for(i=0;i<=k;i++) for(j=0;j<=k;j++) {
《数值分析》大作业（3）
一、 算法设计方案
1．算法的总体设计 根据题目中给定的关于 z、t、u 的数表及方程组，可以得到关于 z、x、
y 的数表，以此数表作为基础，可以得到 Z(x，y)的表达式。
要 在 区 域 D x, y 0 x 0.8, 0.5 y 1.5 上 作 二 元 拟 合 函 数
p x, y
k

crs xr y s
，并使
p x,
y
在某一最小的
k
值达到如下精度，
r ,s 0
10 20

f
xi, y j p xi, y j
2
107 ，其中 xi 0.08i ， y j 0.5 0.05 j 。
i0 j0
拟合节点的 xi ， y j 值已知，要求出相应的 zij g tij,uij f xi, y j 。先
数 {r
(x)
s ( y)(r

0...k; s

0...k )}
，其中 r
(x)

xr ,
s (y)

ys
。拟合函数为
p(x,
y)

k

k

crsr
(x)
s
(
y)
，其中
C

crs
为要求的系数矩阵。实际上就是先
s0 r0
求解系数 C 的矩阵：
C (BT B)1 BTUG(GT G)1
t_0[6],double u_0[6],double z_0[6][6]);
void f(double t_1,double u_1,int flag[2]);
int nihe(double u[11][21],double c[11][11],double delt[11]);
void gauss_1(double b[11][11],double u[11][21],double a[11][21],int l);
double
z_0[6][6]={{-0.5,-0.34,0.14,0.94,2.06,3.5},{-0.42,-0.5,-0.26,0.3,1.18,2.3
8},{-0.18,-0.5,-0.5,-0.18,0.46,1.42},
{0.22,-0.34,-0.58,-0.5,-0.1,0.62},{0.78,-0.02,-0.5,-0.66,-0.5,-0.02},
2．二元拟合选择过程的 k ， 值。
3 ． 达 到 精 度 要 求 的 k ， 的 值 及 p x, y 中 的 系 数 crs r 0,1,, k; s 0,1,, k 。
4．数表：
x
* i
，
y j* ，
f
x
* i
,
y
* jBiblioteka ，pxi*, 解非线性方程组由一组 xi, y j 可得到一组 t ,ij uij ，再根据所给的数表作二
元分片插值得到函数 z g t,u ，最后将 t ,ij uij 代入函数 z g t,u 得到
zij g t ij, uij f xi, y j 。
y
* j