Proof.
Let \((x^1,\ldots ,x^n)\)  be local coordinates on \(M\)  such that
\(\Phi = \mathrm {d}x^1\wedge \ldots \wedge \mathrm {d}x^n.\) 
With respect to these coordinates, the 2-vector \(P\)  is given by
\(P= \frac{1}{2} P_{ij} \frac{\partial }{\partial x^i}\wedge \frac{\partial }{\partial x^j}, \mbox{ with } P_{ij}=-P_{ji},\) 
and the vector field \(X_F^P\)  is given by \(X_F^P=\displaystyle \frac{\partial F}{\partial x^i}P_{ij}\frac{\partial }{\partial x^j}\) . As a consequence,
\((\mbox{div}_\Phi X_F^P)(x) \Phi &=& (\mathcal {L}_{X_F^P}\Phi )(x)=(\mathcal {L}_{\frac{\partial F}{\partial x^i}P_{ij}\frac{\partial }{\partial x^j}}\mathrm {d}x^1\wedge \ldots \mathrm {d}x^n)(x) \\&=& \sum _j\Big ( \mathrm {d}x^1 \wedge \ldots \mathrm {d}\Big (\sum _i \frac{\partial F}{\partial x^i}P_{ij} \Big ) \wedge \ldots \mathrm {d}x^n \Big )(x) \\&=& \sum _{i,j} \Big ( \frac{\partial ^2 F}{\partial x^j \partial x^i}_{| x}P_{ij}(x) + \frac{\partial F}{\partial x^i}_{| x}\frac{\partial P_{ij}}{\partial x^j}_{|x} \Big ) \Big ( \mathrm {d}x^1 \wedge \ldots \wedge \mathrm {d}x^n \Big )(x) \\&=& 0.\) 
Here in the last equality, it has been used that \(\mathrm {d}F(x)=0\) , and that \(P_{ij}=-P_{ji}\)  and
\(\frac{\partial ^2 F}{\partial x^j \partial x^i}_{| x}=\frac{\partial ^2 F}{\partial x^i \partial x^j}_{| x}.\) 
\(\Box \)
Examples
Kozlov's unimodularity result
Let \(\mathfrak {g}\)  be a real Lie algebra of dimension \(n\)  and \(\mathfrak {g}^*\)  be the dual vector space endowed with the Lie-Poisson structure \(\Pi _{LP}\) 
(see Example REF ).
If \(I:\mathfrak {g}^*\times \mathfrak {g}^*\rightarrow \mathbb {R}\)  is a symmetric \(\mathbb {R}\) -bilinear form then we can consider the Hamiltonian function \(H_I:\mathfrak {g}^*\rightarrow \mathbb {R}\)  given by
\(H_I(\mu ) = \frac{1}{2}I(\mu ,\mu ),\mbox{ for } \mu \in \mathfrak {g}^*.\) 
So, if \(\lbrace e_\gamma \rbrace \)  is a basis of \(\mathfrak {g}\)  with dual basis \(\lbrace e^\gamma \rbrace \)  for \(\mathfrak {g}^*\)  and \(\lbrace x_\gamma \rbrace \)  the corresponding global coordinates of \(\mathfrak {g}^*\) , we have that
\(H_I(x) = \frac{1}{2}I^{\alpha \beta }x_\alpha x_\beta .\) 
It is clear that the identity element of \(\mathfrak {g}^*\)  as an abelian Poisson-Lie group (that is, the zero vector 0) is a singular point of \(H_I.\)  Moreover, it is easy to prove that \(H_I\)  is Morse at 0 if and only if \(I\)  is non-degenerate. So, if \(I\)  is non-degenerate then, using Corollary REF  and Theorem REF , we deduce that the Hamiltonian vector field \(X_{H_I}^{\Pi _{LP}}\)  preserves a volume form on \(\mathfrak {g}^*\)  if only if \(\mathfrak {g}\)  is unimodular. Moreover, in such a case, the function \(f_0:\mathfrak {g}^*\rightarrow \mathbb {R}\)  given by (REF ) identically vanishes which implies that the Hamiltonian vector field \(X_{H_I}^{\Pi _{LP}}\)  preserves the volume form \(\Phi =dx_1\wedge \dots \wedge dx_n.\)
This extends a previous result by Kozlov [1]} for the particular case when \(H_I\)  is a Hamiltonian function of kinetic type, that is, the symmetric bilinear form \(I\)  is positive definite.

Contrast functions on Lie groups
Assume that \(\iota :G\rightarrow \mathbb {R}^n\)  is an embedding of the Lie group \(G\)  in \(\mathbb {R}^n\) , with \(\iota (e)=0\) , \(e\)  being the identity element of \(G.\)  If \(\langle \cdot ,\cdot \rangle \)  is a non-degenerate bilinear symmetric form on \(\mathbb {R}^n\) , we can consider the Hamiltonian function \(H:G\rightarrow \mathbb {R}\)  defined by
\(H=H_{\langle \cdot ,\cdot \rangle }\circ \iota :G\rightarrow \mathbb {R},\) 
with \(H_{\langle \cdot ,\cdot \rangle }:\mathbb {R}^n\rightarrow \mathbb {R}\)  given by
\(H_{\langle \cdot ,\cdot \rangle }(x)=\frac{1}{2} \langle x ,x \rangle ,\mbox{ for }x\in \mathbb {R}^n.\) 
Then, it is clear that the unit element \(e\)  is a singular point of \(H\) . In addition, one can show that if the restriction of the symmetric bilinear form \(\langle \cdot ,\cdot \rangle \)  to the subspace \(\mathfrak {g}=T_eG\subseteq T_0\mathbb {R}^n\cong \mathbb {R}^n\)  is nondegenerate then \(e\)  is a non-degenerate singular point for \(H\) . In particular, if \(\langle \cdot ,\cdot \rangle \)  is a definite positive symmetric bilinear form on \(\mathbb {R}^n\)  then \(H\)  is Morse at \(e\) . Thus, if \(G\)  is a matrix Lie group (a Lie subgroup of \(\mathrm {GL}(n,\mathbb {R})\) ) then the smooth function \(H:G\rightarrow \mathbb {R}\)  given by
\(H(A)= \mbox{Tr}\,((Id-A)(Id-A^t)), \qquad A\in G,\) 
is Morse at \(Id\in G.\)  Here \(Id\)  is the identity matrix in \(\mathrm {GL}(n,\mathbb {R}).\)
So, from Corollary REF  and Theorem REF , it follows that, for a Poisson-Lie structure \(\Pi \)  on \(G,\)  the Hamiltonian vector field of \(H\)  preserves a volume form on \(G\)  if and only if \(\Pi \)  is unimodular.
Remark 4.1 In the particular case when \(G=\mathrm {GL}(n,\mathbb {R})\) , the previous function \(H\)  was used in [2]}
as a metric contrast function on \(\mathrm {GL}(n,\mathbb {R})\) . We remark again that this kind of function plays an important role in information geometry (see [2]}, for more details).


