# The variational Lagrangian

Suppose we have a [wave function model](#wave-function-models) $\require{physics} \ket{\Psi(\vb{p})}$ whose the optimal parameters $\vb{p}^\star(\boldsymbol{\eta})$ don't solve any stationary equations, but rather are the solutions of a (not necessarily linear) system of $y$ equations:
\begin{equation}
    \forall \boldsymbol{\eta}: \qquad
    f_a (
        \boldsymbol{\eta},
        \vb{p}^\star(\boldsymbol{\eta})
    )
    = 0
    \thinspace ,
\end{equation}
like the equations for the cluster amplitudes for [coupled-cluster](#the-coupled-cluster-wave-function-model) types of calculations, or any general [projected Schrödinger equations](the-projected-schrodinger-equation-framework.Rmd). These types of methods have a separate energy function
\begin{equation*}
    E(\vb{p}, \boldsymbol{\eta})
    \thinspace ,
\end{equation*}
so solving the equations doesn't necessarily lead to a vanishing gradient of this energy function. This, in turn, means that we cannot ignore the first-order contribution due to the wave function response in the [Hellmann-Feynman principle](first-order-response-derivatives-for-variationally-determined-wave-functions). The corresponding wave functions, energies and parameters $\vb{p}^\star(\boldsymbol{\eta})$ are therefore said to be _non-variational_, short for _non-variationally determined_.

This non-variationality, however, is short-lived. According to Helgaker [@Helgaker2000] [@Saue2002] [@Wilson1992], we can propose the Lagrangian 
\begin{equation}
    \mathscr{L}(\boldsymbol{\eta}, \vb{p}, \boldsymbol{\lambda})
    = E(\boldsymbol{\eta}, \vb{p})
    + \sum_a^l \lambda_a f_a(\boldsymbol{\eta}, \vb{p})
\end{equation}
whose goal is to minimize the energy function $E(\boldsymbol{\eta}, \vb{p})$, subject to the constraint that the optimality equations should be fulfilled. The vector $\boldsymbol{\lambda}$ thus represents $y$ Lagrangian multipliers, one for every (non-)linear equation. The variational conditions are now expressed with respect to the Lagrangian $\mathscr{L}$:
\begin{align}
    & \forall \boldsymbol{\eta}: \qquad
    \eval{
        \pdv{
            \mathscr{L} (
                \boldsymbol{\eta},
                \vb{p}^\star(\boldsymbol{\eta}),
                \boldsymbol{\lambda})
        }{\lambda_a}
    }_{ \boldsymbol{\lambda}^\star(\boldsymbol{\eta}) }
    =
    f_a (
        \boldsymbol{\eta},
        \vb{p}^\star(\boldsymbol{\eta})
    ) = 0 \\
    %
    & \forall \boldsymbol{\eta}: \qquad
    \eval{
        \pdv{
            \mathscr{L} (
                \boldsymbol{\eta},
                \vb{p},
                \boldsymbol{\lambda}^\star(\boldsymbol{\eta})
            )
        }{p_i}
    }_{ \vb{p}^\star(\boldsymbol{\eta}) }
    =
    \eval{
        \pdv{
            E(\boldsymbol{\eta}, \vb{p})
        }{p_i}
    }_{ \vb{p}^\star(\boldsymbol{\eta}) }
    + \sum_a^y
        \lambda_a^\star(\boldsymbol{\eta})
        \qty(
            \eval{
                \pdv{
                    f_a(\boldsymbol{\eta}, \vb{p})
                }{p_i}
            }_{ \vb{p}^\star(\boldsymbol{\eta}) }
        )
    = 0
    \thinspace ,
\end{align}
in which now both sets of optimized parameters $\vb{p}^\star$ and $\boldsymbol{\lambda}^\star$ have implicit dependences on the external perturbation $\boldsymbol{\eta}$:
\begin{equation}
    \vb{p}^\star(\boldsymbol{\eta})
    \qqtext{and}
    \boldsymbol{\lambda}^\star(\boldsymbol{\eta})
    \thinspace .
\end{equation}
We should note that the first set of stationary equations recovers the original equation constraints because the Lagrangian parameters $\boldsymbol{\lambda}$ enter the Lagrangian linearly. The second set of the stationary equations then constitutes a linear system of equations for the Lagrange multipliers, which requires the extra cost of calculating the partial derivatives
\begin{equation}
    \eval{
        \pdv{
            E(\boldsymbol{\eta}, \vb{p})
        }{p_i}
    }_{ \vb{p}^\star(\boldsymbol{\eta}) }
    %
    \qquad \qqtext{and} \qquad
    %
    \eval{
        \pdv{
            f_a(\boldsymbol{\eta}, \vb{p})
        }{p_i}
    }_{ \vb{p}^\star(\boldsymbol{\eta}) }
    \thinspace .
\end{equation}

After the solutions $\vb{p}^\star(\boldsymbol{\eta})$ and $\boldsymbol{\lambda}^\star(\boldsymbol{\eta})$ have been obtained, the energy is straightforwardly calculated as an evaluation of the Lagrangian function:
\begin{equation}
    \mathcal{E}(\boldsymbol{\eta})
    = E (
            \boldsymbol{\eta},
            \vb{p}^\star(\boldsymbol{\eta})
    )
    = \mathscr{L} (
        \boldsymbol{\eta},
        \vb{p}^\star(\boldsymbol{\eta}),
        \boldsymbol{\lambda}^\star(\boldsymbol{\eta})
    )
    \thinspace ,
\end{equation}
because the variational conditions on the Lagrangian hold at the optimal point $(\vb{p}^\star, \boldsymbol{\lambda}^\star)$.

Through the use of the variational Lagrangian, we have succeeded in making a non-variational theory, variational, similarly to [regular variationally determined wave functions](#variationally-determined-wave-functions):

- For first-order response properties, we can [recover _a_ Hellmann-Feynman theorem](#first-order-response-derivatives-for-the-variational-lagrangian);
- For second-order response properties, we require the [first-order parameter and multiplier responses](#second-order-response-derivatives-for-the-variational-Lagrangian).
