2017-05-26 2 views
1

Je suis en train de mettre en œuvre le code Matlab pour résoudre l'équation d'onde, ma fonction ressemble à ceci:Équation d'onde avec FDM, Matlab

function [x,t,w] = wave_eqn(xl,xr,yb,yt,M,N,f,l,r,p) 
% input: space interval [xl,xr], time interval [yb,yt] 
% number of space steps M, number of time steps N 
% output: solution w 
D=2;      % diffusion coefficient 
h=(xr-xl)/M; k=(yt-yb)/N; m=M-1; n=N; 
sigma=D*k/(h*h); 
a=diag(1-2*sigma*ones(m,1))+diag(sigma*ones(m-1,1),1); 
a=a+diag(sigma*ones(m-1,1),-1); % define matrix a 
lside=l(yb+(0:n)*k); rside =r(yb+(0:n)*k); 
w(:,1)=f(xl+(1:m)*h)'; % initial conditions 
for j=1:n 
    w(:,j+1)=a*w(:,j)-w(:,j-1)+sigma^2*[lside(j);zeros(m-2,1);rside(j)]; 
end 
w=[lside;w;rside];  % attach boundary conds 
x=(0:m+1)*h;t=(0:n)*k; 
% view(60,30);axis([xl xr yb yt -1 1]) 
end 

%source: numerical analysis 2nd edition 

Je continue à obtenir une erreur à l'équation dans la boucle avec le w (:, j-1) terme qui dit: Les indices d'indice doivent être soit des entiers positifs réels, soit des logiques.

Je ne suis pas sûr de savoir comment résoudre ce problème. Il faut aussi noter que f, p, l, r sont toutes les fonctions d'entrée de x et t. J'ai utilisé un modèle pour l'équation de la chaleur pour faire ce code, mais je ne suis pas sûr de savoir comment mettre en œuvre la quatrième fonction, p. Merci.

Répondre

0

Vous parcourez 1:n, donc j-1=0 et 0 n'est pas un index valide dans Matlab.

Au lieu de cela, boucle de 2:n-1 pour également prendre en compte le terme .

Vous avez déjà déclaré votre condition initiale w(:,1), mais votre méthode numérique nécessite 2 étapes précédentes, vous aurez également besoin d'assigner une condition initiale à w(:,2), peut-être aussi simple que w(:,2) = w(:,1). Puis boucle en utilisant:

for j=2:n-1 
    w(:,j+1)=a*w(:,j)-w(:,j-1)+sigma^2*[lside(j);zeros(m-2,1);rside(j)]; 
end