## 3 Generation of Equations
### Realizability Analysis
#### 3.2.2 Algorithm
This is achieved by numerical schemes that provide a discrete version of the chain rule for computing the derivative of the composite function \(\mathsf{E}=\mathcal{H}(\mathbf{x})\)This is the case of the forward difference scheme, for which first order approximation of the differential applications \(\mathsf{dx}(t,\mathsf{d}t)=\frac{\mathsf{d}\mathbf{x}}{\mathsf{d}t}(t)\cdot \mathsf{d}t\) and \(\mathsf{d}\mathcal{H}(\mathbf{x},\mathsf{dx})=\
abla\mathcal{H}(\mathbf{x})^{ \intercal}\cdot\mathsf{dx}\) on the sample grid \(t\equiv kT\), \(k\in\mathbb{Z}\) are given by

\[\delta\mathbf{x}(k,T) = \mathbf{x}(k+1)-\mathbf{x}(k), \tag{20}\] \[\delta\mathcal{H}\big{(}\mathbf{x}(k),\delta\mathbf{x}(k,T)\big{)} = \mathcal{H}\big{(}\mathbf{x}(k)+\delta\mathbf{x}(k,T)\big{)}- \mathcal{H}\big{(}\mathbf{x}(k)\big{)}\] \[= \
abla_{d}\mathcal{H}\big{(}\mathbf{x}(k),\mathbf{x}(k)+\delta \mathbf{x}(k,T)\big{)}^{\intercal}\cdot\delta\mathbf{x}(k,T).\]

where, for mono-variate energy storing components (\(\mathcal{H}(\mathbf{x})=\sum_{n=1}^{n_{S}}h_{n}(x_{n})\)), the \(n\)-th coordinate is given by

\[[\
abla_{d}\mathcal{H}\big{(}\mathbf{x},\mathbf{x}+\delta\mathbf{x}\big{)}]_{ n}=\left\{ \begin{array}{ll}\frac{h_{n}\big{(}x_{n}+\delta x_{n}\big{)}-h_{n} \big{(}x_{n}\big{)}}{\delta x_{n}}&\text{if }\delta x_{n}\
eq 0,\\ h^{\prime}_{n}(x_{n})&\text{otherwise}.\end{array} \right. \tag{22}\]

A discrete chain rule is indeed recovered

\[\frac{\delta\mathsf{E}(k,T)}{T}=\
abla_{d}\mathcal{H}\big{(}\mathbf{x}(k), \mathbf{x}(k+1)\big{)}^{\intercal}\cdot\frac{\delta\mathbf{x}(k,T)}{T} \tag{23}\]

so that the following substitution in (5)

\[ \begin{array}{rcl}\frac{\mathsf{dx}}{\mathsf{d}t}(t)&\to&\frac{\delta \mathbf{x}(k,T)}{T}\\ \
abla\mathcal{H}(\mathbf{x})&\to&\
abla_{d}\mathcal{H}\big{(}\mathbf{x}(k), \mathbf{x}(k+1)\big{)}\end{array}  \tag{24}\]

leads to

\[ \begin{array}{rcl}0&=&\mathbf{b}(k)^{\intercal}\cdot\mathbf{J}\cdot\mathbf{b }(k)=\mathbf{b}(k)^{\intercal}\cdot\mathbf{a}(k)\\ &=&\underbrace{\left[\
abla_{d}\mathcal{H}^{\intercal}\cdot\frac{\delta \mathbf{x}}{\delta t}\right](k)}_{\frac{\delta\mathbf{x}(k,T)}{T}}+\underbrace {\mathbf{z}(\mathbf{w}(k))^{\intercal}\cdot\mathbf{w}(k)}_{\mathsf{D}(k)}- \underbrace{\mathbf{u}(k)^{\intercal}\cdot\mathbf{y}(k)}_{\mathsf{S}(k)}.\end{array}  \tag{25}\]

**Remark 4** (Multi-variate components). : _The case of mono-variate energy storing components covers most of the applications in electronics. Additionally, a generalization of the discrete gradient for multi-variate Hamiltonians such that Equations (20) and (21) are satisfied is given in Appendix C._

In this paper, we consider the class of the PHS composed of a collection of linear energy storing components, with quadratic Hamiltonian \(h_{n}(x_{n})=\frac{x_{n}^{2}}{2\mathcal{C}_{n}}\) (\(C_{n}\) is a capacitance or an inductance and we define \(\mathbf{Q}=\text{diag}(C_{1}\cdots C_{n_{S}})^{-1}\)). Then the discrete gradient (22) reads

\[\
abla_{d}\mathcal{H}\big{(}\mathbf{x},\mathbf{x}+\delta\mathbf{x}\big{)} = \mathbf{Q}\left(\mathbf{x}(k)+\frac{\delta\mathbf{x}(k)}{2}\right), \tag{26}\]

which restores the midpoint rule that coincides in this case with the trapezoidal rule. For nonlinear cases, (22) leads to another numerical scheme depending on the nonlinearity, still preserving passivity (see (25) and SS4.3). ### Solving the Implicit Equations

Injecting the numerical scheme (26) in (5) and solving for the quantity \(\delta\mathbf{x}(k)=\mathbf{x}(k+1)-\mathbf{x}(k)\) leads to the following energy-preserving numerical system:

\[\left( \begin{array}{c}\delta\mathbf{x}(k)\\ \mathbf{w}(k)\\ \mathbf{y}(k)\end{array} \right)=\left( \begin{array}{ccc}\mathbf{A}_{ \mathbf{x}}&\mathbf{B}_{\mathbf{x}}&\mathbf{C}_{\mathbf{x}}\\ \mathbf{A}_{\mathbf{w}}&\mathbf{B}_{\mathbf{w}}&\mathbf{C}_{\mathbf{w}}\\ \mathbf{A}_{\mathbf{y}}&\mathbf{B}_{\mathbf{y}}&\mathbf{C}_{\mathbf{y}}\end{array}  \right)\cdot\left( \begin{array}{c}\mathbf{x}(k)\\ \mathbf{z}\big{(}\mathbf{w}(k)\big{)}\\ \mathbf{u}(k)\end{array}