how can i run this code need help?

%its about Solving PDEs by Asymmetric MQ Collocation the above code is already given in the book which %i have
%pasted below i havn't much experience in MATALAB thats why i shared it here...
%Neumann Boundary Conditions
%The function f is set to f(x, y) = -2(2y3 - 3y2 + 1) + 6(1 - x2)(2y - 1) and Dirichlet boundary conditions
%of u(0, y) = 2y2 - 3y2 + 1 and u(1, y) = 0 are applied as well as Neumann boundary conditions of
%∂u/∂y = 0, along y = 0 and y = 1.
%A shape parameter of ε = 2.5 results in an evaluation matrix to be inverted that has a very high condition %number
%κ(H) ≈ 4.6e19.
clc
Nb = 22; % number of boundary point s
Np = 60; % t o t a l number of c ent e r , i n t e r i o r and boundary
centers = dlmread(centersUnitCircle60.txt ,'_');
x = centers(:,1);y = centers(:,2);
u = 65./(65+(x-0.2).^2+(y+0.1).^2 );% e x ac t s o l u t i o n
f = 130./(65+(x-0.2).^2+(y+0.1).^2 ).^3.*(2.*x-0.4).^2 - ...
260./(65+(x-0.2).^2+(y+0.1).^2).^2+...
130./(65+(x-0.2).^2+(y+0.1).^2).^3.*(2.*y+0.2).^2;
H = zeros (Np,Np ); rx = zeros (Np,Np ); ry = zeros (Np,Np); r = zeros (Np,Np );
f ( 1 :Nb) = u ( 1 :Nb ) ; % Di r i c h l e t Boundary c ond i t i ons
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
index = 1 ;
for shape =1.0: -0.01: 0.1;
H( 1 :Nb , : ) = mq( r ( 1 :Nb , : ) , shape ) ; % enf or c e boudary c ond i t i ons
Hxx = mqDerivatives( r (Nb+1:Np , : ) , rx (Nb+1:Np , : ) , shape , 2 ) ; % enf or c e PDE
Hyy = mqDerivatives( r (Nb+1:Np , : ) , ry (Nb+1:Np , : ) , shape , 2 ) ;
H(Nb+1:Np , : ) = Hxx + Hyy ; % e v a l u a t i on matrix
kappa ( index ) = cond(H) ;
warning off
lambda = H\ f ; % expansion c o e f f i c i e n t s v i a Gaussian e l imi na t i on
warning on
B = mq( r , shape ) ; % system matrix
uh = B*lambda ; % approximate s o l u t i o n
e r ( index ) = norm(u?uh , i n f ) ;
sh ( index ) = shape ;
index = index + 1 ;
end
semilogy( sh , er , 'b' ) , xlabel ’ shape parameter ’ , ylabel ’max e r r o r ’
figure , semilogy( sh , kappa ) , xlabel ’ shape parameter ’ , ylabel ’ \kappa (H) ’

댓글 수: 14

Walter Roberson
Walter Roberson 2022년 3월 9일
what problems are you having? We do not have your file so we cannot test your code.
it ita can not run
its a poissonNeumannQM code
shaficode
Unrecognized function or variable 'linespace'.

Error in shaficode (line 4)
[X,Y] = meshgrid(linespace(0,1,N),linespace(0,1,N));
it ita can not run
its a poissonNeumannQM code
neman
Unable to resolve the name centersUnitCircle60.txt.

