Inequality constraints in calculus of variations

$\def\d{\mathrm{d}}$It turns out that Yuri's answer to my earlier question, whilst correct (and I thank him for his effort), was not quite what I desired. I had not posed the question properly, so I have chosen to re-ask as I am still struggling to implement the strict inequality constraints.

Let me begin by reposing the question. I wish to extremise $$Q = \int_0^h u \, \, \d y$$ subject to the constraints$$ B = B_0,\ M = M_0,$$ where $$B = \int_0^h ug \, \d y, \ M = \int_0^h \left( u^2 + \int_0^y g(y^*) \, \d y^*\right) \, \d y.$$ An important new difference is that we have that $h$, the upper limit, can vary, and crucially we need

$$u(h) = 0, \ g(h) = 0, \\ \int_0^h u \, \d y> 0, \ \int_0^h g \, \d y = G(h) > 0.$$

Here is my attempt (mostly following @Yuri). We set $$ G(y) = \int_0^y g(y^*) \, \d y^*$$

such that $$B = \int_0^h uG' \, \d y, \ M = \int_0^h u^2 + G \, \d y$$ and hence the Lagrangian becomes

$$L(u,G,G',\lambda,\mu) = Q + \lambda(B-B_0) + \mu(M-M_0).$$ The variations of this (using Euler-Lagrange $\frac{\partial L}{\partial F} - \frac{\d}{\d y}\frac{\partial L}{\partial F'} = 0$) become

$$1 + \lambda G' + 2\mu u = 0, \ \mu - \lambda u' = 0, \\ \int_0^h uG' = B_0, \ \int_0^h u^2 + G = M_0.$$

The second equation here tells us

$$u' = \frac{\lambda}{\mu} \Longrightarrow u = \frac{\lambda}{\mu}y + u_0.$$

Now since $h$ can vary we have two further transversality conditions

