How do I get this code to converge?

% This code runs, but theta_hat should converge to 1, and e which is the
% error should converge to 0.


clc; 
clear all; 
close all; 

T=120;
%Simulation step size 
dt = 0.0001; 
t(1)=0;
N=T/dt;

% Initial conditions
psi = [1/sqrt(2) 1/sqrt(2)]';
rho(:,:,1) = psi*psi';          
psi_hat = [1/sqrt(5) 1/sqrt(5)]';
rho_hat(:,:,1) = psi_hat*psi_hat';
y(1) = 0;
y_hat(1) =  0;
theta_hat(1) = 1.5;  

% Parameters definitions
A = 1;
w = 1;
omega=1;
theta = 1;
gamma = 0.1; 
GAMMA = 1;
mu = [0 1;1 0];
P = [1 0;0 0];
H = [w/2 0;0 -w/2];


rho11(1)=rho(1,1,1);
rho12(1)=rho(1,2,1);
rho21(1)=rho(2,1,1);
rho22(1)=rho(2,2,1);

rho_hat11(1)=rho_hat(1,1,1);
rho_hat12(1)=rho_hat(1,2,1);
rho_hat21(1)=rho_hat(2,1,1);
rho_hat22(1)=rho_hat(2,2,1);

for k=1:N
    t(k+1)=t(k)+dt;
    
    u=A*cos(omega*t(k));
    rho(:,:,k+1)=rho(:,:,k)+dt*(-1i)*((H+u*theta*mu)*rho(:,:,k)-rho(:,:,k)*(H+u*theta*mu));
       rho11(k+1)=rho(1,1,k);
       rho12(k+1)=rho(1,2,k);
       rho21(k+1)=rho(2,1,k);
       rho22(k+1)=rho(2,2,k);
       

  y(k+1) = trace(P*rho(:,:,k+1));

% Solving for rho_hat %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

  drho_hat_dt = -1i*((H+u*theta*mu)*rho_hat(:,:,k) - rho_hat(:,:,k)*(H+u*theta*mu))*GAMMA*(y(k)-y_hat(k))*(P*rho_hat(:,:,k) ...
  +rho_hat(:,:,k)*P - 2*trace(P*rho_hat(:,:,k))*rho_hat(:,:,k));

    rho_hat(:,:,k+1)=rho_hat(:,:,k)+dt*drho_hat_dt ;
    
% extracting matrix elements
       rho_hat11(k+1)=rho_hat(1,1,k);
       rho_hat12(k+1)=rho_hat(1,2,k);
       rho_hat21(k+1)=rho_hat(2,1,k);
       rho_hat22(k+1)=rho_hat(2,2,k);

% Solve for y_hat   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
    y_hat(k+1) = trace(P*rho_hat(:,:,k));

% solving for theta_hat
       
   dtheta_hat_dt = -1i*gamma*u*trace(P*(mu*rho_hat(:,:,k)-rho_hat(:,:,k)*mu))*(y(k)-y_hat(k));  
   theta_hat(k+1) = theta_hat(k)  +  dt*dtheta_hat_dt;
   
% Evaluating the error

   e(k+1) = y(k) - y_hat(k);       
       
end

subplot(223)
plot(t,abs(theta_hat)); 
grid on; 
xlabel('Time [sec]','fontsize',20), 
ylabel('$\hat{\theta}$','interpreter','latex','fontsize',20);
title("Plots",'fontsize',16);

subplot(224)
plot(t,abs(e)); 
grid on; 
xlabel('Time [sec]','fontsize',20), 
ylabel('$e$','interpreter','latex','fontsize',20);
title("Plots",'fontsize',16);

Answer 1

[Vaibhav]

Hi Kaydian

The issue with the convergence of "theta_hat" to 1 and e to 0 in your code may be due to several factors, including the learning rate "gamma", the initial conditions, the model used for the system and observer, or numerical stability.

Here are a few suggestions that might help improve convergence:

1. Check the Adaptation Gain (gamma): The value of "gamma" significantly affects the speed of convergence and stability. If "gamma" is too small, convergence will be slow, and if it's too large, the system might become unstable or oscillate.
2. Initial Conditions: The choice of initial conditions, especially for "theta_hat", can affect the convergence. You've already set theta_hat(1) = 1.5;, which is close to the true value of theta = 1. However, depending on the system dynamics and adaptation law, starting closer or exactly at the true value might not always be necessary or beneficial. Experiment with different initial conditions for "theta_hat" and "psi_hat".
3. Numerical Stability: The choice of "dt" (time step) can impact numerical stability. While your dt = 0.0001 seems small, ensuring numerical stability involves more than just the time step size; it also includes how the differential equations are solved. You might want to experiment with slightly larger "dt" values to see if it affects convergence without sacrificing too much precision.

Hope this helps!