A Hamiltonian system on the Poisson-Lie group \(\mathrm {SL}(2,\mathbb {R})\)
Let \(\mathrm {SL}(2,\mathbb {R})\)  be the special linear group of order 2, that is, the Lie subgroup of \(\mathrm {GL}(2,\mathbb {R})\)  defined by
\(\mathrm {SL(2,\mathbb {R})} = \left\lbrace  A =\begin{pmatrix}a_{11} && a_{12} \\a_{21} && a_{22} \\\end{pmatrix}\in \mathrm {GL}(2,\mathbb {R}) \,|\, \det A = 1\right\rbrace .\) 
The vector space \(\mathfrak {gl}(2,\mathbb {R})\)  of matrices of order 2 with real coefficients is the Lie algebra of \(\mathrm {GL}(2,\mathbb {R})\)  and
\(\mathfrak {sl}(2,\mathbb {R})=\lbrace  A\in \mathfrak {gl}(2,\mathbb {R})\, |\mbox{ Tr}\,A=0\rbrace \) 
is the Lie algebra of \(\mathrm {SL}(2,\mathbb {R}).\)
We can consider the basis \(\mathcal {B}=\lbrace  J_3, J_+,J_-,\displaystyle \frac{1}{2}Id\rbrace \)  of \(\mathfrak {gl}(2,\mathbb {R})\) ,
given by
\(J_3 =\begin{pmatrix}1 && 0 \\0 && -1 \\\end{pmatrix}, \qquad J_+ =\begin{pmatrix}0 && 1 \\0 && 0 \\\end{pmatrix}, \qquad J_- =\begin{pmatrix}0 && 0 \\1 && 0 \\\end{pmatrix}.\) 
Note that \(\lbrace J_3,J_{+},J_{-}\rbrace \)  is a basis of \({\mathfrak {g}}=\mathfrak {sl}(2,\mathbb {R})\)  and, moreover, the non-zero commutation relations for the basis \({\mathcal {B}}\)  read
\([ J_3,J_{-} ]_{\mathfrak {g}} = -2 J_{-}, \qquad [ J_3,J_{+} ]_{\mathfrak {g}} = 2 J_{+}, \quad [ J_+,J_- ]_\mathfrak {g}=J_3.\) 

The basis of left-invariant vector fields on \(\mathrm {GL}(2,\mathbb {R})\)  associated with \(\mathcal {B}\)  is just
\({\frac{1}{2}Id^l} &=& \frac{1}{2}\Big ( a_{11} \frac{\partial }{\partial {a_{11}}}+a_{12} \frac{\partial }{\partial {a_{12}}} + a_{21} \frac{\partial }{\partial {a_{21}}} + a_{22} \frac{\partial }{\partial {a_{22}}} \Big ) ,\\J_3^l&=& a_{11} \frac{\partial }{\partial {a_{11}}}- a_{12} \frac{\partial }{\partial {a_{12}}} + a_{21} \frac{\partial }{\partial {a_{21}}} - a_{22} \frac{\partial }{\partial {a_{22}}} , \\J_+^l &= & a_{11} \frac{\partial }{\partial {a_{12}}} + a_{21} \frac{\partial }{\partial {a_{22}}} ,\\J_-^l &=& a_{12} \frac{\partial }{\partial {a_{11}}} + a_{22} \frac{\partial }{\partial {a_{21}}}.\) 
The dual basis of \(\mathcal {B}\)  is
\(\lbrace (J^3)^l, (J^+)^l, (J^-)^l ,\mathrm {d} (\det A) \rbrace \) , \(\lbrace  J^3, J^+, J^- \rbrace \)  being the dual basis of \(\lbrace J_3,J_+,J_- \rbrace \) .
Next, consider the Poisson-Lie structure on \(\mathrm {GL} (2,\mathbb {R})\)  given by
\(\begin{split}&\lbrace  a_{11}, a_{12} \rbrace  = a_{11} a_{12}, \quad \lbrace  a_{11}, a_{21} \rbrace  = a_{11} a_{21} , \quad \lbrace  a_{11}, a_{22} \rbrace  = 2 a_{12} a_{21} , \\&\lbrace  a_{12}, a_{21} \rbrace  = 0 , \quad \lbrace  a_{12}, a_{22} \rbrace  = a_{12} a_{22} , \quad \lbrace  a_{21}, a_{22} \rbrace  = a_{21} a_{22} .\end{split}\) 
This is indeed a Poisson structure on \(\mathrm {SL} (2,\mathbb {R}),\)  since \(\det A = a_{11} a_{22} - a_{12} a_{21}\)  is a Casimir for this structure (see, for instance,
[4]}). For this Poisson-Lie group structure, we have that
\(A=\begin{pmatrix} a_{11}& 0\\ 0 &a_{22} \end{pmatrix}\in GL(2,\mathbb {R})\) 
is a singular point of \(\Pi \) .
The previous Poisson-Lie structure is the one defined by the so-called (standard) Drinfel'd-Jimbo \(r\) -matrix
and, interestingly, it is the same as the Sklyanin bracket appeared in the singular value decomposition in [5]} for Toda-SVD in [6]}.
With respect to the basis \(\lbrace J_3,J_+,J_- \rbrace \) , the (standard) Drinfel'd-Jimbo \(r\) -matrix is given by
\(r = J_- \wedge J_+ \in \wedge ^2 \mathfrak {sl}(2,\mathbb {R}).\) 
In order to compute the tangent Lie algebra to this Poisson-Lie structure using that, for coboundary structures, we have
\(\delta (\xi ) = \mathrm {ad}_\xi (r),\) 
and, from (REF ) and (REF ),
\(\delta (J_3) = 0, \qquad \delta (J_\pm ) = J_3 \wedge J_\pm .\) 
Therefore,
\([ J^3,J^\pm ]_{\mathfrak {g}^*} = J^\pm , \qquad [ J^+,J^- ]_{\mathfrak {g}^*} = 0 .\) 
So, the modular character \(\mathcal {M}_{\mathfrak {g}^*}\)  of the dual Lie algebra can be computed from (REF ) and reads
\(\mathcal {M}_{\mathfrak {g}^*} =2 J_3\ne 0 .\) 

