Finding Maximum of Symbolic Function of Fourth Derivative Using MATLAB

Question

I am trying to find the maximum $y$-value of the fourth derivative of the function $f(x) = \frac{1}{1.1+cos(x)}$.

I am limited to real numbers. I know that the answer should be 6100 and occurs at $x_0=-\pi$ and $x_1=\pi$.

I am having problems trying to to write code in MATLAB to solve the problem. I am still relatively new to MATLAB, so it might be the way I go about trying to find multiple critical points.

Logically, what my program does is that it looks at the function $f(x)$, finds the fourth derivative, finds the critical points looking at the fifth derivative, and then substitutes those values into the fourth derivative to find the maximum $y$-value.

Here is what I have so far:

syms x
format long

fx = 1/(1.1 + cos(x));

for i = 1:4 % Computes fourth derivative of fx
    fx = diff(fx); 
end

fourth_derivative = fx;

fx = diff(fx);
fifth_derivative = fx;

% Set fifth_derivative = 0 and solve
crit_pts = real(solve(fifth_derivative))

% This is the y value that my code thinks is the maximum
max_y = eval(subs(fourth_derivative,crit_pts))

% This is the correct value of the maximum of y
correct_max_y = eval(subs(fourth_derivative,pi))

% Plotting to show that my code can't find the max
ezplot(fourth_derivative)
hold on
plot(crit_pts, subs(fifth_derivative,crit_pts), 'ro')
title('Max and Min of fourth derivative')
axis([-6 6 -100 6200])
hold off

Answer 1

Reliably finding multiple roots is challenging. As many of these roots can't be found analytically for your system, they must be found numerically.

Matlab's solve function is a general method. It tries many things in order to return solutions, but it does not guarantee that it will return all solutions (indeed, your system, being periodic, has an infinite number). I recommend reading the documentation for solve. The first thing it tries is to find an exact symbolic solution. In your case, it appears that it found one root at zero, but it did not find any other so it stopped. If it had not found any solution it would have attempted to find a root numerically (similar to calling vpasolve directly).

What can you do? In Matlab you can take advantage of some of the feature of MuPAD that have not yet been fully-integrated into the Symbolic Math toolbox. In particular, the numeric::realroots (documentation) and numeric::solve (documentation)functions are of interest. Here's an example of how you might use the latter for your system:

syms x
fx = 1/(1.1 + cos(x));
fourth_derivative = diff(fx,4);
fifth_derivative = diff(fourth_derivative);
max_y = feval(symengine,'numeric::solve',fifth_derivative,'x=-pi..pi','AllRealRoots')

which returns roots at

[ -3.1415926535897932384626433832795, 0, 3.1415926535897932384626433832795, ...
  -2.8854984092570082133339706132691, -2.3733079921147789584596786124712, ...
   2.3733079921147789584596786124712, 2.8854984092570082133339706132691]

Also, if you can bracket the root you're interested in, you can use vpasolve to find it, e.g.:

max_y = vpasolve(fifth_derivative,x,[3 3.2])

You can use a similar approach to find the root numerically via fzero.