Chebyshev smoother, I cannot get this to run, I dont understand why.
1 vue (au cours des 30 derniers jours)
Afficher commentaires plus anciens
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1363733/image.jpeg)
dtVx = rho*resolx.^2./(4*eta_nbv*(1.+eta_b)*4.1); %dtVx is a very large vector
A = (4*eta_nbv*(1.+eta_b)*4.1); %Here I am using Ax=b for the above dtVx vector
A = diag(A);
[V,D] = eigs(A);
lam1 = min(diag(D))
lam2 = max(diag(D))
%END OF ARNOLDI
S = sparse(A);
S = S \ A;
b = rho*resolx.^2; %This b is from the dtVx vector above
r = S.*(b - A.*dtVx); %redidual
thet = 0.5*(lam2+lam1); %Theta value
delt = 0.5*(lam2-lam1); %Delta value
sig = thet/delt; %sigma value
rho1 = 1/sig;
dk = r/thet;
xk = 0;
nbv = 44229;
for k = 1:nbv
xk = xk + dk;
r = r - S*A*dk;
rhok = 1/(2*sig-rhok);
dk = (rhok*rho1*dk + 2*rhok*r)/(delt+1);
rho1 = rhok
end
xk = xk + dk;
0 commentaires
Réponse acceptée
Paul
le 22 Avr 2023
Hi Nathan,
Can't run the code because not all of the information needed to run it is provided.
dtVx = rho*resolx.^2./(4*eta_nbv*(1.+eta_b)*4.1); %dtVx is a very large vector
Comparing the code to the Algorithm, it looks like there may be a few issues.
A = (4*eta_nbv*(1.+eta_b)*4.1); %Here I am using Ax=b for the above dtVx vector
A = diag(A);
[V,D] = eigs(A);
lam1 = min(diag(D))
lam2 = max(diag(D))
%END OF ARNOLDI
S = sparse(A);
This line should result in S being the identity matrix, to within numerical precision.
S = S \ A;
b = rho*resolx.^2; %This b is from the dtVx vector above
This line using element-wise multiplication and not matrix multiplication. Is thar correct? Also, this value of r is used and updated in the loop, even though the Algorithm doesn't initialize it as r_1. Is that a typo in the Algorithm?
r = S.*(b - A.*dtVx); %redidual
thet = 0.5*(lam2+lam1); %Theta value
delt = 0.5*(lam2-lam1); %Delta value
sig = thet/delt; %sigma value
rho1 = 1/sig;
dk = r/thet;
xk = 0;
nbv = 44229;
for k = 1:nbv
xk = xk + dk;
r = r - S*A*dk;
rhok = 1/(2*sig-rhok);
The next line doesn't look correct. Accoding to the algorithm it should be
dk = (rhok*rho1*dk + 2*rhok*r)/(delt+1);
% dk = rhok*rho1*dk + 2*rhok*r/delt % should be this?
rho1 = rhok;
end
xk = xk + dk;
0 commentaires
Plus de réponses (0)
Voir également
Catégories
En savoir plus sur Sparse Matrices dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!