Remez exchange Algorithm for approximation
41 次查看(过去 30 天)
显示 更早的评论
Hello,
I tried to write Remez excahnge algorithm function but I don't know how to improve this code:
function [out] = remez(f,a,b) %% we want at the end of the function x, p et e
%% f(xi) - p(xi) = ((-1)^i)*e for i=0,...,n+1
%% we work on the interval [a,b]
%% a < x0 < ... < x(n) < x(n+1) < b : these are the points we work with in the function
x = a:0.05:b
y = f(x);
p = polyval(poly_alt_e(f,x)); %% poly_alt_e function return the alternating coefficients of the polynom p
%% and it gives also the value of e
v = y - p;
M = max(abs(v)); %% M is the maximum point of abs(f-p)
%% we want to construct new (n+2) points including M and the (n+1) points we already have
%% we have to make sure that the signs of (f-p) alternate and I don't kwow how to make sure of that.
if a <= M < x(2)
if sign(v(M)) == sign(v(x(2)))
x(x == x(2)) = [];
else x = [M x];
x(end)=[];
end
%%% if M is between two points xj and xj+1 :
%% we reject the point where the signs of (v(xj) or v(xj+1) and (v(M)) are equals
for i = 1:length(x)
M = max(abs(v));
while sign(v(M))~=sign(v(i))
if i < M < i+1
x(x == i) = [];
else x(x == i+1) = [];
end
end
end
if v(end-1)< M <= b
if sign(v(M)) == sign(v(x(end-1)))
x(x == x(end-1)) = [];
else x = [x M];
x(end)=[];
end
for i = 1:length(x)
M = max(abs(v));
while sign(v(M))~=sign(v(i))
if i < M < i+1
x(x == i) = [];
else x(x == i+1) = [];
end
end
end
%% we have to stop when the value of e is not changing anymore
%% or when x (here M) is not changing.
end
Can anyone tell me how to have the polynom p, the value e and x (here M in my function) ??
Thank you so much!
0 个评论
采纳的回答
Jan
2021-10-13
编辑:Jan
2021-10-13
if a <= M < x(2)
if v(end-1)< M <= b
if i < M < i+1
These commands will not do, what you want. The comparison is evaluated from left to right. The first expression a <= M replies TRUE or FALSE, which are treated as 1 or 0. This is the input for the 2nd comparison: 1 < x(2) or 0 < x(2) respectively. I assume you want:
if a <= M && M < x(2)
if v(end-1) < M && M <= b
if i < M && M < i+1
This looks strange:
for i = 1:length(x)
M = max(abs(v));
while sign(v(M))~=sign(v(i))
if i < M < i+1
x(x == i) = [];
else x(x == i+1) = [];
end
end
end
M = max(abs(v)) is calculated in each iteration with the same input. Move the call before the loop and do it once only to save time. But
M = max(abs(v));
while sign(v(M)) ~= sign(v(i))
might be a bug already: M is the maximum absolute value and it is used as index of v. Perhaps you want the index of the maximum absolute value?
[~, M] = max(abs(v));
The missing comments make it hard to guess, what each line of code should do. Code without clear comments cannot be maintained or debugged efficiently.
0 个评论
更多回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Creating and Concatenating Matrices 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!