# NSolve/FindRoot matrix equation and conditions

I have an eigenvalue (vector) problem (please see the theory, maybe I missed something in the post or did not understand Eqs (3.5.24) ) how to obtain two unknown matrices (36 unknowns) to satisfy all given conditions. Unknown matrices are

 XX = Table[X[10 i + j], {i, 6}, {j, 3}];
YY = Table[Y[10 i + j], {i, 6}, {j, 3}];


which should be joined and present the 36 unknowns

 unknown = ArrayFlatten[({{XX, YY}})];


Let's suppose that we have two matrices in the following form:

JJ = {{0., 0., 0., 1., 0., 0.}, {0., 0., 0., 0., 1.,
0.}, {0., 0., 0., 0., 0., 1.}, {-1., 0., 0., 0., 0.,
0.}, {0., -1., 0., 0., 0., 0.}, {0., 0., -1., 0., 0.,
0.}};

HH={{5.33449*10^7, -3.14159*10^7, 3.99684*10^6, 0., 199.842,
0.}, {-3.14159*10^7, 2.82358*10^7, -200000., -200.,
0., -199.842}, {3.99684*10^6, -200000., 4.23581*10^6, 0., 199.842,
0.}, {0., -200., 0., 0.01, 0., 0.}, {199.842, 0., 199.842, 0.,
0.00999211, 0.}, {0., -199.842, 0., 0., 0., 0.00999211}};


We obtained eigenvalues s1, s2 and s3, (s4=-s1, s5=-s2, s6=-s3) and take theirs imaginary parts

q0 = Eigenvalues[N[{HH, JJ}, 5]];

s1 = N[Im[q0]][[5]];
s2 = N[Im[q0]][[3]];
s3 = N[Im[q0]][[1]];


to construct the diagonal matrix s in the form

s = {{s1, 0, 0}, {0, s2, 0}, {0, 0, s3}};


Now, we need to satisfy (find unknowns) with NSolve or FindRoot all conditions where

 zero = Table[0, {i, 3}, {j, 3}];
onem = IdentityMatrix[3];

one = Transpose[XX].JJ.XX;
two = Transpose[YY].JJ.YY;
three = Transpose[YY].HH.XX;
four = Transpose[XX].HH.YY;

five = Transpose[XX].HH.XX;
six = Transpose[YY].HH.YY;
seven = Transpose[XX].JJ.YY;
eight = -Transpose[YY].JJ.XX;

FindRoot[{one == zero, two == zero, three == zero, four == zero,
five == s, six == s, seven == onem,
eight == onem}, {{X[11], 1}, {X[12], 1}, {X[13], 1}, {X[21],
1}, {X[22], 1}, {X[23], 1}, {X[31], 1}, {X[32], 1}, {X[33],
1}, {X[41], 1}, {X[42], 1}, {X[43], 1}, {X[51], 1}, {X[52],
1}, {X[53], 1}, {X[61], 1}, {X[62], 1}, {X[63], 1}, {Y[11],
1}, {Y[12], 1}, {Y[13], 1}, {Y[21], 1}, {Y[22], 1}, {Y[23],
1}, {Y[31], 1}, {Y[32], 1}, {Y[33], 1}, {Y[41], 1}, {Y[42],
1}, {Y[43], 1}, {Y[51], 1}, {Y[52], 1}, {Y[53], 1}, {Y[61],
1}, {Y[62], 1}, {Y[63], 1}}]


but we have more equations then unknowns. If I fill use only first four equations, I can get the solution, but how to satisfy the rest number of equations?

Also, I choose roots around 1. How to get one set of solutions to satisfy all 8 matrix equations?

• FindInstance might be helpful. Anyways, if you end up with more equations than unknowns, then you might have redundancies in your equations or your problem is not well-posed in the first place. Are you sure that your problem should have solutions? Can you tell in advance how many solutions there might be? Commented May 15, 2018 at 9:29
• there are a more solutions... it is posted well for sure, these conditions are simple eigen vector problem and they are well known as a conditions of symplectic matrix... if you know a way to find only one solution, please write Commented May 15, 2018 at 10:38
• If it is really only an eigenvector problem, then it should be solvable with Eigensystem. You need to Flatten your matrix variables and rephrase left and right multiplication by other matices as linear mapping on the flattened matrix variables. Commented May 15, 2018 at 10:41
• thanks for the comment, can you please check what you wrote because i am not sure that you explained well Commented May 15, 2018 at 10:47
• Before I start with that, would you please be so kind and write down the mathematica eigenvalue problem at the beginning your question (in $\LaTeX$ if possible)? I am not sure if I understood you correctly. I might write an extended answer afterwards. Commented May 15, 2018 at 10:49

