Finding coefficients of a 3rd degree polynomial while i impose a numerical value to its maximum value

Hi, I have a polynomial of 3rd degree in time. a0 + a1×t + a2×t^2 + a3×t^3. I want to impose a numerical value to its maximum (say 50) in a known time interval and find out the corresponding values of the coefficients a0 to a3. Is there any direct way ? I have a tendency towards some optimization scheme. Any ideas ??

Answer 1

[John D'Errico]

Solving for a cubic polynomial such that the maximum value is constrained over an interval is not simply done, although I do it (approximately) in my SLM tool box.

The trick is to understand the difference between necessary and sufficient constraints on your goal.

Suppose that you picked a finite set of points to test the cubic polynomial at within your interval? If the chosen polynomial satisfies the constraints at each point, then you will accept that this is a necessary condition for the polynomial to be less than your goal everywhere in the interval, recognizing that this is not a sufficient condition. There may still be some points at other locations in the interval where the polynomial exceeds that maximum. Of course your function is only a cubic polynomial, so it must be smooth and continuous. Therefore if you have chosen enough (well placed) points it will not exceed the goal by a lot.

You can define these necessary constraint locations using the linear constraints found in lsqlin, so lsqlin is a viable solver for this. You need not use fmincon at all.

Or, you could just use my slmengine to solve your problem, it defines the set of points as a list of 11 values in the interval, plus it will add the endpoints of the interval to that list. Thus, this call will be all it takes to fit a single cubic polynomial over the interval [0,1], such that it (hopefully) never passes above the value 0.5. It uses the necessary set of constraints I have described above, so it may go slightly above that value though.

x = rand(10,1); y = rand(10,1); slm = slmengine(x,y,'knots',[0 1],'maxvalue',0.5,'result','pp');

This will return a spline in a pp form. Note that the single cubic segment it generates will be relative to the left hand end point, effectively evaluating P(x-x_left).

You can find SLM on the File Exchange:

http://www.mathworks.com/matlabcentral/fileexchange/24443-slm-shape-language-modeling

If you wanted to use more/different point locations, you would need to change the code, but that point is rather easy to find, and easy to modify. (I suppose I could have given that as an option. Ugh.) Or, as I said, just use lsqlin, defining the constraints yourself.

You really cannot easily do any better unless you are willing to get your hands a bit dirty. You would then use fmincon, defining a nonlinear constraint as the maximum value of the cubic over that interval. This is easily found for a cubic. Simply test the points where the cubic polynomial has a zero slope, IF those points lie within the interval of interest. Compare the value found to the value of the cubic at the end points of the interval, returning the overall maximum value as the constraint value. Use of this form of constraint with fmincon will yield a result that is both necessary and sufficient for your goal.

Finally, a flaw with the use of fmincon here is the resulting constraint will not be a differentiable function of the parameters of the polynomial. (The constraint value will be continuous, but it can under some circumstances be non-differentiable.) fmincon will probably survive it, although it may leave you at a sub-optimal location.

Answer 2

[Walter Roberson]

Calculus.

The maximum of a cubic must occur either at one of the two interval end-points, or at one of the two places where the differential of the function is 0.

The differential is a1+2*a2*t+3*a3*t^2 which has 0's at t =

(1/3)*(-a2+sqrt(a2^2-3*a3*a1))/a3

-(1/3)*(a2+sqrt(a2^2-3*a3*a1))/a3

Differentiate the original polynomial twice and look at the sign to determine whether the point will be a maximum or minimum:

sign(2*a2+6*a3*t)

Substitute the two potential solutions for the differential being 0 back in to the original expression and simplify,

(1/27)*(27*a0*a3^2 - 9*a3*a1*a2 + 6*a3*a1*sqrt(a2^2-3*a3*a1) + 2*a2^3 - 2*a2^2*sqrt(a2^2-3*a3*a1)) / a3^2

-(1/27)*(-27*a0*a3^2 + 9*a3*a1*a2 + 6*a3*a1*sqrt(a2^2-3*a3*a1) - 2*a2^3 - 2*a2^2*sqrt(a2^2-3*a3*a1)) / a3^2

To constrain the maximum to be C, pick one of those two and solve for it being equal to C subject to the constraint that the t so implied is within the target interval. Solving it for a2 will get you a longish cubic; solving for a3 gets you two more compact solutions. Pick one of the two and substitute it back in to the double-differentiated polynomial and simplify. Now impose the condition that the sign of that should be negative (for a maximum). The easiest way to approach that is to solve for the t that makes the expression 0. The expression will be linear in t so everything to one side will be negative and everything to the other side will be positive. You will end up with something like

t = 9*a2*(-a0+C)^2*a1 / (9*C*a2*a1^2 + 6*C*a2 * sqrt(a1^2*(3*C*a2-3*a2*a0+a1^2)) - 9*a2*a0*a1^2 - 6*a2*a0 * sqrt(a1^2*(3*C*a2-3*a2*a0+a1^2)) + 2*a1^4+2*a1^2 * sqrt(a1^2*(3*C*a2-3*a2*a0+a1^2)))

I have not done a formal examination of this expression, but it appears to me to have enough free play to be able to have it come out to pretty much any real value, and thus to be within whatever interval you imposed.

The implication of this approach is that you can impose your maximum, C, and impose your interval, and that you can possibly do so for any arbitrary a1 or a2 (the other would be determined by the location within the interval that you wish to be the peak.)

Minimization would be a waste of time as there are probably infinite solutions. Pick a random t in your interval to be the peak, pick a random a1 or a2, and the rest can be solved for algebraically.

Answer 3

[Bruno Luong]

Sure, a0=-Inf, a1=a2=a3=0

Bruno

Answer 4

[aft]

thanks all for your solutions. these helped me understand better my problem.