How do I solve the Riccati Differential Equation with fsolve in MATLAB / Octave?

Question

Let's say that I have this equation. This is the Riccati Differential Equation:

$$P A + A^T P - P B R^{-1} B^T P + Q = 0$$

I know $ A, B, Q, R$. My goal is to find $P$. How do I use the MATLAB / Octave command fsolve to find $P$ ?

EDIT: $$ $$

Here is the answer - An example to get LQR matrix.

>> fun = @(P, A, B, Q, R) P*A+A'*P-P*B*inv(R)*B'*P + Q
fun =

@(P, A, B, Q, R) P * A + A' * P - P * B * inv (R) * B' * P + Q
>> p0 = 0.1*ones(2)
p0 =

   0.10000   0.10000
   0.10000   0.10000

>> P = fsolve(@(P) fun(P, Amat, Bmat, Q, R), p0)
P =

   55.67637    0.57241
    0.57241    2.29935

>>

>> inv(R)*Bmat'*P
ans =

   2.2269e-02   2.2895e-04
  -1.4350e-03  -5.7645e-03

>> lqr(Amat, Bmat, Q, R)
ans =

   2.2269e-02   2.2899e-04
  -1.4353e-03  -5.7644e-03

>>

Answer 1

Define a function, demofun:

function err = demofun(p, A, B, Q, R)
P = reshape(p, [3 3]); %% Reshape p to your desired matrix size
F = P*A+A'*P-P*B*inv(R)*B'*P + Q;

% Calculate the error value
err = sum(F(:).*F(:));

return

Then write the script that optimizes the function

A = randn(3);
B = randn(3);
Q = randn(3);
R = randn(3);

p0 = ones(9, 1)*0.1; % Initial guess for P, vectorized to p

tmp = fsolve(@(p) demofun(P, A, B, Q, R), p0);

P = reshape(tmp, [3 3]); % Reshape vectorized answer to matrix format