Calculating second centred difference for FDM - indices for iteration
7 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
I have an NxN matrix 'phi' that is a function of x and y, and is constantly getting updated over iteration (is a function of t as well). To calculate the Laplacian (sum of spatial second derivatives), I calculate the second centred difference, defined as :
I am attaching the codes I have been given in two cases: first when phi is just a function of x (1D) and second when it is a function of x and y (2D). In both cases, my doubt is why pind and nind are defined this way (pind and nind refer to the positive index - forward step in x or y dirn - and negative index - backward step - respectively). I believe it is just a question of me not understanding the syntax used.
pind=[2:N,N];
nind=[1,1:N-1];
for i=1:length(c)
dphidx = (phi(pind)-phi(nind))/(2*h);
d2phi = (phi(pind)+phi(nind)-2*phi)/(h^2);
end
code i) 1 dimensional example
pind=[2:N,1];
nind=[N,1:N-1];
%iteration begins
for iter = 1:Niter,
dphidx = (phi(:,pind)-phi(:,nind))/(2*h);
dphidy = (phi(pind,:)-phi(nind,:))/(2*h);
d2phi = (phi(:,pind)+phi(:,nind)+phi(pind,:)+phi(nind,:)-4*phi)/(h^2);
end
code ii) 2 dimensional example
I have checked test outputs to see how pind and nind change over the iterations, but I still don't understand why the definition is done like this. Can attach screenshots of that test output if required.
0 Commenti
Risposte (1)
darova
il 14 Giu 2020
i think something is wrong here
pind = 2 3 4 5
nind = 1 1 2 3
correct way for 1D would be:
for i = 2:length(phi)-1
dphidx = (phi(i-1)-phi(i+1))/2/h;
d2phidx = (phi(i-1)+2*phi(i)-phi(i+1))/h^2;
end
2d:
for i = 2:size(phi,1)-1
for j = 2:size(phi,2)-1
dphidx = (phi(i-1,j)-phi(i+1,j))/2/h; % dphi/dx
dphidy = (phi(i,j-1)-phi(i,j+1))/2/h; % dphi/dy
jj = j-1:j+1;
dphidx1 = (phi(i-1,jj)-phi(i+1,jj))/2/h;
d2phi = (dphidx1(3) - dphidx1(1))/2/h; % d2phi/dx/dy
end
end
2 Commenti
darova
il 16 Giu 2020
Here is the formula
dphidx1 is . It has length of 3 elements.
Here is a simple scheme:
Vedere anche
Categorie
Scopri di più su Directed Graphs in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!