I want to run this code?
显示 更早的评论
shape = 2.5;
N = 30;
Np = N^2;
[X,Y] = meshgrid(linespace(0,1,N),linespace(0,1,N));
x = reshape(X,N^2,1);
y = reshape(Y,N^2,1);
f = -2*(2*y.^3-3*y.^2+1)+6*(1-x.^2).*(2*y-1);
u = (1-x.^2).*(2*y.^3-3*y.^2+1); % exact solution
H = zeros(Np,Np);
rx = zeros(Np,Np);
ry = zeros(Np,Np);
r = zeros(Np,Np);
dirichletBCs = find(x==0 | x==1); % identify boundry and interior centers
neumannBCs = find ((y==0 | y==1)& ~(x==0 | x==1));
interior = find('x~=0' & 'x~ = 1' & 'y~ = 0' & 'y~ = 1' );
f (dirichletBCs)=u(dirichletBCs);
f (neumannBCs)= 0;
for i =1:Np
for j =1:Np
rx ( i , j ) = x ( i )-x ( j ) ;
ry ( i , j ) = y ( i )-y ( j ) ;
r ( i , j ) = sqrt ( rx ( i , j )^2 + ry ( i , j )^2 );
end
end
H(dirichletBCs,:)=mq(r(dirichletBCs,:),shape); % e v a l u a t i on matrix
H(neumannBCs,:)= mqDerivatives(r(neumannBCs,:),ry(neumannBCs,:),shape,1);
Hxx = mqDerivatives(r(interior,:),rx(interior,:),shape,2);
Hyy = mqDerivatives(r(interior,:),ry(interior,:),shape,2);
H(interior,:)=Hxx+Hyy;
kappa=cond(H);
alpha=H\f;
B = mq( r , shape ) ; % system matrix
uh = B*alpha;
error = norm(uh-u,inf);
t = delaunay(x,y);trisurf(t,x,y,abs(uh-u));
xlabel 'x' , ylabel 'y' colormap('Summer')
回答（1 个）
DGM
2022-3-9
I can't test anything, since there are missing functions, but:
% there is no linespace()
%[X,Y] = meshgrid(linespace(0,1,N),linespace(0,1,N));
[X,Y] = meshgrid(linspace(0,1,N),linspace(0,1,N));
% this is an inexplicable misuse of char vectors and invalid operator spacing
%interior = find('x~=0' & 'x~ = 1' & 'y~ = 0' & 'y~ = 1' );
interior = find(x~=0 & x~=1 & y~=0 & y~=1 );
% using command syntax for xlabel,ylabel is technically valid
% but the line continuation will cause an error at colormap()
%xlabel 'x' , ylabel 'y' colormap('Summer')
xlabel('x')
ylabel('y')
colormap('Summer')
3 个评论
Akhtar Jan
2022-3-9
DGM
2022-3-9
Unless you're declaring that you've given up, you'll have to explain what exactly you did and what the observed outcomes were. If there are errors, show them. Do you have the functions mq() and mqDerivatives()? I don't.
You are the only person who has access to the information necessary to troubleshoot the problem.
Akhtar Jan
2022-3-9
类别
在 帮助中心 和 File Exchange 中查找有关 Graphics Objects 的更多信息
2022-3-9
2022-3-9
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