This implies that the Poisson manifold \((\mathrm {SL(2,\mathbb {R})},\lbrace  \cdot ,\cdot \rbrace )\)  is not unimodular. In fact,
the modular vector field of the Poisson manifold \((\mathrm {SL(2,\mathbb {R})},\lbrace  \cdot ,\cdot \rbrace )\)  is
\(\frac{1}{2} (\mathcal {M}_{\mathfrak {g}^*}^r \mathcal {+} {\mathcal {M}}_{\mathfrak {g}^*}^l ) = 2 a_{11}\frac{\partial }{\partial a_{11}} -2 a_{22}\frac{\partial }{\partial a_{22}},\) 
which once more proves the Poisson structure is not unimodular. Note that modular character \(\mathcal {M}_{\mathfrak {g}}=0\) .
Now, we consider the Hamiltonian function \(H\in C^\infty (\mathrm {GL}(2,\mathbb {R}))\)  given by
\(H = \frac{1}{2} \operatorname{Tr}(A^T A) = \frac{1}{2} (a_{11}^2 + a_{12}^2 + a_{21}^2 + a_{22}^2) = \frac{1}{2} \sum _{i,j=1}^2 a_{ij}^2,\) 
Note that this Hamiltonian function generates the Toda-SVD flow defined by the Poisson bracket (REF ) (see [5]}).
Since \(\mathrm {d} H = a_{11} \mathrm {d} a_{11} + a_{12} \mathrm {d} a_{12} + a_{21} \mathrm {d} a_{21} + a_{22} \mathrm {d} a_{22}\) , it follows that
\(\mathrm {d} H (Id) \ne  0.\) 
Therefore, we cannot use Theorem REF  but nevertheless we can use Remark REF .
Indeed, the Hamiltonian vector field of \(H\)  is given by
\(X_H^{\Pi } & = & \left( -a_{11} \left(a_{12}^2+a_{21}^2\right)-2(a_{22}-1)a_{12}a_{21} \right) \frac{\partial }{\partial a_{11}} +\left(( a_{11}-1)a_{11}a_{22}-( a_{22}-1)a_{12} a_{22} \right) \frac{\partial }{\partial a_{12}} \\& &+\left( (a_{11}-1)a_{11} a_{12}-(a_{22}-1)a_{12} a_{22} \right) \frac{\partial }{\partial a_{21}} +\left( 2( a_{11}-1) a_{12} a_{21}+a_{22} \left(a_{12}^2+a_{21}^2\right) \right) \frac{\partial }{\partial a_{22}}\\&=&X^{11}\frac{\partial }{\partial a_{11}} + X^{12}\frac{\partial }{\partial a_{12}}+X^{21}\frac{\partial }{\partial a_{21}}+X^{22}\frac{\partial }{\partial a_{22}},\) 
which is the flow of Toda-SVD in dimension two.
In addition, if \(a\ne 0\)  then, from (REF ),
\(A=\left( \begin{array}{cc} a & 0 \\0 &\displaystyle \frac{1}{a} \end{array}\right)\in \mathrm {SL}(2,\mathbb {R})\) 
is an equilibrium point for the Hamiltonian vector field \(X_H^{\Pi }\)  on \(\mathrm {SL}(2,\mathbb {R})\) .
Moreover, using (REF ) and (REF ), we deduce that
\({\mathcal {M}}^l_{\mathfrak {g}^*}(A)(H)+{\mathcal {M}}^r_{\mathfrak {g}^*}(A)(H)=\frac{4}{a^2}(a^2-1)(a^2-a+1).\) 
In particular, if \(a\ne \pm 1\)  it follows that
\({\mathcal {M}}^l_{\mathfrak {g}^*}(A)(H)+{\mathcal {M}}^r_{\mathfrak {g}^*}(A)(H)\ne 0,\) 
which, from Remark REF , implies that \(X_H^\Pi \)  does not preserve a global volume form on \(\mathrm {SL}(2,\mathbb {R})\) .

