Wolfram Function Repository

Instant-use add-on functions for the Wolfram Language

Function Repository Resource:

 BialternateSum

Source Notebook

Get the bialternate sum matrix of a square matrix

Contributed by: E. Chan-López & Víctor Castellanos

ResourceFunction["BialternateSum"][mat]

gives the bialternate sum of the square matrix mat.

Details

The bialternate sum of a matrix allows the study of the local stability of any nonlinear system of ordinary differential equations.
The bialternate sum matrix is a useful criterion for studying supercritical Hopf bifurcations without using eigenvalues and is equivalent to the Liu criterion using the Hurwitz matrix.
If the square matrix mat is of length n, then the bialternate sum is a square matrix of length n(n-1)/2.
The bialternate sum of a matrix is a special case of the bialternate product (also known as the biproduct), where one of the matrices is twice the identity matrix.
ResourceFunction["BialternateSum"] can be used on SparseArray objects, returning a SparseArray object when possible.

Examples

Basic Examples (1) 

Bialternate sum of a matrix of dimension 2:

In[1]:=
ResourceFunction["BialternateSum"][Array[Subscript[a, ##] &, {2, 2}]]
Out[1]=

Scope (3) 

Bialternate sum of a matrix of dimension 3:

In[2]:=
ResourceFunction["BialternateSum"][ Array[Subscript[b, ##] &, {3, 3}]] // MatrixForm
Out[2]=

Bialternate sum of a matrix of dimension 4:

In[3]:=
ResourceFunction["BialternateSum"][ Array[Subscript[c, ##] &, {4, 4}]] // MatrixForm
Out[3]=

The bialternate sum of a sparse matrix is also a sparse matrix:

In[4]:=
ResourceFunction["BialternateSum"][ SparseArray[{{i_, i_} -> -2, {i_, j_} /; Abs[i - j] == 1 -> 1}, {9, 9}]]
Out[4]=

Applications (3) 

Define the Brusselator system:

In[5]:=
Clear[x, y, \[Alpha], \[Beta]]; Brusselator[{x_, y_}][{\[Alpha]_, \[Beta]_}] := {\[Alpha] - (\[Beta] + 1) x + x^2 y, \[Beta] x - x^2 y}; X = {x, y}; \[Mu] = {\[Alpha], \[Beta]}; Brusselator[X][\[Mu]] // Column
Out[9]=

Generate the Jacobian matrix for the system:

In[10]:=
J[{x_, y_}][{\[Alpha]_, \[Beta]_}] = Simplify[D[Brusselator[X][\[Mu]], {X}]]; J[X][\[Mu]] // MatrixForm
Out[11]=

Solve for the equilibrium points:

In[12]:=
X0[{\[Alpha]_, \[Beta]_}] := SolveValues[Brusselator[X][\[Mu]] == 0, X] Map[Composition[MatrixForm, List], X0[\[Mu]]] // Column
Out[13]=

Generate the linear approximation at the equilibrium point (α,β/α):

In[14]:=
J0[{\[Alpha]_, \[Beta]_}] = Simplify[J[X0[\[Mu]][[1]]][\[Mu]]]; J0[\[Mu]] // MatrixForm
Out[15]=

Evaluate the bialternate sum of the linear approximation at (α,β/α):

In[16]:=
A0[{\[Alpha]_, \[Beta]_}] = Simplify[ResourceFunction["BialternateSum"][J0[\[Mu]]]]; A0[\[Mu]] // MatrixForm
Out[17]=

Generate the local stability criteria depending on the control parameter β:

In[18]:=
Reduce[Det[A0[\[Mu]]] < 0 && \[Alpha] > 0 && \[Beta] > 0, \[Beta]]
Out[18]=

Compute the critical Hopf bifurcation value for the control parameter β:

In[19]:=
\[Beta]0 = Normal@SolveValues[ Det[A0[\[Mu]]] == 0 && \[Alpha] > 0 && \[Beta] > 0, \[Beta]][[1]]
Out[19]=

Define a system representing a nonlinear feedback-control system of Lur'e type:

In[20]:=
Clear[x, y, z, \[Alpha], \[Beta]]; Lure[{x_, y_, z_}][{\[Alpha]_, \[Beta]_}] := {y, z, -x + x^2 - \[Beta] y - \[Alpha] z}; X = {x, y, z}; \[Mu] = {\[Alpha], \[Beta]}; Lure[X][\[Mu]] // Column
Out[24]=

Generate the Jacobian matrix:

In[25]:=
J[{x_, y_, z_}][{\[Alpha]_, \[Beta]_}] = Simplify[D[Lure[X][\[Mu]], {X}]]; J[X][\[Mu]] // MatrixForm
Out[26]=

Solve for the equilibrium points:

In[27]:=
X0[{\[Alpha]_, \[Beta]_}] := SolveValues[Lure[X][\[Mu]] == 0, X] Map[Composition[MatrixForm, List], X0[\[Mu]]] // Column
Out[28]=

Generate the linear approximation at the origin:

In[29]:=
J0[{\[Alpha]_, \[Beta]_}] = Simplify[J[X0[\[Mu]][[1]]][\[Mu]]]; J0[\[Mu]] // MatrixForm
Out[30]=

Evaluate the bialternate sum of the linear approximation at origin:

In[31]:=
A0[{\[Alpha]_, \[Beta]_}] = Simplify[ResourceFunction["BialternateSum"][J0[\[Mu]]]]; A0[\[Mu]] // MatrixForm
Out[32]=

Generate the local stability criteria depending on the control parameter α:

In[33]:=
Reduce[Det[A0[\[Mu]]] < 0 && \[Alpha] > 0 && \[Beta] > 0, \[Alpha]]
Out[33]=

Give the critical Hopf bifurcation value for the control parameter α:

In[34]:=
\[Alpha]0 = SolveValues[Det[A0[\[Mu]]] == 0, \[Alpha]][[1, 1]]
Out[34]=

The Lorenz system:

In[35]:=
Clear[x, y, z, \[Sigma], r, \[Beta]]; Lorenz[{x_, y_, z_}][{\[Sigma]_, r_, \[Beta]_}] := {\[Sigma] (y - x), -x z + r x - y, x y - \[Beta] z}; X = {x, y, z}; \[Mu] = {\[Sigma], r, \[Beta]}; Lorenz[X][\[Mu]] // Column
Out[39]=

Generate the Jacobian matrix:

In[40]:=
J[{x_, y_, z_}][{\[Sigma]_, r_, \[Beta]_}] = Simplify[D[Lorenz[X][\[Mu]], {X}]]; J[X][\[Mu]] // MatrixForm
Out[41]=

Solve for the equilibrium points:

In[42]:=
X0[{\[Sigma]_, r_, \[Beta]_}] := SolveValues[Lorenz[X][\[Mu]] == 0, X] Map[Composition[MatrixForm, List], X0[\[Mu]]] // Column
Out[43]=

Linear approximation at the origin:

In[44]:=
J0[{\[Sigma]_, r_, \[Beta]_}] = Simplify[J[X0[\[Mu]][[1]]][\[Mu]]]; J0[\[Mu]] // MatrixForm
Out[45]=

Evaluate the bialternate sum of the linear approximation at the origin:

In[46]:=
A0[{\[Sigma]_, r_, \[Beta]_}] = Simplify[ResourceFunction["BialternateSum"][J0[\[Mu]]]]; A0[\[Mu]] // MatrixForm
Out[47]=

Generate the local stability criteria depending on the control parameter r:

In[48]:=
Simplify@ Reduce[Det[A0[\[Mu]]] < 0 && \[Sigma] > 0 && r > 1 && \[Beta] > 0, r]
Out[48]=

Give the critical Hopf bifurcation value for the control parameter r:

In[49]:=
r0 = Simplify@ SolveValues[ Det[A0[\[Mu]]] == 0 && \[Sigma] > 0 && r > 1 && \[Beta] > 0, r][[ 1]]
Out[49]=

Properties and Relations (5) 

The bialternate sum of a transpose is equal to the transpose of the bialternate sum:

In[50]:=
n = 7; ResourceFunction["BialternateSum"][ Transpose[Array[\[FormalA], {n, n}]]] == Transpose[ ResourceFunction["BialternateSum"][ Array[\[FormalA], {n, n}]]] // Simplify
Out[51]=

The bialternate sum of a sum of matrices is equal to the sum of the individual bialternate sums:

In[52]:=
n = 7; ResourceFunction["BialternateSum"][ Array[\[FormalA], {n, n}] + Array[\[FormalB], {n, n}]] == ResourceFunction["BialternateSum"][Array[\[FormalA], {n, n}]] + ResourceFunction["BialternateSum"][ Array[\[FormalB], {n, n}]] // Simplify

The bialternate sum of a symmetric matrix is also symmetric:

In[53]:=
n = 7; SymmetricMatrixQ[ResourceFunction["BialternateSum"][HilbertMatrix[n]]]
Out[54]=

A matrix and its eigenvalues:

In[55]:=
m = \!\(\* TagBox[ RowBox[{"(", "", GridBox[{ {"18", "52", RowBox[{"-", "148"}], "148"}, { RowBox[{"-", "10"}], RowBox[{"-", "32"}], "100", RowBox[{"-", "100"}]}, { RowBox[{"-", "13"}], RowBox[{"-", "56"}], "173", RowBox[{"-", "168"}]}, { RowBox[{"-", "11"}], RowBox[{"-", "49"}], "150", RowBox[{"-", "145"}]} }, GridBoxAlignment->{"Columns" -> {{Center}}, "Rows" -> {{Baseline}}}, GridBoxSpacings->{"Columns" -> { Offset[0.27999999999999997`], { Offset[0.7]}, Offset[0.27999999999999997`]}, "Rows" -> { Offset[0.2], { Offset[0.4]}, Offset[0.2]}}], "", ")"}], Function[BoxForm`e$, MatrixForm[BoxForm`e$]]]\); eigs = Eigenvalues[m]
Out[56]=

Compute its bialternate sum:

In[57]:=
MatrixForm[bs = ResourceFunction["BialternateSum"][m]]
Out[57]=

The eigenvalues of the bialternate sum are equal to the sums of the eigenvalues of the original matrix, taken in pairs:

In[58]:=
Eigenvalues[bs]
Out[58]=
In[59]:=
Total[Subsets[eigs, {2}], {2}]
Out[59]=

In analyzing the local stability and determining the critical value of the Hopf bifurcation, the method using the bialternate sum is equivalent to using the Liu minor of the corresponding Hurwitz matrix. Demonstrate this using the Lorenz equations:

In[60]:=
Clear[x, y, z, \[Sigma], r, \[Beta]]; Lorenz[{x_, y_, z_}][{\[Sigma]_, r_, \[Beta]_}] := {\[Sigma] (y - x), -x z + r x - y, x y - \[Beta] z}; X = {x, y, z}; \[Mu] = {\[Sigma], r, \[Beta]}; Lorenz[X][\[Mu]] // Column
Out[61]=

Evaluate the Jacobian matrix:

In[62]:=
J[{x_, y_, z_}][{\[Sigma]_, r_, \[Beta]_}] = Simplify[D[Lorenz[X][\[Mu]], {X}]]; J[X][\[Mu]] // MatrixForm
Out[63]=

Give the equilibrium points:

In[64]:=
X0[{\[Sigma]_, r_, \[Beta]_}] := SolveValues[Lorenz[X][\[Mu]] == 0, X] Map[Composition[MatrixForm, List], X0[\[Mu]]] // Column
Out[65]=

Linear approximation at the equilibrium point :

In[66]:=
J0[{\[Sigma]_, r_, \[Beta]_}] = Simplify[J[X0[\[Mu]][[3]]][\[Mu]]]; J0[\[Mu]] // MatrixForm
Out[67]=

Here is where the two methods start to differ. Evaluate the bialternate sum of the linear approximation at :

In[68]:=
A0[{\[Sigma]_, r_, \[Beta]_}] = Simplify[ResourceFunction["BialternateSum"][J0[\[Mu]]]]; A0[\[Mu]] // MatrixForm
Out[69]=

Generate the local stability criteria depending on the control parameter r:

In[70]:=
Reduce[Det[A0[\[Mu]]] < 0 && \[Sigma] > 0 && r > 1 && \[Beta] > 0, r]
Out[70]=

Give the critical Hopf bifurcation value for the control parameter r:

In[71]:=
r0 = Simplify@ SolveValues[ Det[A0[\[Mu]]] == 0 && \[Sigma] > 0 && r > 1 && \[Beta] > 0, r][[ 1]]
Out[71]=

For the method based on the Liu minor, start by computing the characteristic polynomial of the linear approximation at :

In[72]:=
polJ0 = -Collect[ CharacteristicPolynomial[J0[\[Mu]], \[Lambda]], \[Lambda], Simplify]
Out[72]=

Generate the Hurwitz matrix for the characteristic polynomial:

In[73]:=
HMlorenz = ResourceFunction["HurwitzMatrix"][polJ0, \[Lambda]]; HMlorenz // MatrixForm
Out[74]=

Generate the Liu minor:

In[75]:=
MatrixForm[ LiuSubmatrix = Take[HMlorenz, Length[HMlorenz] - 1, Length[HMlorenz] - 1]]
Out[76]=
In[77]:=
LiuMinor = Simplify[Det[LiuSubmatrix]]
Out[78]=

Show the equivalence between the two methods:

In[79]:=
Reduce[Det[A0[\[Mu]]] < 0 && \[Sigma] > 0 && r > 1 && \[Beta] > 0, r] === Reduce[ Reduce[Thread[ DeleteCases[CoefficientList[polJ0, \[Lambda]], _Integer] > 0]] && LiuMinor > 0 && \[Sigma] > 0 && r > 1 && \[Beta] > 0, r]
Out[79]=
In[80]:=
Simplify@ SolveValues[ Det[A0[\[Mu]]] == 0 && \[Sigma] > 0 && r > 1 && \[Beta] > 0, r][[ 1]] === Simplify@ SolveValues[LiuMinor == 0 && \[Sigma] > 0 && r > 1 && \[Beta] > 0, r][[1]]
Out[80]=

Neat Examples (1) 

Visualize the structure of a bialternate sum:

In[81]:=
MatrixPlot[ResourceFunction["BialternateSum"][ToeplitzMatrix[16]]]
Out[81]=

Publisher

Ramón Eduardo Chan López

Version History

  • 1.0.0 – 28 February 2022

Source Metadata

  • Citation:
    • Asada, T., Ishikawa, T., Time and Space in Economics. Springer Science & Business Media, 2007.
    • Fuller, A. T., "Conditions for a matrix to have only characteristic roots with negative real parts." Journal of Mathematical Analysis and Applications, vol. 23, no. 1, 71–98, 1968. DOI: 10.1016/0022-247X(68)90116-9
    • Govaerts, W. J., Numerical methods for bifurcations of dynamical equilibria. Society for Industrial and Applied Mathematics, 2000.
    • Guckenheimer, J., Myers, M., Sturmfels, B., "Computing Hopf bifurcations I." SIAM Journal on Numerical Analysis, vol. 34, no. 1, 1–21, 1997. DOI: 10.1137/S0036142993253461
    • Jury, E., Gutman, S., "On the stability of the A matrix inside the unit circle." IEEE Transactions on Automatic Control, vol. 20, no. 4, 533–535, 1975. DOI: 10.1109/TAC.1975.1100995
    • Kuznetsov, Y. A., Elements of Applied Bifurcation Theory. Applied Mathematical Sciences, Vol. 112. New York: Springer Science & Business Media, 2013.
    • Liu, W. M., "Criterion of Hopf bifurcations without using eigenvalues." Journal of Mathematical Analysis and Applications, vol. 182, no. 1, 250–256, 1994. DOI: 10.1006/jmaa.1994.1079
    • Stéphanos, C. "Sur une extension du calcul des substitutions linéaires." Journal de Mathématiques pures et Appliquées, vol. 6, 73–128, 1900. http://sites.mathdoc.fr/JMPA/PDF/JMPA_1900_5_6_A4_0.pdf

Related Resources

License Information