function interpolation and root finding

2021 7 25
1 个回答
回答已采纳
更新时间：2021 7 25
11 次查看（30 天）

0 个投票

Hello community,
i have 37 points in in the interval [-.5 .5]. and i'm looking for a minimum. I have written a script which uses polyfit for the polynomial and polyder with fzero to find the minium of the interpolated function. but I get some warnings..
points =[...
0.500000000000000 0.500000000000000;...
0.498097349045873 0.494452298903742;...
0.492403876506104 0.477882019169956;...
0.482962913144534 0.450598639000461;...
0.469846310392954 0.413355783194506;...
0.453153893518325 0.367259042487410;...
0.433012701892219 0.313952225075596;...
0.409576022144496 0.255588017440073;...
0.383022221559489 0.194787975761779;...
0.353553390593274 0.134485572953294;...
0.321393804843270 0.077602649602957;...
0.286788218175523 0.026729010250023;...
0.250000000000000 -0.016230934343725;...
0.211309130870350 -0.050225749043852;...
0.171010071662834 -0.075118583665429;...
0.129409522551260 -0.091527009978884;...
0.086824088833465 -0.100617488371580;...
0.043577871373829 -0.103815701636865;...
0.000000000000000 -0.102599285091917;...
-0.043577871373829 -0.098307643482893;...
-0.086824088833465 -0.092067023864118;...
-0.129409522551260 -0.084735622572768;...
-0.171010071662834 -0.076929476902457;...
-0.211309130870350 -0.069040079308481;...
-0.250000000000000 -0.061295542051084;...
-0.286788218175523 -0.053796351887970;...
-0.321393804843270 -0.046573711545795;...
-0.353553390593274 -0.039619979531691;...
-0.383022221559489 -0.032931563548278;...
-0.409576022144496 -0.026527753568514;...
-0.433012701892219 -0.020473888300108;...
-0.453153893518325 -0.014883984708655;...
-0.469846310392954 -0.009918994629169;...
-0.482962913144534 -0.005769046336078;...
-0.492403876506104 -0.002628877217708;...
-0.498097349045873 -0.000667477816139;...
-0.500000000000000 0];
I have written a script which give me the interpolated polynomial, but i get the Warning: Polynomial is badly conditioned. Add points with distinct X values or reduce the degree
I tested the function polyfit also with: [P,S,MU] = polyfit(X,Y,N) but i did not helped me.
Maybe the condition number of the Vandermod matrix is not low enough for the interpolation.
Does somebody know some alternatives? Maybe with an Aitken's Algorithm or Lagrange Interpolation ?
It would be nice to find a numerical accurate solution to interpolate to a polynomial degree of numel(x) -1.
Here is the code is used:
x=points(:,1);
y=points(:,2);
[P, S, mu]=polyfit(x,y,numel(x)-1);
[~,i_xmin]=min(x);
x=fzero(@(x) polyval(polyder(P),x,S,mu),x(i_xmin));
f=polyval(P,x,S,mu);
Best regards
emjay

请先登录，再回答此问题。

 采纳的回答

Star Strider
Star Strider 2021-7-25
编辑：Star Strider 2021-7-25
Ran in:

1 个投票

Ran in:
Try this —
points = [0.500000000000000 0.500000000000000
0.498097349045873 0.494452298903742
0.492403876506104 0.477882019169956
0.482962913144534 0.450598639000461
0.469846310392954 0.413355783194506
0.453153893518325 0.367259042487410
0.433012701892219 0.313952225075596
0.409576022144496 0.255588017440073
0.383022221559489 0.194787975761779
0.353553390593274 0.134485572953294
0.321393804843270 0.077602649602957
0.286788218175523 0.026729010250023
0.250000000000000 -0.016230934343725
0.211309130870350 -0.050225749043852
0.171010071662834 -0.075118583665429
0.129409522551260 -0.091527009978884
0.086824088833465 -0.100617488371580
0.043577871373829 -0.103815701636865
0.000000000000000 -0.102599285091917
-0.043577871373829 -0.098307643482893
-0.086824088833465 -0.092067023864118
-0.129409522551260 -0.084735622572768
-0.171010071662834 -0.076929476902457
-0.211309130870350 -0.069040079308481
-0.250000000000000 -0.061295542051084
-0.286788218175523 -0.053796351887970
-0.321393804843270 -0.046573711545795
-0.353553390593274 -0.039619979531691
-0.383022221559489 -0.032931563548278
-0.409576022144496 -0.026527753568514
-0.433012701892219 -0.020473888300108
-0.453153893518325 -0.014883984708655
-0.469846310392954 -0.009918994629169
-0.482962913144534 -0.005769046336078
-0.492403876506104 -0.002628877217708
-0.498097349045873 -0.000667477816139
-0.500000000000000 0];
p = polyfit(points(:,1), points(:,2), 7) % Fit 7th-Degree Polynomial
p = 1×8
-7.5119 -6.2764 0.8450 4.0360 2.5354 0.7934 -0.0694 -0.1026
dp = polyder(p) % Calculate Derivative
dp = 1×7
-52.5832 -37.6587 4.2248 16.1441 7.6062 1.5868 -0.0694
dpr = roots(dp) % Roots Of Derivative Polynomial
dpr =
0.6763 + 0.0000i -0.4904 + 0.2286i -0.4904 - 0.2286i -0.2242 + 0.3621i -0.2242 - 0.3621i 0.0368 + 0.0000i
dpr_real = dpr(imag(dpr)==0) % Return Real Roots
dpr_real = 2×1
0.6763 0.0368
dpr_in_interval = dpr_real(dpr_real>=min(points(:,1)) & dpr_real<=max(points(:,1))) % Return Real Roots Within Interval
dpr_in_interval = 0.0368
dpryval = polyval(p, dpr_in_interval) % Evaluate Original Polynomial At That Value
dpryval = -0.1039
dpv = polyval(dp,points(:,1));
v = polyval(p, points(:,1));
figure
plot(points(:,1), points(:,2),'.')
hold on
plot(points(:,1), v)
plot(dpr_in_interval, dpryval, 'sg')
hold off
grid
legend('Original Data',' Fitted Polynomial', 'Calculated Minimum', 'Location','best')
EDIT — (25 Jun 2021 at 17:55)
It is not necessary to use fzero with polyder. Just use the roots function to find the roots of the derivative of the polynomial. It is likely more accurate than fzero, and likely much more efficient.
If you absolutely must use fzero:
pointsmin = fzero(@(x)polyval(dp,x), -0.5)
pointsmin = 0.0368
With the identical result.
.

