Second-order properties
In order to calculate second-order properties, we calculate the second total derivative of the energy with respect to the corresponding perturbation:
\[\begin{equation}
\begin{split}
\frac{ %
\dd^2{E(\boldsymbol{\eta})} %
}{\dd{\eta_m} \dd{\eta_n}} = %
& %
\pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}^\star, \boldsymbol{\lambda}^\star) %
}{\eta_m}{\eta_n} \\
& + \sum_i^m %
\eval{ %
\pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}, \boldsymbol{\lambda}^\star) %
}{p_i}{\eta_n} %
}_{\vb{p}^\star} %
\qty( %
\pdv{ %
p_i^\star(\boldsymbol{\eta}) %
}{\eta_m} %
) \\ %
&+ %
\sum_a^l %
\eval{ %
\pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}^\star, \boldsymbol{\lambda}) %
}{\lambda_a}{\eta_n} %
}_{\boldsymbol{\lambda}^\star} %
\qty( %
\pdv{ %
\lambda_a^\star(\boldsymbol{\eta}) %
}{\eta_m} %
) \thinspace ,
\end{split}
\label{eq:molecular_Hessians}
\end{equation}\]
from which we see that we require the first-order response of the wave function parameters
\[\begin{equation}
\vb{x} = %
\pdv{ %
\vb{p}^\star(\boldsymbol{\eta}) %
}{\boldsymbol{\eta}}
\end{equation}\]
and of the Lagrangian multipliers:
\[\begin{equation}
\vb{y} = %
\pdv{ %
\boldsymbol{\lambda}^\star(\boldsymbol{\eta}) %
}{\boldsymbol{\eta}} %
\thinspace ,
\end{equation}\]
where we have introduced the short-hand symbols
\(\vb{x}\) and
\(\vb{y}\). In order to determine the linear wave function response
\(\vb{x}\), we take the total derivative of equation
\(\eqref{eq:Lagrangian_var_lambda}\), yielding:
\[\begin{equation}
\vb{k}_{\vb{p}} \thinspace \vb{x} %
= - \vb{F}_{\vb{p}} %
\thinspace ,
\end{equation}\]
in which we will call
\(\vb{k}_{\vb{p}}\) the parameter response constant:
\[\begin{equation}
\qty( %
\vb{k}_{\vb{p}}
)_{ai} %
= %
\eval{ %
\pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}, \boldsymbol{\lambda}) %
}{\lambda_a}{p_i} %
}_{\vb{p}^\star, \boldsymbol{\lambda}^\star} %
= \eval{ %
\pdv{ %
f_a(\boldsymbol{\eta}, \vb{p}) %
}{p_i} %
}_{\vb{p}^\star}
\end{equation}\]
and
\(\vb{F}_{\vb{p}}\) the parameter response force:
\[\begin{equation}
\qty( %
\vb{F}_{\vb{p}}
)_{am} %
= %
\eval{ %
\pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}^\star, \boldsymbol{\lambda}) %
}{\lambda_a}{\eta_m} %
}_{\boldsymbol{\lambda}^\star} %
= %
\pdv{ %
f_a(\boldsymbol{\eta}, \vb{p}^\star) %
}{\eta_m} %
\thinspace ,
\end{equation}\]
in analogy with Hooke’s law. We note that the parameter response constant
\(\vb{k}_{\vb{p}}\) is related to the Jacobian of the PSEs. After having solved the linear equation for the first-order wave function response
\(\vb{x}\), we are able to solve a linear equation for the first-order Lagrange multiplier response
\(\vb{y}\):
\[\begin{align}
\vb{k}_{\boldsymbol{\lambda}} \thinspace \vb{y} %
&= - \vb{F}_{\boldsymbol{\lambda}} \\
%
&= - \vb{A}_{\boldsymbol{\lambda}} %
- \vb{B}_{\boldsymbol{\lambda}} \vb{x} %
\thinspace ,
\end{align}\]
where
\(\vb{k}_{\boldsymbol{\lambda}}\) is then accordingly called the multiplier response constant:
\[\begin{equation}
\qty( %
\vb{k}_{\boldsymbol{\lambda}}
)_{ia} %
=
\qty( %
\vb{k}_{\vb{p}}
)_{ai}
\end{equation}\]
and
\(\vb{F}_{\boldsymbol{\lambda}}\) the multiplier response force:
\[\begin{equation}
\qty( %
\vb{F}_{\vb{p}}
)_{am} %
= %
\eval{ %
\pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}^\star, \boldsymbol{\lambda}) %
}{\lambda_a}{\eta_m} %
}_{\boldsymbol{\lambda}^\star} %
= %
\pdv{ %
f_a(\boldsymbol{\eta}, \vb{p}^\star) %
}{\eta_m}
\end{equation}\]
that can be calculated using the intermediaries
\[\begin{align}
\qty( %
\vb{A}_{\boldsymbol{\lambda}}
)_{im} %
& = %
\eval{ %
\pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}, \boldsymbol{\lambda}^\star) %
}{p_i}{\eta_m} %
}_{\vb{p}^\star} \\
& = %
\eval{ %
\pdv{ %
E(\boldsymbol{\eta}, \vb{p}) %
}{p_i}{\eta_m} %
}_{\vb{p}^\star} %
+ \sum_a^l \lambda_a^\star %
\eval{ %
\pdv{ %
f_a(\boldsymbol{\eta}, \vb{p}) %
}{p_i}{\eta_m} %
}_{\vb{p}^\star}
\end{align}\]
and
\[\begin{align}
\qty( %
\vb{B}_{\boldsymbol{\lambda}}
)_{ij} %
&= %
\eval{ %
\pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}, \boldsymbol{\lambda}) %
}{p_i}{p_j} %
}_{\vb{p}^\star, \boldsymbol{\lambda}^\star} \\
&= %
\eval{ %
\pdv{ %
E(\boldsymbol{\eta}, \vb{p}) %
}{p_i}{p_j} %
}_{\vb{p}^\star} %
+ \sum_a^l \lambda_a^\star %
\eval{ %
\pdv{ %
f_a(\boldsymbol{\eta}, \vb{p}) %
}{p_i}{p_j} %
}_{\vb{p}^\star} %
\thinspace .
\end{align}\]
After having solved these two sets of linear response equations, we can calculate second-order properties as:
\[\begin{equation}
\begin{split}
\frac{ %
\dd^2{E(\boldsymbol{\eta})} %
}{\dd{\eta_m} \dd{\eta_n}} %
= %
& \pdv{ %
\mathscr{L}(\boldsymbol{\eta}, \vb{p}^\star, \boldsymbol{\lambda}^\star) %
}{\eta_m}{\eta_n} %
+ \qty( %
\vb{x}^{\text{T}} \vb{A}_{\boldsymbol{\lambda}} %
)_{mn} \\
& + \qty( %
\vb{y}^{\text{T}} \vb{F}_{\vb{p}} %
)_{mn} %
\thinspace .
\end{split}
\end{equation}\]
For the PSE framework, we list here the following energy derivatives:
\[\begin{equation}
\eval{ %
\pdv{ %
E(\boldsymbol{\eta}, \vb{p}) %
}{\eta_m}{\eta_n} %
}_{\vb{p}^\star} %
= %
\frac{ %
\matrixel**{0}{ %
\pdv[2]{ %
\hat{\mathcal{H}}(\boldsymbol{\eta}) %
}{\eta_m}{\eta_n} %
}{\Psi(\vb{p}^\star)} %
}{ \braket{0}{\Psi(\vb{p}^\star)} }
%
\end{equation}\]