A Hamiltonian system on the Poisson-Lie group \(S^3\cong \mathrm {SU}(2,\mathbb {C})\)
Let \(S^3\)  be the unit sphere in \(\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace \) . In \(\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace \)  one can consider standard global coordinates
\((x,y,z,t)\)  and the Lie group structure with multiplication
\((x,y,z,t)(x^{\prime },y^{\prime },z^{\prime },t^{\prime })=(xx^{\prime }-yy^{\prime }-zz^{\prime }-tt^{\prime }, xy^{\prime }+yx^{\prime }-zt^{\prime }+tz^{\prime },zx^{\prime }-ty^{\prime }+xz^{\prime }+yt^{\prime }, zy^{\prime }+tx^{\prime }+xt^{\prime }-yz^{\prime }),\) 
which is just the Lie group structure identifying the quaternions \(\mathbb {H}\)  with \(\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace \) . Note that the identity element is \((1,0,0,0)\)  and that \(S^3\)  is a closed normal Lie subgroup of \(\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace \) . If \(\lbrace  e_1,e_2,e_3,e_4\rbrace \)  is the canonical basis in \(\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace \) , then the associated left invariant vector fields are
\(e_1^l&=& x\frac{\partial }{\partial x}+y\frac{\partial }{\partial y}+z\frac{\partial }{\partial z}+t\frac{\partial }{\partial t}\\e_2^l&=& -y \frac{\partial }{\partial x}+x\frac{\partial }{\partial y}-t\frac{\partial }{\partial z}+z\frac{\partial }{\partial t}\\e_3^l&=& -z\frac{\partial }{\partial x}+t\frac{\partial }{\partial y}+x\frac{\partial }{\partial z}-y\frac{\partial }{\partial t}\\e_4^l&=& -t\frac{\partial }{\partial x}-z\frac{\partial }{\partial y}+y\frac{\partial }{\partial z}+x\frac{\partial }{\partial t}.\) 
Then, the non-zero commuting relations for this basis read
\([e_2,e_3]_\mathfrak {g}=-[e_3,e_2]_\mathfrak {g}=-2e_4, \qquad [e_2,e_4]_\mathfrak {g}=-[e_4,e_2]_\mathfrak {g}=2e_3, \qquad [e_3,e_4]_\mathfrak {g}=[e_4,e_3]_\mathfrak {g}=-2e_1.\) 
In fact, \(\lbrace e_2,e_3,e_4\rbrace \)  is a basis of the Lie algebra of the Lie subgroup \(S^3\cong \mathrm {SU}(2,\mathbb {C})\) . On the other hand, there is a Poisson structure \(\Pi \)  on \(\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace \)  whose Poisson bracket \(\lbrace \cdot ,\cdot \rbrace \)  is characterized by the following conditions
\(\lbrace  x,y\rbrace =-(z^2+t^2),\\\lbrace  x,z \rbrace  =yz, \quad \lbrace  x,t \rbrace  =yt,\\\lbrace  y,z \rbrace  =-xz, \quad \lbrace  y,t \rbrace  =-xt.\\\) 
Indeed, \(\Pi \)  is a Poisson-Lie structure on \(\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace \)  and
\(\Pi = -(z^2+t^2) \frac{\partial }{\partial x}\wedge \frac{\partial }{\partial y}+yz\frac{\partial }{\partial x}\wedge \frac{\partial }{\partial z}+yt\frac{\partial }{\partial x}\wedge \frac{\partial }{\partial t}-xz\frac{\partial }{\partial y}\wedge \frac{\partial }{\partial z}-xt\frac{\partial }{\partial y}\wedge \frac{\partial }{\partial t}.\) 
We remark that
\( \Vert \cdot \Vert :\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace \rightarrow \mathbb {R},\;\;\; (x,y,z,t)\rightarrow \Vert (x,y,z,t)\Vert ^2=x^2+y^2+z^2+t^2\) 
is a Casimir function and, thus, \(\Pi \)  induces a Poisson-Lie structure on \(S^3,\)  which we also denote by \(\Pi .\)  It is the standard Poisson-Lie structure on \(S^3\cong \mathrm {SU}(2,\mathbb {C})\)  (see [4]}). Moreover, using
the fact that the adjoint 1-cocycle \(\delta \)  satisfies
\([\xi ^l, \Pi ] = (\delta \xi )^l, \mbox{ for }\xi \in \mathfrak {su}(2),\) 
we deduce that
\(\delta (e_1)=\delta (e_2)=0,\quad \delta (e_3)=-e_2\wedge e_3,\quad \delta (e_4)=-e_2\wedge e_4,\) 
and for the Lie algebra structure on \(\mathfrak {g}^*\)  the non-zero commutation relations read
\(\,[e^2,e^3]_{\mathfrak {g}^*}=-[e^3,e^2]_{\mathfrak {g}^*}=-e^3, \\\,[e^2,e^4]_{\mathfrak {g}^*}=-[e^4,e^2]_{\mathfrak {g}^*}=-e^4.\) 
So, the modular character of the dual of the Lie algebra, \(\mathcal {M}_{\mathfrak {g}^*}\) , is simply \(-2e_2\) . Therefore,
the modular vector field of the Poisson manifold \((\mathbb {R}^4-\lbrace (0,0,0,0\rbrace , \Pi )\)  is given by
\(\frac{1}{2} ( \mathcal {M}_{\mathfrak {g}^*}^l+ \mathcal {M}_{\mathfrak {g}^*}^r)= -2(x\frac{\partial }{\partial y}-y\frac{\partial }{\partial x} ).\) 
Thus, it follows that \((S^3,\lbrace  \cdot ,\cdot \rbrace )\)  is not unimodular.
Now, it is easy to show that for any smooth functions \(F\colon \mathbb {R}\rightarrow \mathbb {R}\)  and \(G\colon \mathbb {R}^2\rightarrow \mathbb {R}\) , the function
\(\begin{array}{rccl}H\colon &\mathbb {R}^4-\lbrace (0,0,0,0)\rbrace & \rightarrow & \mathbb {R} \\& (x,y,z,t) & \mapsto & F(x^2+y^2)+G(z,t)\end{array}\) 
is a smooth first integral of \(\frac{1}{2} ( \mathcal {M}_{\mathfrak {g}^*}^l+ \mathcal {M}_{\mathfrak {g}^*}^r)\) . In particular, the Hamiltonian function
\(H\colon S^3\rightarrow \mathbb {R},\quad (x,y,z,t)\mapsto P(z,t),\) 
where \(P\)  is an arbitrary smooth function on \(\mathbb {R}^2\) , is a first integral of
\(\frac{1}{2} ( \mathcal {M}_{\mathfrak {g}^*}^l+ \mathcal {M}_{\mathfrak {g}^*}^r)_{|S^3}\) . Therefore, using Theorem REF , we deduce that the Hamiltonian vector field \(X^{\Pi }_H\)  preserves any left-invariant volume form on \(S^3\) .

A deformation of a conservative Lorenz system
In [9]}, the authors presented a method for obtaining integrable deformations of Lie-Poisson bi-Hamiltonian systems. The initial data is a dynamical system \(\mathcal {D}\)  on \(\mathbb {R}^n\) , which is bi-Hamiltonian with respect to two compatible linear Poisson structures. Then, they construct a bi-Hamiltonian deformation \(\mathcal {D}_\eta \)  of the dynamical system \(\mathcal {D}\) , whose bi-Hamiltonian structure is provided by a pair of Poisson-Lie structures on a non-abelian Lie group \(G_\eta \) . Under certain conditions, the limit \(\eta \rightarrow 0\)  of \(\mathcal {D}_\eta \)  is just the initial dynamical system \(\mathcal {D}\) . The Poisson coalgebra approach to Hamiltonian systems was discussed in [10]}.
In what follows, we illustrate our results with two examples coming from the Poisson-Lie deformation theory of Lie-Poisson bi-Hamiltonian systems.
As a first example, consider the integrable deformation of a conservative Lorenz system studied in [11]}, whose equations of motion read
\(\begin{split}\dot{x} &= y, \\\dot{y} &= \frac{x}{2} \left( 4 + \frac{e^{-2 \eta (z + w)} - 1}{\eta } + \eta (2 x^2 - y^2) \right), \\\dot{z} &= x y , \\\dot{w} &= 0.\end{split}\) 
This dynamical system can be written as a Hamiltonian system on the 4-dimensional Lie group \(G_{\eta } \simeq B_2 \times B_2\) , where \(B_2\)  is the unique connected simply connected non-abelian Lie group of dimension 2. We start by introducing exponential coordinates of the second kind on \(G_{\eta }\) , defined by
\(g (x,y,z,w) = \exp ( x \, X) \exp ( y \, Y) \exp ( z \, Z) \exp ( w \, W) .\) 
In terms of these coordinates the group law reads