2 个评论
显示 无 隐藏 无

Marko
Marko 2021-7-25
Hello Star Strider, thank you for the information to use roots, instead of fzero.
Do you know if there is possible increase the polynomial order to numel(x)-1 without this warning:
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
> In polyfit (line 83)
I belive that polyfit is not capable of this high order, because of the high condition number of the Vandermonde matrix. So an another algorithm is needed.
Star Strider
Star Strider 2021-7-25
Ran in:
Ran in:
My pleasure!
The problem with the polynomial order is that any order more than 7 will not produce a more accurate fit. That is simply the nature of polynomilal approximations. (Even increasing the order from 5 originally to 7 in the edited version did not change the result beyond a few non-signifcant decimal places.) I also douobt that another algorithm is necessary, simply because a higer order will not increase the accuracy of the fit.
It is relatively easy to do that experiment.
Example —
points = [0.500000000000000 0.500000000000000
0.498097349045873 0.494452298903742
0.492403876506104 0.477882019169956
0.482962913144534 0.450598639000461
0.469846310392954 0.413355783194506
0.453153893518325 0.367259042487410
0.433012701892219 0.313952225075596
0.409576022144496 0.255588017440073
0.383022221559489 0.194787975761779
0.353553390593274 0.134485572953294
0.321393804843270 0.077602649602957
0.286788218175523 0.026729010250023
0.250000000000000 -0.016230934343725
0.211309130870350 -0.050225749043852
0.171010071662834 -0.075118583665429
0.129409522551260 -0.091527009978884
0.086824088833465 -0.100617488371580
0.043577871373829 -0.103815701636865
0.000000000000000 -0.102599285091917
-0.043577871373829 -0.098307643482893
-0.086824088833465 -0.092067023864118
-0.129409522551260 -0.084735622572768
-0.171010071662834 -0.076929476902457
-0.211309130870350 -0.069040079308481
-0.250000000000000 -0.061295542051084
-0.286788218175523 -0.053796351887970
-0.321393804843270 -0.046573711545795
-0.353553390593274 -0.039619979531691
-0.383022221559489 -0.032931563548278
-0.409576022144496 -0.026527753568514
-0.433012701892219 -0.020473888300108
-0.453153893518325 -0.014883984708655
-0.469846310392954 -0.009918994629169
-0.482962913144534 -0.005769046336078
-0.492403876506104 -0.002628877217708
-0.498097349045873 -0.000667477816139
-0.500000000000000 0];
size(points,1)
ans = 37
dpr_in_interval = num2cell(NaN(size(points,1)-1,1));
for k = 1:size(points,1)-1
order(k,:) = k;
p = polyfit(points(:,1), points(:,2), k); % Fit 7th-Degree Polynomial
dp = polyder(p); % Calculate Derivative
dpr = roots(dp); % Roots Of Derivative Polynomial
dpr_real = dpr(imag(dpr)==0); % Return Real Roots
select = dpr_real(dpr_real>=min(points(:,1)) & dpr_real<=max(points(:,1))) % Return Real Roots Within Interval
if ~isempty(select)
dpr_in_interval{k,:} = select;
end
end
select = 0×1 empty double column vector
select = -0.1325
select = 2×1
-0.4525 0.0113
select = 0.0156
select = 0.0351
select = 0.0416
select = 0.0368
select = 0.0368
select = 0.0362
select = 0.0359
select = 0.0359
select = 0.0359
select = 0.0359
select = 0.0360
select = 0.0359
select = 0.0359
select = 0.0359
select = 0.0360
select = 0.0359
select = 0.0359
select = 0.0359
select = 0.0359
select = 0.0360
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0360
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT.
select = 0.0359
Results = table(order,dpr_in_interval)
Results = 36×2 table
order dpr_in_interval _____ _______________ 1 {[ NaN]} 2 {[ -0.1325]} 3 {2×1 double} 4 {[ 0.0156]} 5 {[ 0.0351]} 6 {[ 0.0416]} 7 {[ 0.0368]} 8 {[ 0.0368]} 9 {[ 0.0362]} 10 {[ 0.0359]} 11 {[ 0.0359]} 12 {[ 0.0359]} 13 {[ 0.0359]} 14 {[ 0.0360]} 15 {[ 0.0359]} 16 {[ 0.0359]}
So the result does not change significantly (and then only in the third decimal place) beyond the 7th-order polynomial.
.

请先登录，再进行评论。

更多回答（0 个）

类别

在 帮助中心 和 File Exchange 中查找有关 Polynomials 的更多信息

产品

版本

R2021a

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by