How to use iterations in solving system of equations? want to update the vector Y at each iteration but ending up with same Y defined in the beginning.
1 visualizzazione (ultimi 30 giorni)
Mostra commenti meno recenti
Tanya Sharma
il 2 Feb 2021
Commentato: Alan Stevens
il 3 Feb 2021
Hi,
I am runnig the code for an iterative solution.
The problem is that I want to update the vector Y at each iteration , and solve the system ynew=A\R with new values of Y .
But in my case I am unable to update the Y after the iteration and ending up using the same Y as defined in the beginning of the code.
I have attached the cheb function too.
clear;
clc;
%%%%----------diff matrix---------------%%
N=5;
[D,x]=cheb(N);
D1=2*D; %%derivative matrix
D2=D1*D1;
D3=D1*D1*D1;
I=eye(N+1);
%%%%---parameters-----------------%%%
Re=5;
alpha=-1;
eta= (1+x)/2;
%-------------initial guess--------%%%
change =1;
it=0;
while change >1e-15
Y = (-eta.^2+1); %%this is the initial approximation for Y %%----NEED to UPDATE AT EACH ITERATION LEVEL
dY = D1*Y; %%derivative for Y
a0 = diag(2*alpha*Re*Y+4*alpha*alpha); %%%----a0 and a1 contain Y and dY
a1 = diag(2*alpha*Re*dY);
A = D3 + a0.*D1 + a1.*I;
A(1,:)=0;
A(1,1)=1;
A(N+1,:)=0;
A(N+1,N+1)=1;
A(N,:)=D(N+1,:);
R = zeros(N+1,1)+2*alpha*Re*Y.*dY; %%----R contains Y and dY
R(1)=0;
R(N+1)=1;
R(N)=0;
ynew = A\R;
change = norm(ynew-Y,inf);
Y = ynew; %%updated Y, it must update the Y so that all the calculation can run with new Y at each iteration
it = it+1;
end
figure(1)
plot(eta,ynew,'r-');
hold on
0 Commenti
Risposta accettata
Alan Stevens
il 2 Feb 2021
Try taking
Y = (-eta.^2+1);
outside the while command.. i.e. have
Y = (-eta.^2+1);
while change >1e-15
...etc
4 Commenti
Alan Stevens
il 3 Feb 2021
Ok. It works just fine (though I suggest you reduce your change test to, say, 10^-8, as this makes little difference to the accuracy, but speeds things up enormously as you increase the size of N). However, I can't comment on the validity of the results as I know nothing of the background to the problem you are trying to solve.
Più risposte (0)
Vedere anche
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!