\(\begin{array}{rcl}g (x,y,z,w) \cdot g ^{\prime }(x^{\prime },y^{\prime },z^{\prime },w^{\prime }) &=&g \bigg (x+ \displaystyle \frac{e^{ -\eta (z + w)}}{2} \bigg ( 2 x^{\prime } \cosh \left( \displaystyle \frac{\eta w}{\sqrt{2}} \right) - \sqrt{2} y^{\prime } \sinh \left( \displaystyle \frac{\eta w}{\sqrt{2}} \right) \bigg ),\\[8pt]&&y+e^{ -\eta (z + w)} \bigg ( y^{\prime } \cosh \left( \displaystyle \frac{\eta w}{\sqrt{2}} \right) - \sqrt{2} x^{\prime } \sinh \left( \displaystyle \frac{\eta w}{\sqrt{2}} \right) \bigg ), z+z^{\prime }, w+w^{\prime } \bigg )\end{array}\) 

where \(\eta > 0\)  and the Poisson bivector field is given by
\(\Pi _{B_2 \times B_2} = \frac{e^{-2 \eta (z + w)} - 1 + \eta ^2 (2 x^2 - y^2)}{4 \eta } \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial y} +\frac{y}{2} \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial z} + x \frac{\partial }{\partial y} \wedge \frac{\partial }{\partial z} .\) 
\(\Pi _{B_2 \times B_2}\)  defines a Poisson-Lie structure on \(G_\eta \)  and the dynamical system (REF ) can be recovered by means of the Hamiltonian function
\(\mathcal {H} = 2 (z-w) - x^2 .\) 
The Lie algebra \(\mathfrak {g}_{\eta }\)  associated to the Lie group \(G_{\eta }\)  is characterized by the commutators
\([X,Z]_\eta =\eta X, \quad [X,W]_\eta =\eta (X +Y),\quad [Y,Z]_\eta =\eta Y, \quad [Y,W]_\eta =\eta \left( \frac{1}{2} X + Y \right), \quad [X,Y]_\eta =[Z,W]_\eta =0.\) 
To prove that the Poisson-Lie structure (REF ) is unimodular and therefore the system (REF ) admits invariant volume forms, we need to compute the commutator map \(\delta _{B_2 \times B_2} = d_e \Pi _{B_2 \times B_2} : {\mathfrak {g}}\rightarrow \wedge ^2 {\mathfrak {g}}\) . Explicitly, we have
\(\delta _{B_2 \times B_2} (X) = Y \wedge Z, \qquad \delta _{B_2 \times B_2} (Y) = \frac{1}{2} X \wedge Z, \qquad \delta _{B_2 \times B_2} (Z) = - \frac{1}{2} X \wedge Y, \qquad \delta _{B_2 \times B_2} (W) = - \frac{1}{2} X \wedge Y .\) 
From this cocommutator, we directly obtain that the Lie algebra \(\mathfrak {g}^*\)  of the dual Poisson-Lie group is
\([\bar{X}, \bar{Y}]_{\mathfrak {g}^*} = -\frac{1}{2} (\bar{Z} + \bar{W}), \qquad [\bar{X}, \bar{Z}]_{\mathfrak {g}^*} = \frac{1}{2} \bar{Y}, \qquad [\bar{Y}, \bar{Z}]_{\mathfrak {g}^*} = \bar{X}, \qquad [\bar{W}, \bar{X}]_{\mathfrak {g}^*} = [\bar{W}, \bar{Y}]_{\mathfrak {g}^*} = [\bar{W}, \bar{Z}]_{\mathfrak {g}^*} = 0 .\) 
Therefore, \(\mathfrak {g}^* \simeq \mathfrak {sl}(2, \mathbb {R}) \oplus \mathbb {R}\) , which is clearly unimodular since it is a direct sum of two unimodular Lie algebras.
To have an invariant volume form, we obtain the left invariant vector fields as
\(\begin{array}{rcl}X^l&= & e^{ -\eta (z + w)} \cosh ( \displaystyle \frac{\eta w}{\sqrt{2}})\frac{\partial }{\partial x}-\sqrt{2}e^{ -\eta (z + w)} \sinh ( \displaystyle \frac{\eta w}{\sqrt{2}})\frac{\partial }{\partial y}\\[6pt]\end{array}Y^l&=& - \displaystyle \frac{\sqrt{2}}{2}e^{ -\eta (z + w)} \sinh ( \displaystyle \frac{\eta w}{\sqrt{2}}) \displaystyle \frac{\partial }{\partial x}+ e^{ -\eta (z + w)} \cosh \displaystyle (\frac{\eta w}{\sqrt{2}}) \displaystyle \frac{\partial }{\partial y}\\[6pt]Z^l&=& \displaystyle \frac{\partial }{\partial z}\\[6pt]W^l&=& \displaystyle \frac{\partial }{\partial w}\\\) 