The book's author's primary goal seems to be his self-entertainment.

He wants to solve the first order ODE

$$\dot{\mathbf{u}}(t) = \mathbf{J} \, \mathbf{H} \, \mathbf{u}(t)$$

and applies an exponential ansatz

$$\mathbf{u} = \sum_{r=1}^{2 n} \hat{\mathbf{u}}_r \operatorname{e}^{\lambda_r \,t + c_r}$$

with integration constants $c_r$ and solutions $(\hat{\mathbf{u}}_r,\lambda_r)$ of the eigenvalue problem

$$\mathbf{J} \, \mathbf{H} \, \mathbf{u}_r = \lambda_r \, \mathbf{u}_r.$$

So. Let's. Just. Solve. That.

n = 3;
JJ = ArrayFlatten[{{ConstantArray[0., {n, n}],
N@IdentityMatrix[n]}, {-N@IdentityMatrix[n],
ConstantArray[0., {n, n}]}}];
HH = N[{{4.93480220054467930941724555.*^7, -3.1415926535897937*^7,
0, 0, 0, 0}, {-3.1415926535897937*^7,
2.0238963636413604*^7, -200000.5., 0, 0, 0}, {0, -200000.5.,
238963.63641360097489457615., 0, 0, 0}, {0, 0, 0, 0.015., 0,
0}, {0, 0, 0, 0, 0.00999211054574252555., 0}, {0, 0, 0, 0, 0,
0.00999211054574252555.}}];

{λ, vecs} = Eigensystem[JJ.HH];


Bamm. That's it.

Let's synthesize the general solution $\mathbf{u}(t)$ of the ODE and test it:

u = t \[Function] Evaluate[Simplify[Total[Exp[t λ + Array[C, Length[λ]]] vecs]]];
u'[t] - JJ.HH.u[t] // Simplify // Chop


{0, 0, 0, 0, 0, 0}

Done.

Afterwards, the author means to say (but he did't):

$$\mathbf{x}_r = \Re(\hat{\mathbf{u}}_e) \quad \text{and} \quad \mathbf{y}_r = \Im(\hat{\mathbf{u}}_r).$$

In Mathematica terms, this is

x = Re[vecs];
y = Im[vecs];


Afterwards, he collects a bunch of equations that are automatically satisfied. Here are some of them (3.5.12 and 3.5.13):

ppos = Flatten[Position[λ/I // Chop, _?Positive]];
npos = Flatten[Position[λ/I // Chop, _?Negative]];
ω = 0.5 (λ[[ppos]]/I - λ[[npos]]/I) // Chop;

Table[HH.x[[ppos[[i]]]] - (ω[[i]] JJ.y[[ppos[[i]]]]), {i, 1, n}] // Chop
Table[HH.x[[npos[[i]]]] - (-ω[[i]] JJ.y[[npos[[i]]]]), {i, 1, n}] // Chop


{{0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0}}

{{0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0}}

The matrices $\mathbf{X}$ and $\mathbf{Y}$ can be obtained from that by

XX = Transpose[x[[1 ;; n]]];
YY = Transpose[y[[1 ;; n]]];


I leave checking all the other conditions to you.

• Dear Henrik, thank you for your post, but the XX and YY matrices should be 3x6 that T=(XX YY) be 2nx2n (6x6) so that the conditions 3.5.25 and 3.5.26 should be satisfied, but there are not. Commented May 15, 2018 at 22:59
• try with this matrix HH which I provided, because this HH matrix satisfy condition to have real and imaginary part of eigenvector. With code up something must be different that we can satisfy conditions 3.5.25 and 3.5.26 Commented May 15, 2018 at 23:02
• for sure XX and YY are not good Commented May 15, 2018 at 23:10
• Henrik can you please just extract XX and YY which I need? Commented May 16, 2018 at 2:01
• Henrik can you please take a look on Eq. 3.5.12. It is written that H (6x6) matrix multiple vector xr(3x1). This is strange I can not form the T matrix. I determined the eigenvector what you wrote, it is simple. Matrix T is a point of the question Commented May 17, 2018 at 0:18