Error in neman (line 4)
centers = dlmread(centersUnitCircle60.txt ,'_');
This line:
centers = dlmread(centersUnitCircle60.txt ,'_');
should be:
centers = dlmread('centersUnitCircle60.txt');
Akhtar Jan
Akhtar Jan 2022년 3월 9일
still not run
Walter Roberson
Walter Roberson 2022년 3월 9일
We do not have your txt file, so we cannot run the code to test it. You will need to show us the error message instead of just telling us that it does not run.
Jan
Jan 2022년 3월 9일
편집: Jan 2022년 3월 9일
@Akhtar Jan: Please do do explain, what is not happening, but explain, what you observe. "still not run" does not explain anything. You do have a useful error message on the screen already, so share it with the readers.
The code is full of typos. Did you copy&paste it from a PDF? E.g.
e r ( index ) = norm(u?uh , i n f ) ;
% should be:
er(index) = norm(uh , inf);
Akhtar Jan
Akhtar Jan 2022년 3월 9일
Pls check this
Check what? It is your turn either to post some code and explain clearly, which problem you have with it, or to explain clearly, what you want to achieve and showing, what you have tried so far.
Currently we see a massively distorted text, which looks almost like Matlab, but the readers have to guess many details. Lines like:
e r ( index ) = norm(u?uh , i n f ) ;
are pure nonsense and you cannot expect, that readers in the internet can magically guess, what its intention is.
Walter Roberson
Walter Roberson 2022년 3월 9일
편집: Walter Roberson 2022년 3월 9일
clc
Nb = 22; % number of boundary point s
Np = 60; % t o t a l number of c ent e r , i n t e r i o r and boundary
centers = dlmread('centersUnitCircle60.txt' ,'_');
x = centers(:,1);y = centers(:,2);
u = 65./(65+(x-0.2).^2+(y+0.1).^2 );% e x ac t s o l u t i o n
f = 130./(65+(x-0.2).^2+(y+0.1).^2 ).^3.*(2.*x-0.4).^2 - ...
260./(65+(x-0.2).^2+(y+0.1).^2).^2+...
130./(65+(x-0.2).^2+(y+0.1).^2).^3.*(2.*y+0.2).^2;
H = zeros (Np,Np ); rx = zeros (Np,Np ); ry = zeros (Np,Np); r = zeros (Np,Np );
f ( 1 :Nb) = u ( 1 :Nb ) ; % Di r i c h l e t Boundary c ond i t i ons
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
index = 1 ;
for shape =1.0: -0.01: 0.1;
H( 1 :Nb , : ) = mq( r ( 1 :Nb , : ) , shape ) ; % enf or c e boudary c ond i t i ons
Hxx = mqDerivatives( r (Nb+1:Np , : ) , rx (Nb+1:Np , : ) , shape , 2 ) ; % enf or c e PDE
Hyy = mqDerivatives( r (Nb+1:Np , : ) , ry (Nb+1:Np , : ) , shape , 2 ) ;
H(Nb+1:Np , : ) = Hxx + Hyy ; % e v a l u a t i on matrix
kappa ( index ) = cond(H) ;
warning off
lambda = H\f ; % expansion c o e f f i c i e n t s v i a Gaussian e l imi na t i on
warning on
B = mq( r , shape ) ; % system matrix
uh = B*lambda ; % approximate s o l u t i o n
er( index ) = norm(u - uh , inf ) ;
sh( index ) = shape ;
index = index + 1 ;
end
semilogy( sh , er , 'b' ) , xlabel 'shape parameter' , ylabel 'max error'
figure , semilogy( sh , kappa ) , xlabel 'shape parameter', ylabel '\kappa (H)'
Walter Roberson
Walter Roberson 2022년 3월 9일
편집: Walter Roberson 2022년 3월 9일
shape = 2.5;
N = 30;
Np = N^2;
[X,Y] = meshgrid(linspace(0,1,N),linspace(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')
Walter Roberson
Walter Roberson 2022년 3월 9일
Neither of those are tested, as we do not have the data file and we do not have functions mq or mqDerivative
Akhtar Jan
Akhtar Jan 2022년 3월 10일
편집: Akhtar Jan 2022년 3월 10일
%its about Solving PDEs by Asymmetric MQ Collocation the above code is already given in the book which %i have
%pasted below i havn't much experience in MATALAB thats why i shared it here...
%Neumann Boundary Conditions
%The function f is set to f(x, y) = -2(2y3 - 3y2 + 1) + 6(1 - x2)(2y - 1) and Dirichlet boundary conditions
%of u(0, y) = 2y2 - 3y2 + 1 and u(1, y) = 0 are applied as well as Neumann boundary conditions of
%∂u/∂y = 0, along y = 0 and y = 1.
%A shape parameter of ε = 2.5 results in an evaluation matrix to be inverted that has a very high condition %number
%κ(H) ≈ 4.6e19.
Jan
Jan 2022년 3월 10일
@Walter Roberson: Although the OP does not want to do it, I want to say: Thank you for spending your time for fixing the code.

댓글을 달려면 로그인하십시오.

답변 (1개)

Image Analyst
Image Analyst 2022년 3월 9일

0 개 추천

The function is called 'linspace', not 'linespace'.

댓글 수: 4

And many of your spaces are not allowed, like
close a l l
and
0 . 2
and
i n f
etc.
Akhtar Jan
Akhtar Jan 2022년 3월 10일
i have tried it without spaces i mean after modifications i tried but not.....run
Akhtar Jan
Akhtar Jan 2022년 3월 10일
%its about Solving PDEs by Asymmetric MQ Collocation the above code is already given in the book which %i have
%pasted below i havn't much experience in MATALAB thats why i shared it here...
%Neumann Boundary Conditions
%The function f is set to f(x, y) = -2(2y3 - 3y2 + 1) + 6(1 - x2)(2y - 1) and Dirichlet boundary conditions
%of u(0, y) = 2y2 - 3y2 + 1 and u(1, y) = 0 are applied as well as Neumann boundary conditions of
%∂u/∂y = 0, along y = 0 and y = 1.
%A shape parameter of ε = 2.5 results in an evaluation matrix to be inverted that has a very high condition %number
%κ(H) ≈ 4.6e19.
Jan
Jan 2022년 3월 10일
편집: Jan 2022년 3월 10일
@Akhtar Jan: Your code contains to many spaces. Your description contains a lot of % characters. It is a good idea to spend the time to post clean messages, if you want to motivate others to help you.
You have been told already, that statements as "but not.....run" are not useful. Post the error message you see, or even better: read the message and try to fix the specific problem.

댓글을 달려면 로그인하십시오.

카테고리

도움말 센터File Exchange에서 Partial Differential Equation Toolbox에 대해 자세히 알아보기

제품

릴리스

R2015a

태그

질문:

2022년 3월 9일

댓글:

Jan
2022년 3월 10일

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by