%βֵ
len=length(x);
alpha=0.3*ones(len,len);
N=length(x);
M=length(y);
xx=0;
yy=0;
fzz=0;
xp=NaN;
yp=NaN;
zp=NaN;
s=1;
s1=1;
s2=1;
s3=0;
s4=0;
for i=1:N-1
 for j=1:N
    xx=x(i)+(x(i+1)-x(i))*(x(j)-x(1))/(x(N)-x(1));
 xp(s1)=xx;sl=s1+1;
 end
end
for i=1:M-1
 for j=1:M
    yy=y(i)+(y(i+1)-y(i))*(y(j)-y(1))/(y(N)-y(1));
 yp(s2)=yy;s2=s2+1;
 end
end
KK=[];
for i=2:N
for j=2:M
c=(z(i-1,j-1)-z(i-1,j)-z(i,j-1)+z(i,j)-alpha(i,j)*(z(1,1)-z(N,1)-z(1,M)+z(N,M)))/(x(1)*y(1)-x(N)*y(1)-x(1)*y(M)+x(N)*y(M));
b=(z(i-1,j-1)-z(i,j-1)-alpha(i,j)*(z(1,1)-z(N,1))-c*(x(1)*y(1)-x(N)*y(1)))/(x(1)-x(N));
d=(z(i-1,j-1)-z(i-1,j)-alpha(i,j)*(z(1,1)-z(1,M))-c*(x(1)*y(1)-X(1)*y(M)))/(y(1)-y(M));
k=z(i,j)-b*x(N)-d*y(M)-alpha(i,j)*z(N,M)-c*x(N)*y(M);
[cbdk];
i;
j;
for t1=1:N
for t2=l:M
KK=[KK;[(i-2)*M+tl,(j-2)*N+t2]];  
zz((i-2)*N+tl,(j-2)*M+t2)=b*x(t1)+d*y(t2)+c*x(t1)*y(t2)+alpha(i,j)*z(tl,t2)+k;  
end  
end  
end  
end  
surf(xp,yp,zz)  
x=xp;  
y=yp;  
z=zz;