# FullInformationOutputRegulator: Unable to determine the composite center manifold

Suppose we have the following example taken by Isidori's book on Nonlinear Control Systems. I've asked some days ago here how the function works. Back then I have been addressed to this very book for the understanding of the function FullInformationOutputRegulator implemented in Mathematica.

h = x1 - w1;
f = {x2, 0, x3 + x1 + x2^2};
x = {x1, x2, x3};
g = {{0}, {0}, {1}};
asys = AffineStateSpaceModel[{f, g, {h}}, x];


Where the exosystem has been defined as:

tmodel = NonlinearStateSpaceModel[{{w2, -w1}, {}}, {w1, w2}, {}];
sys = SystemsModelMerge[{asys, tmodel}];


Following an example in the Documentation I've merged the two systems, yet Mathematica won't find the feedback law.

fb = FullInformationOutputRegulator[
sys, {"Poles", {-1, -1.5 + I, -1.5 - I}}];


The error is the following:

FullInformationOutputRegulator::csm Unable to determine the composite center manifold.


Would someone kindly help me out? Have I misunderstood how to introduce the exosystem in the equations?

Assuming I've correctly written my problem, it could be a bug.

EDIT

I've calculated the feedback law following Isidori's exercise.

    Clear[l];
h = x1 - w1;
f = {x2, 0, x3 + x1 + x2^2};
x = {x1, x2, x3};
g = {{0}, {1}, {0}};
s = {a w2, -a w1};
w = {w1, w2};
k1 = w1;
k2 = DLie[k1, s, w, 1];
deqn = D[l[w1, w2], w1] a w2 - D[l[w1, w2], w2] a w1 ==
l[w1, w2] + w1 + a^2 w2^2;
aSoln = l -> (a1 #1 + a2 #2 + a11 #1^2 + a12 #1 #2 + a22 #2^2 &);
eqn = deqn /. aSoln;
coefEqn = Thread[(CoefficientList[#, {w1, w2}] // Flatten) & /@ eqn];
coefSoln =
Solve[coefEqn, {a1, a2, a11, a12, a22}, Reals][] // Simplify;
soln = aSoln /. coefSoln;
l = l[w1, w2] /. soln;
pi = {{k1}, {k2}, {l}};
c = DLie[k2, s, w, 1];
k = Ackermann[f, g u, h, x][];
Clear[t]
alpha = c + k.(x - pi) /. a -> 1 /. x1 -> x1[t] /. x2 -> x2[t] /.
x3 -> x3[t] /. w1 -> w1[t] /. w2 -> w2[t];


So u=alpha. I have defined DLie here. As for Ackermann it is a simple function:

Ackermann[f_, g_, h_, x_] := With[{},
matrixA =
DJacobian[f, x] /. Table[x[[i]] -> 0, {i, 1, Dimensions[x][]}];
vectorB =
D[g, u] /. Table[x[[i]] -> 0, {i, 1, Dimensions[x][]}];
vectorC = Table[D[h, x[[i]]], {i, 1, Dimensions[x][]}];
sys = StateSpaceModel[{matrixA, vectorB, {vectorC}}];
vectorK =
StateFeedbackGains[sys, {-40, -0.4, -0.4}, Method -> "Ackermann"]]


A possible reason or bug could be in the calculation of l(w1,w2). Indeed as pointed out here Mathematica finds some trouble solving the proper partial differential equation with DSolve.

• One of the assumptions for the theory to work is that the exosystem be neutrally stable, i.e. some poles on the imaginary axis and any remaining must be in the r.h.p. There should probably be an example in possible issues highlighting this assumption. Nov 9 '16 at 17:32
• You are right. I wrote wrong the exosystem, missed a minus. Yet it still doesn't work. @SubaThomas Nov 9 '16 at 17:38
• @SubaThomas Might it be possible that Mathematica has had trouble solving the center manifold equations because it doesn't have a unique solution? I think I'm going to write a possible answer to this. Nov 9 '16 at 18:13
• This needs to be investigated further. Probably a bug. Nov 9 '16 at 18:24
• @SubaThomas I have updated my question adding the code I wrote for the calculation of the feedback law, following Isidori's path. Nov 9 '16 at 19:03

This is not a bug.

With the updated equations it works.

h = x1 - w1;
f = {x2, 0, x3 + x1 + x2^2};
x = {x1, x2, x3};
g = {{0}, {1}, {0}};

sys = SystemsModelMerge[{AffineStateSpaceModel[{f, g, {h}}, x],
NonlinearStateSpaceModel[{{ w2, - w1}, {}}, {w1, w2}, {}]}];

FullInformationOutputRegulator[sys, {"Poles", {-1, -2, -3}}]


{w1 + 18 (-w1 + x1) + 7 (-w2 + x2) + 24 (w1/2 + (2 w1^2)/5 + w2/2 - (2 w1 w2)/5 + (3 w2^2)/5 + x3)}

However, when $g$ is $\left( \begin{array}{ccc} 0 & 0 & 1 \\ \end{array} \right)^T$ as you originally have it cannot be solved because the conditions for regulation cannot be satisfied even for the linear case. (Refer to Isidori for the equations and notation.)

g = {{0}, {0}, {1}};
w = {w1, w2};
A = D[f, {x}];
P = D[f, {w}];
B = g;
S = D[{w2, -w1}, {w}];
CC = D[{h}, {x}];
Q = D[{h}, {{w1, w2}}];
\[CapitalPi] = Array[Subscript[p, ##] &, {3, 2}];
\[CapitalGamma] = Array[Subscript[m, ##] &, {1, 2}];