How do I optimize this problem where the constraints and objective are variable?

Question

Problem Definition:

Pa = Constant
Pb = Constant
Vmax_a = Constant
Vmax_b = constant

Objective Function:

maximize: PaVa + PbVb - PoVo
Where: Vo = f(Po) ~ C  + e*Po

Constraints:

0 <= Va <= Vmax_a
0 <= Vb <= Vmax_b
Va + Vb = Vo
P_min < Po < P_max

Variables to Optimize:

Va, Vb, Po

Answer 1

Your merit function is nonlinear, continuous and differentiable.

So long as $C$ is unknown, variable and unbounded, the fact that $V_o$ may be expressed as a function of $P_o$ is immaterial. If $C$ is bounded or fixed then constraints need to be added to the demonstration below, but I assume that that is not the case.

When optimizing, $V_a$ and $V_b$ being bounded implies that $V_o$ is also bounded, but this bound can stay implicit; and the simpler form for the merit function is one where $V_o$ is substituted with $V_a + V_b$. As a cost (negative merit),

$$ c = P_o (V_a + V_b) - P_a V_a - P_b V_b $$

The gradient of the cost function is easy and can help many optimizers, so include it in your problem construction:

$$ \nabla c = \begin{bmatrix} P_o - P_a \\ P_o - P_b \\ V_a + V_b \end{bmatrix} $$

All together, an implementation could look like

import numpy as np
import pandas as pd
from scipy.optimize import check_grad, minimize

Pa = 0.3
Pb = 0.6
Pmin = 0.2
Pmax = 0.7
Vmax_a = 12
Vmax_b = 9


def cost(params: np.ndarray) -> float:
    Va, Vb, Po = params
    return Po*(Va + Vb) - Pa*Va - Pb*Vb


def cost_jacobian(params: np.ndarray) -> tuple[float, float, float]:
    Va, Vb, Po = params
    return (
        Po - Pa,
        Po - Pb,
        Va + Vb,
    )


error = check_grad(func=cost, grad=cost_jacobian, x0=(0.9, 1.1, 1.4))
assert error < 1e-6

result = minimize(
    fun=cost, jac=cost_jacobian,
    x0=(1, 1, 1),
    bounds=(
        (0, Vmax_a),
        (0, Vmax_b),
        (Pmin, Pmax),
    ),
)
assert result.success, result.message

Va, Vb, Po = result.x
Vo = Va + Vb
C = Vo - np.e*Po
df = pd.DataFrame(
    {
        'V': (Va, Vb, Vo),
        'P': (Pa, Pb, Po),
    },
    index=pd.Index(['a', 'b', 'o']),
)
print(df)

print(f'\nC = {C:.3f}')

      V    P
a  12.0  0.3
b   9.0  0.6
o  21.0  0.2

C = 20.456

Answer 2

I am not 100% sure, but I am thinking like this: I suppose you can rewrite the objective to: $$\max V_a(P_a-P_o)+V_b(P_b-P_o)$$ since $V_a + V_b = V_o$.

And since that sum is also equal to whatever is the RHS of that strange equation @prubin asked about and since either $P_a-P_o$ or $P_b-P_o$ will always be bigger depends only on $P_a$ and $P_b$ and not on $P_o$, you can assign both $V_a$ and $V_b$ from the RHS. If $P_a > P_b$ you will max out $V_a$ before assigning the remaining to $V_b$, otherwise it's just the opposite. So you can essentially search only for $P_o$ in the range $P_\min$ and $\min(P_\max, \max(P_a, P_b))$. Especially because the RHS has to be non-negative due to being equal to the sum of two non-negative numbers -> if $P_o$ would surpass $P_a$ and $P_b$ will make the objective negative, while it may be non-negative if $P_o$ is smaller than the larger one.

The problem becomes (assuming $P_a > P_b$):

\begin{align}\max&\quad {\rm RHS} \cdot (P_a - P_o)\\ \text{s.t.}&\quad 0 \le{\rm RHS}\le {V_a}_\max + {V_b}_\max\\ &\quad P_\min \le P_o \le \min(P_\max, P_a)\end{align}

And RHS being the function of $P_o$.

Not sure how to solve that optimally though. Maybe use heuristics?