6 views (last 30 days)

Show older comments

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
on 25 Jul 2021

Edited: Star Strider
on 25 Jul 2021

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

dp = polyder(p) % Calculate Derivative

dpr = roots(dp) % Roots Of Derivative Polynomial

dpr_real = dpr(imag(dpr)==0) % Return Real Roots

dpr_in_interval = dpr_real(dpr_real>=min(points(:,1)) & dpr_real<=max(points(:,1))) % Return Real Roots Within Interval

dpryval = polyval(p, dpr_in_interval) % Evaluate Original Polynomial At That Value

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)

With the identical result.

.

Star Strider
on 25 Jul 2021

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)

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

Results = table(order,dpr_in_interval)

So the result does not change significantly (and then only in the third decimal place) beyond the 7th-order polynomial.

.

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

Start Hunting!