and the right invariant vector fields as
\(\begin{array}{rcl}X^r&= & \displaystyle \frac{\partial }{\partial x}\\[6pt]\end{array}Y^r&=& \displaystyle \frac{\partial }{\partial y}\\[6pt]Z^r&=&-\eta x \displaystyle \frac{\partial }{\partial x}-\eta y\frac{\partial }{\partial y}+ \displaystyle \frac{\partial }{\partial z} \\[6pt]W^r&=&-(\eta x+ \displaystyle \frac{\eta }{2} y)\frac{\partial }{\partial x}- (\eta x+\eta y) \displaystyle \frac{\partial }{\partial y}+ \displaystyle \frac{\partial }{\partial w}.\\\) 
For invariant 1-forms we have that
\(\begin{array}{rcl}\bar{X}^l&= & e^{ \eta (z + w)} \cosh ( \displaystyle \frac{\eta w}{\sqrt{2}})\mathrm {d}x+\frac{\sqrt{2}}{2} e^{ \eta (z + w)} \sinh ( \displaystyle \frac{\eta w}{\sqrt{2}})\mathrm {d}y\\[10pt]\bar{Y}^l&=& \sqrt{2} e^{ \eta (z + w)} \sinh (\frac{\eta w}{\sqrt{2}})\mathrm {d}x+ e^{ \eta (z + w)} \cosh ( \displaystyle \frac{\eta w}{\sqrt{2}})\mathrm {d}y\\[6pt]\bar{Z}^l&=& \mathrm {d}z\\[6pt]\bar{W}^l&=&\mathrm {d}w\\\end{array}\) 
and
\(\begin{array}{rcl}\bar{X}^r&= & \mathrm {d}x-\eta x \mathrm {d}z +(\eta x+\eta y)\mathrm {d}w\\[6pt]\end{array}\bar{Y}^r&=& \mathrm {d}y-\eta y \mathrm {d}z+(\eta x +\eta y)\mathrm {d}w\\[6pt]\bar{Z}^r&=&\mathrm {d}z \\[6pt]\bar{W}^r&=&\mathrm {d}w,\\\) 
The modular character of Lie algebra \(\mathfrak {g}_{\eta }\)  is \(\mathcal {M}_{\mathfrak {g}_{\eta }}=-2\eta (\bar{Z}+\bar{W})\)  and so \(\mathcal {M}_{\mathfrak {g}_{\eta }}^r=-2\eta (\bar{Z}^r+\bar{W}^r)=-2\eta (\mathrm {d}z+\mathrm {d}w)=\mathrm {d}(\log f_0)\) . Therefore, \(f_0=e^{-2\eta (z+w)}\) . For \(v=\bar{X}\wedge \bar{Y}\wedge \bar{Z}\wedge \bar{W}\in \Lambda ^4\mathfrak {g}^*\) , we get the invariant volume form
\(\Phi =\sqrt{f}_0 v^l=e^{ \eta (z+w)}\mathrm {d}x\wedge \mathrm {d}y\wedge \mathrm {d}z\wedge \mathrm {d}w.\) 
 It is easy to see that \(\mathcal {L}_{{X_{\mathcal {H}}^\Pi }} \Phi = \mathrm {d}(i_{{X_{\mathcal {H}}^\Pi } }\Phi ) =0\) .
In fact, it is easy to show that \(\mathcal {L}_{X_{\mathcal {F}}^{\Pi }} \Phi = 0\) , for an arbitrary smooth function \(\mathcal {F} : G_{\eta } \rightarrow \mathbb {R}\) , in agreement with Corollary REF .
Remark 4.2 As stated above, (REF ) is an integrable deformation of a different system, which we call a conservative Lorenz system (see [11]} for details), with equations
\(\dot{x} = y, \qquad \dot{y} = x(1-z), \qquad \dot{z} = x y .\) 
Since \(B_2 \times B_2\)  is diffeomorphic to \(\mathbb {R}^4\)  as a smooth manifold, this implies that we can recover the previous dynamical system (REF ) starting from (REF ) and taking the \(\eta \rightarrow 0\)  limit. Moreover, we can also take the \(\eta \rightarrow 0\)  limit in the group law (REF ) obtaining the abelian Lie group structure on \(\mathbb {R}^4\)
\(g (x,y,z,w) \cdot g ^{\prime }(x^{\prime },y^{\prime },z^{\prime },w^{\prime }) = g \big (x + x^{\prime }, y + y^{\prime }, z+z^{\prime }, w+w^{\prime } \big )\) 
and in the Poisson-Lie bivector
\(\Pi _{\mathbb {R}^4} = -\frac{1}{2} (x+w) \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial y} +\frac{y}{2} \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial z} + x \frac{\partial }{\partial y} \wedge \frac{\partial }{\partial z} ,\) 
if we take the \(\eta \rightarrow 0\)  limit in (REF ).
Note that this Poisson-Lie structure is nothing more than the Lie-Poisson structure associated with \(\mathfrak {g}^*\) , i.e. the linear Poisson structure on the vector space \((\mathfrak {g}^*)^* \simeq \mathfrak {g}\)  associated with the Lie algebra \(\mathfrak {g}^* \simeq \mathfrak {sl}(2, \mathbb {R}) \oplus \mathbb {R}\) . Therefore, using the Hamiltonian (REF ) we obtain the system
\(\dot{x} = y, \qquad \dot{y} = x(2-z-w), \qquad \dot{z} = x y, \qquad \dot{w} = 0 .\) 
From this system, in order to recover (REF ), we simply restrict to the submanifold \(\mathcal {S} = \lbrace  (x,y,z,w) \in \mathbb {R}^4 | w=1\rbrace \) , taking into account that (REF ) is tangent to \(\mathcal {S}\) .


