I am trying to code the following maximization program that I am not able to solve algebraically (Of course I manage to get the first derivative, but struggle finding its root):

$\begin{equation} Max_{\{F\}} F \cdot (1+\alpha \cdot ln(1-\frac{F}{\alpha})) \end{equation} \\ \text{Subject to: }F>0, F<\alpha \\ \text{where } \alpha > 0 \text{ is a parameter (real).}$

The function is indeed concave and has indeed a solution (unique), based on the graph:

f[a_Real, x_Real] := x*(1 + a*Log[1 - x/a])
Manipulate[Plot[f[a, x], {x, 0, a}], {a, 0, 5, 0.01}]

So far I've tried the brute force (no restriction on domain(s)) following code (solving for first derivative equal zero):

g[a_Real, x_Real] := D[x*(1 + a*Log[1 - x/a]), x]
Solve[g[a, x] == 0,x]

Unfortunately, I get the following error "Solve: Inverse functions are being used by Solve​, so some solutions may not be found; use Reduce for complete solution information."

I am VERY new to Wolfram Mathematica, and don't really understand how to do this...

Would there be a more efficient code that would allow me to get a closed form solution for $F$ (as a function of $\alpha$) solving this program ?

  • $\begingroup$ Similar questions were asked and answered at this forum a lot. The first step is g[a_?NumericQ] := Maximize[{x*(1 + a*Log[1 - x/a]), x >= 0 && x <= a}, x]. For example, g[1] performs {((E^2 - E^ProductLog[E^2]) (1 + Log[1 - (E^2 - E^ProductLog[E^2])/E^2]))/E^2, {x -> ( E^2 - E^ProductLog[E^2])/E^2}}. You may prefer NMaximize instead of. $\endgroup$
    – user64494
    Jul 19 at 17:40

You can get an analytical solution for the maximum with a litte trick (here called sol3).

I know some people don't like this way, but if you test the result numericaly (here sol4), i see no obstacle.

Since Solve does not work for general a, but for given numbers, insert an otherwise never used constant (here EulerGamma), solve and reinsert a to find the x-coordinate for maximum.

f[a_, x_] = x*(1 + a*Log[1 - x/a]);

sol3[a_] = 
  First@Solve[{D[f[EulerGamma, x], x] == 0, 
  D[f[ EulerGamma, x], {x, 2}] < 0, x < EulerGamma}, x, Reals] /. 
  EulerGamma -> a // Expand

(*   {x -> a - a E^(-1 - 1/a + ProductLog[E^(1 + 1/a)])}   *)

Test with explicit numerical solution gives the same.

sol4[a_] := 
  First@Solve[{D[f[a, x], x] == 0, D[f[ a, x], {x, 2}] < 0, x < a}, x, 

Table[{(x /. sol3[a]), (x /. sol4[a])}, {a, 10^-2, 5, Pi/5}] // 
  N // TableForm

Define point coordinates for the found maximum.

pt[a_] = Point[{x /. sol3[a], f[a, x /. sol3[a]]}] // Simplify

  Plot[f[a, x], {x, 0, a}, 
    Epilog :> {Red, PointSize[.03], N[pt[a]]}], {{a, 1}, 0, 5}]

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.