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 visualización (últimos 30 días)
Mostrar comentarios más antiguos
Tanya Sharma
el 2 de Feb. de 2021
Comentada: Alan Stevens
el 3 de Feb. de 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 comentarios
Respuesta aceptada
Alan Stevens
el 2 de Feb. de 2021
Try taking
Y = (-eta.^2+1);
outside the while command.. i.e. have
Y = (-eta.^2+1);
while change >1e-15
...etc
4 comentarios
Alan Stevens
el 3 de Feb. de 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.
Más respuestas (0)
Ver también
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!