A deformed Euler top
Following a similar procedure as above we can consider the integrable deformation of the Euler top (see [9]} for more details), given by
\(\begin{split}\dot{x} &= e^{\eta x} (y^2 - z^2), \\\dot{y} &= \eta e^{\eta x} y z^2 - \frac{1}{2} \eta e^{\eta x} y (y^2+z^2) + \frac{\sinh (\eta x) }{\eta } (2 z - y), \\\dot{z} &= -\eta e^{\eta x} y^2 z + \frac{1}{2} \eta e^{\eta x} z (y^2 + z^2) + \frac{\sinh (\eta x)}{\eta } (z - 2 y) .\end{split}\) 
This system can be written as a bi-Hamiltonian system on the so-called “book” Lie group \(G_{\eta }\) , with the Lie algebra \(\mathfrak {g}\)  defined by
\([X,Y]=-\eta Y, \quad [X,Z]=-\eta Z, \quad [Y,Z]= 0,\) 
where \(\eta \in \mathbb {R}\)  is the deformation parameter. Since \(G_{\eta }\)  is diffeomorphic to \(\mathbb {R}^3\)  we can choose global coordinates such that the group law reads
\(g (x,y,z) \cdot g ^{\prime }(x^{\prime },y^{\prime },z^{\prime }) = g \big (x + x^{\prime }, y + y^{\prime } e^{- \eta x}, z + z^{\prime } e^{- \eta x} \big ) .\) 
The system (REF ) is bi-Hamiltonian since we can find two different Poisson-Lie structures on \(G_{\eta }\) , namely
\(\Pi _{G,0} = - z \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial y} + y \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial z} + \frac{1}{2} \bigg ( - \eta (y^2 + z^2) + \frac{e^{- 2 \eta x} - 1}{\eta } \bigg ) \frac{\partial }{\partial y} \wedge \frac{\partial }{\partial z}\) 
and
\(\Pi _{G,1} = - y \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial y} + z \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial z} + \bigg ( - \eta y z + \frac{e^{- 2 \eta x} - 1}{\eta } \bigg ) \frac{\partial }{\partial y} \wedge \frac{\partial }{\partial z} ,\) 
and two different functions
\(\mathcal {H}_0 = y z e^{\eta x} + 2 \left( \frac{\cosh (\eta x) - 1}{\eta ^2} \right),\) 
and
\(\mathcal {H}_1 = -\frac{1}{2} (y^2 + z^2) e^{\eta x} + \frac{\cosh (\eta x) - 1}{\eta ^2} ,\) 

such that (REF ) is the dynamical system defined by the Hamiltonian vector field \(\Pi _{G,0}^\sharp (d \mathcal {H}_0) = \Pi _{G,1}^\sharp (d \mathcal {H}_1)\) . Note that both Poisson-Lie structures form a Poisson pencil on \(G_{\eta }\) , i.e. the convex combination \(\lambda \Pi _0 + (1 - \lambda ) \Pi _1\)  is also a Poisson structure, for \(\lambda \in \mathbb {R}.\)
In this case, the commutator maps \(\delta _{G,i} = d_e \Pi _{G,i}\)  associated to (REF ) and (REF ), respectively, read
\(\delta _{G,0} (X) = - Y \wedge Z, \qquad \delta _{G,0} (Y) = X \wedge Z, \qquad \delta _{G,0} (Z) = - X \wedge Y ,\) 
and
\(\delta _{G,1} (X) = - 2 Y \wedge Z, \qquad \delta _{G,1} (Y) = - X \wedge Y, \qquad \delta _{G,1} (Z) = X \wedge Z .\) 
and therefore the dual Lie algebras read
\(\begin{split}&[\bar{X}, \bar{Y}]_{\mathfrak {g}^*_0} = - \bar{Z}, \qquad [\bar{X}, \bar{Z}]_{\mathfrak {g}^*_0} = \bar{Y}, \qquad [\bar{Y}, \bar{Z}]_{\mathfrak {g}^*_0} = - \bar{X}, \\&[\bar{X}, \bar{Y}]_{\mathfrak {g}^*_1} = - \bar{Y}, \qquad [\bar{X}, \bar{Z}]_{\mathfrak {g}^*_1} = \bar{Z}, \qquad [\bar{Y}, \bar{Z}]_{\mathfrak {g}^*_1} = -2 \bar{X}, \\\end{split}\) 

