How to treat sequential numerical optimization (bilevel optimization)?

Question

I am trying to understand how to set up sequential optimization problems. I first want to understand how to treat a basic problem like

\begin{equation} \max_{x \in \Omega} f(x,m(x)) \quad \text{with} \quad m(x) = \max_{y}g(x,y) \end{equation}

In my actual problem, I have very complicated nonlinear vector functions and several constraints building the set $\Omega$.

But, let us consider an example with \begin{eqnarray} f(x,m(x)) &=& \sin(x)m(x) \ , \\ m(x) &=& \max_{y} g(x,y) \ , \\ g(x,y) &=& -(y-\sin(x))^2+\cos(x) \ , \\ \Omega &=& \{x \ | \ 0 \leq x \leq \pi \} \end{eqnarray} Naturally, $m(x) = \cos(x)$ and $\max_{x \in \Omega}f = 1/2$ for $x=\pi/4 \approx 0.785398$ holds. But suppose you dont have an analytical solution for the maximum of $g$ in respect to $y$ (due to high nonlinearity) and you want to treat the problem numerically.

How would you approach this problem numerically?

My attempt is to define all functions with numerical arguments and use NMaxValue and NMaximize, but this takes forever (around 17 mins on my laptop) already for this 1D problem.

g[x_?NumericQ, y_?NumericQ] := -(y - Sin[x])^2 + Cos[x];
m[x_?NumericQ] := NMaxValue[g1[x, y], y];
Omega = ImplicitRegion[0 <= x <= Pi, x];
f[x_?NumericQ] := Sin[x]*m[x];
DateString[]
NMaximize[{f[x], 0 <= x <= Pi}, x]
DateString[]
DateDifference[t1, t2, {"Minute", "Second"}]

"Sat 30 Dec 2017 15:17:01"

{0.5, {x -> 0.785398}}

"Sat 30 Dec 2017 15:33:54"

16 min 53 s

EDIT: I just learned, that this kind of problems is referred to as bilevel optimization, see Wikipedia. How do you treat these problems with Mathematica?

Answer 1

Assuming that the maxima of g is in the interior and not on a boundary, one can use the constraint $\frac{\partial g}{\partial y}=0$. Then, the maxima can be obtained with:

Maximize[{f[x, g[x, y]], Derivative[0, 1][g][x, y] == 0}, {x, y}] //TeXForm

$\left\{\frac{1}{2},\left\{x\to -2 \left(24 \pi -\tan ^{-1}\left(\frac{4+3 \sqrt{2}}{-2-\sqrt{2}}\right)\right),y\to \frac{2 \left(-1-\sqrt{2}\right)}{1+\left(-1-\sqrt{2}\right)^2}\right\}\right\}$