$$\left.\frac{\partial L}{\partial G'}\right|_{y=h} = 0, \ \left.L - G' \frac{\partial L}{\partial G'}\right|_{y=h} = 0.$$

The first transversality condition tells us that

$$\lambda u(h) = 0 \Longrightarrow u(h) = -\frac{\lambda}{\mu} \Longrightarrow u = -\frac{\mu}{\lambda}\left(h-y\right).$$

This is great as it satisifes $u(h) = 0$.

The second transversality condition tells us (using the first variation equation to give us $G$) that

$$\mu G(h) = 0 \Longrightarrow h = \frac{\lambda}{\mu^2}.$$

(However I think this condition violates the requirement $G(h) > 0$)

We can now calculate everything we need,

$$B_0 = -\frac{1}{6\mu^3}, \ M_0 = \frac{\lambda}{2\mu^4} \\ \lambda = \frac{M_0}{3\left(6^{1/3}B_0^{4/3}\right)}, \ \mu = -\frac{1}{6^{1/3}B_0^{1/3}}, \ h = \frac{2^{1/3}M_0}{3^{2/3}B_0^{2/3}} \\ Q = \frac{M_0}{6^{1/3}B_0^{1/3}}.$$

If we now sub in some values $B_0 = 6, M_0 = 36.5$ I am especially happy with the answer for $Q \approx 11.05, h \approx 6.696$ as it verifies a separate result I have reached. The problem comes with the profile for $g$. If look at the profile for $u$ we have everything behaving as expected. $u(h) = 0$ and $\int_0^h u \, \, dy >0$.

However, the issue comes when looking at the profile for $g$. We can immediately see that $g(h) \neq 0$ (strangely $g(h/2) = 0$), and also we have $\int_0^h g \, \d y = 0$ which is not strictly bigger than 0.

So my question is, how can I build the machinery for the constraints on $g$ into the Lagrangian. I have read a little online about possibly using slack variables, but I am not sure how to implement them. Any pointers would be greatly appreciated.

Interestingly, if I were to transform the profile of $g$ to $\hat{g}$ such that $g(h/2) = \hat{g}(h) = 0$ and $g(0) = g^*, \hat{g}(0) = g^* / 2$ the solution for $\hat{g}$ would match the separate result I have reached identically.

EDIT

So I have a found a paper in which the author uses slack variables to impose a condition $\theta_z \geq 0$, my thinking is that this condition is equivalent to $G' >0$. Though in his derivation, his multiplier becomes a function of $z$. Something that I have never seen before.

EDIT 2

I am now not convinced the second transversality condition is correct. This condition $G(h) = \int_0^h g = 0$ is what is forcing the profile to be negative above $h/2$ (to keep the area 0). This actually contradicts the essential requirement that $\int_0^h g > 0$.

EDIT 3

I thought it would be helpful to add some further motivation. This method of calculus of variations as described above finds two solution profiles for $u$ and $g$ that are given by $$ u_1 = \frac{3 B_0 (h - y)}{M_0} \\ g_1 = \frac{3 \left(6 B_0^2 (h-y)-6^{1/3} B_0^{4/3} M_0\right)}{M_0^2} $$

which we can see satisfy the equations and maximise $Q$.

I believe that the true solution I am seeking is given by profiles

$$ u_2 = \frac{3 B_0 (h - y)}{M_0}, \ g_2 = \frac{9 B_0^2 (h-y)}{2 M_0^2}, \ h = \frac{2^{1/3}M_0}{3^{2/3}B_0^{2/3}}.$$

which can be shown to satisfy all the equations for $M,B$ and also give the same maximal $Q$. Furthermore, the profile for $g_2$ in particular satisfies the requirements that $g_2(h) = 0$ and $\int_0^h g_2 > 0$.

It would appear to me that the $g_1, g_2$ profiles are equivalent solutions and I would expect that the calculus of variations to find both? I find it confusing that it doesn't.

I suppose one could frame the question in the way. What constraints needed to be added to the Lagrangian such that the calculus of variations gives the me the desired $g_2$ profile, rather than the unphysical $g_1$ profile?

Solution 1:

So, the task is to maximize the value of $$Q=\int\limits_0^h udy,\qquad(1)$$ taking in account the integral relations of $$\begin{cases} B(uG') = \int\limits_0^h uG'dy = B_0\\ M(u^2+G) = \int\limits_0^h (u^2 + G)dy = M_0\\ G(y) = \int\limits_0^y g(y^*)dy^*\\ \int\limits_0^h g dy > 0, \end{cases}\qquad(2)$$ and boundary conditions of $$\begin{cases} u(h) = 0\\ g(h) = 0. \end{cases}\qquad(3)$$ That task can be solved by using of Euler-Lagrange theorem of variations $$\dfrac{\delta \mathcal L(y,f,f')}{\delta f} = \mathcal L'_f - \dfrac{d}{dy}\mathcal L'_{f'}.$$ Then, for the Lagrangian in the form of $$\mathcal L(u, G, \lambda, \mu) = Q(u) + \lambda(B(uG')-B_0) + \mu(M(u^2+G)-M_0)$$ the requirements of zero variations with attention to $u$ and $G$ give the next complement equations: $$\begin{cases} 1 + \lambda g + 2\mu u = 0\\ -\lambda u' + \mu = 0.\\ \end{cases}\qquad(4)$$ $(4.2)$ and $(3.1)$ leads to $$u= - \frac\mu\lambda(h-y),\qquad(5)$$ then from $(4.1)$ $$g(y) = \dfrac{2\mu^2}{\lambda^2}(h-y) - \dfrac1\lambda.$$ At this point it turns out that condition $(3.2)$ leads to $$\dfrac1\lambda = 0,\quad u=0,\quad g=0,$$ so it contradicts with $(2.4).$

At that time, the physical representation of the task means that $h$ isn't a constant and must be used as optimization variable in the form $$\dfrac{\partial \mathcal L}{\partial h} = 0,$$ or $$u(h) + \lambda u(h)g(h) + \mu (u^2(h)+G(h)) = 0,$$ $$G(h) = 0.$$ Then it's easy to get $$G(y) = \dfrac{\mu^2}{\lambda^2}(2hy - y^2) - \dfrac y\lambda,$$ $$\dfrac{\mu^2}{\lambda^2}h^2 - \dfrac h\lambda = 0,$$ and that gives $$\begin{cases} h = \dfrac\lambda{\mu^2}\\ g(y) = \dfrac{\mu^2}{\lambda^2}(h-2y)\\[4pt] G(y) = \dfrac{\mu^2}{\lambda^2}y(h-y). \end{cases}\qquad(6)$$ Then $$B_0 = -\,\dfrac{\mu^3}{\lambda^3}\int\limits_0^h(h-y)(h-2y)\,dy = -\,\dfrac{\mu^3}{\lambda^3}\left(h^2y - \dfrac32hy^2 + \dfrac23y^3\right)\biggr|_0^h = -\,\dfrac{\mu^3h^3}{6\lambda^3},$$ $$M_0 = \,\dfrac{\mu^2}{\lambda^2}\int\limits_0^h\left((h-y)^2 + y(h-y)\right)\,dy = \,\dfrac{\mu^2}{\lambda^2}h\int\limits_0^h(h-y)\,dy = \,\dfrac{\mu^2h^2}{2\lambda^2}.$$ This leads to $$\dfrac\mu\lambda h = -\,\sqrt[3]{6B_0\,}$$ $$\,\sqrt[3]{6B_0\,} = \,\sqrt{2M_0\,},$$ $$Q = -\,\dfrac\mu\lambda\int\limits_0^h (h-y)\,dy = -\,\dfrac\mu\lambda\dfrac{ h^2}2 = \sqrt[3]{6B_0}\dfrac h2.$$

And the source data is not sufficient for an unambiguous answer