Since \(\mathfrak {g}^*_0 \simeq \mathfrak {so}(3)\)  and \(\mathfrak {g}^*_1 \simeq \mathfrak {sl}(2, \mathbb {R})\)  are both unimodular Lie algebras, both Poisson-Lie structures (REF ) and (REF ) are unimodular, and therefore the dynamical system (REF ) admits invariant volume forms.
Taking the \(\eta \rightarrow 0\)  limit in the equations (REF ), (REF ) and (REF ) we obtain the 3-dimensional abelian Lie group \(\mathbb {R}^3\)  together with the Poisson bivectors
\(\begin{split}\Pi _{\mathbb {R}^3,0} &= - z \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial y} + y \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial z} - x \frac{\partial }{\partial y} \wedge \frac{\partial }{\partial z} , \\\Pi _{\mathbb {R}^3,1} &= - y \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial y} + z \frac{\partial }{\partial x} \wedge \frac{\partial }{\partial z} - 2 x \frac{\partial }{\partial y} \wedge \frac{\partial }{\partial z} ,\end{split}\) 
which are precisely the Lie-Poisson structures associated with the Lie algebras \(\mathfrak {g}^* \simeq \mathfrak {so}(3)\)  and \(\mathfrak {g}^* \simeq \mathfrak {sl}(2, \mathbb {R})\) , respectively. Also, the Hamiltonian functions (REF ) and (REF ) go to
\(\mathcal {H}_0^{\eta = 0} = x^2 + y z ,\) 
and
\(\mathcal {H}_1^{\eta = 0} = -\frac{1}{2} (x^2 + y^2 + z^2) .\) 
The dynamical system defined by the Hamiltonian vector field \( \Pi _{\mathbb {R}^3,0}^\sharp (d \mathcal {H}_0^{\eta = 0}) = \Pi _{\mathbb {R}^3,1}^\sharp (d \mathcal {H}_1^{\eta = 0})\)  reads
\(\dot{x} = y^2 - z^2, \qquad \dot{y} = x(2z-y), \qquad \dot{z} = x (z-2y) ,\) 
which is equivalent to a particular case of the Euler top.
Now, we obtain the volume form preserved by the flow of the Hamiltonian vector field. A basis of left and right-invariant vector fields for the basis \(\lbrace \partial _x,\partial _y,\partial _z\rbrace \)  on \({\mathfrak {g}}_\eta ,\)  are
\(X^l=\partial _x,\quad Y^l=e^{-\eta x}\partial _y,\quad Z^l=e^{-\eta x}\partial _z,\) 
\(X^r= \partial _x-\eta y \partial _y-\eta z \partial _z,\quad Y^r= \partial _y,\quad Z^r= \partial _z,\) 
and a basis of the left and right invariant one forms for the dual basis \(\lbrace dx,dy,dz\rbrace \)  on \({\mathfrak {g}}^*_\eta \) , are
\(\bar{X}^l=dx,\quad \bar{Y}^l=e^{\eta x}dy,\quad \bar{Z}^l=e^{\eta x}dz.\) 
\(\bar{X}^r=dx,\quad \bar{Y}^r=dy+\eta y dx,\quad \bar{Z}^r=dz+\eta z dx.\) 
The modular character of Lie algebra \(\mathfrak {g}\)  is \(\mathcal {M}_{\mathfrak {g}}=-2\eta \bar{X}\)  and so \(\mathcal {M}_{\mathfrak {g}}^r =-2\eta \mathrm {d}x= \mathrm {d}(\log f_0)\) . Therefore, \(f_0=e^{-2\eta x}\) . Then, by the Corollary REF , for \(v=\bar{X}\wedge \bar{Y}\wedge \bar{Z}\in \Lambda ^3\mathfrak {g}^*\) , the volume form
\(\Phi =\sqrt{f_0} \nu ^l=e^{\eta x}dx\wedge dy\wedge dz\) 
is preserved by all Hamiltonian vector fields.
In fact, it is easy to show that \(\mathcal {L}_{X_{\mathcal {G}}^{\Pi _{G,0}}} \Phi = L_{X_{\mathcal {F}}^{\Pi _{G,1}}} \Phi = 0\) , for arbitrary functions \(\mathcal {G}, \mathcal {F} : G_{\eta } \rightarrow \mathbb {R}\) , in agreement with Corollary REF .

Conclusions and future work
We have discussed the existence of invariant volume forms for Hamiltonian systems on Poisson-Lie groups. In particular, if the identity element of the Lie group is a non-degenerate singularity of the Hamiltonian function, the Hamiltonian system admits an invariant volume form if and only if the Poisson-Lie structure is unimodular. Several examples which illustrate our theoretical results are also exhibited.
We also present some examples of Hamiltonian functions on Poisson-Lie groups with non-degenerate singularity at the identity element which appear, in a natural way, as metric contrast functions (a special class of divergence functions in information geometry) on matrix groups (see Section REF ).
On the other hand, in [14]} the authors discuss a connection between information geometry and discrete geometric mechanics. More precisely, they present the relationship between divergence functions (in information geometry) and (regular) discrete Lagrangian functions (in discrete geometric mechanics).
So, following the previous ideas, it would be interesting to study the relation between Hamiltonian dynamics on Poisson-Lie groups, discrete geometric mechanics and information geometry. This will be postponed for a future research.
In another direction, an interesting class of Poisson manifolds closely related to Poisson-Lie groups are Poisson homogeneous spaces, given by the quotient \(G/H\)  of a Poisson-Lie group \(G\)  with a closed Lie subgroup \(H\)  of \(G\)  (see [15]}). As an extension of the results in this paper, it would be interesting to discuss the existence of invariant volume forms for Hamiltonian systems on \(G/H\) . In this setting we remark the following facts:

As in the case of Poisson-Lie groups, the explicit description of the symplectic leaves and the symplectic structure on them for the Poisson structure on \(G/H\)  is, in general, a hard problem.

The Poisson homogeneous space \(G/H\)  is integrable and an explicit symplectic groupoid \({\mathcal {G}}(G/H)\)  integrating \(G/H\)  may be constructed (see [16]} and also [17]}, [18]}). But, in general, the dimension of \({\mathcal {G}}(G/H)\)  is \(2\dim (G/H)\) .

So, in order to study the existence of invariant volume forms, for a Hamiltonian system on \(G/H\) , it seems reasonable to discuss the problem directly in the homogeneous Poisson space \(G/H\) . In fact, attempts to give a description of the modular class of the Poisson structure on \(G/H\)  has been given in [19]}, [20]}, and therefore, we could discuss the relation between the unimodularity of the Poisson structure and the existence of invariant volume forms for the Hamiltonian system.
Acknowledgments: D. Iglesias, JC Marrero and E. Padrón acknowledge financial support from the Spanish Ministry of
Science and Innovation under grant PGC2018-098265-B-C32. I. Gutiérrez-Sagredo acknowledge financial support from Agencia Estatal de Investigación (Spain) under grant PID2019-106802GB-I00/AEI/10.13039/50110 0011033.
