How to solution of linear-time varying system correctly

Question

I tried the given below to get the solution of the linear-time varying equation using for-loop, but the solution is not correct as I want to. I don’t know where I am getting wrong while implementing this.

r0 = 0.05;
L = 0.1;
d = 0.005;
w0 = 1.5;

Ts=10;
t=[0:Ts:800];
x0=[0 0 1]';
y0=1;
x_value=[];

for k=1:(length(t)) % Number of Iterations
x_value=[x_value x0];
g(k) =  (2*r0*L*sinh((d*(k))/2))./(d*cosh((d*(k))/2)+L*sinh((d*(k))/2));
    A0 = [ -0.5*g(k)   -w0    0 ; 
       w0     -g(k)*0.5   0 ; 
       0        0     -g(k)];
   B0 =[0; 0; -g(k)];
Q=integral(@(u) expm(A0*((k+1)*Ts)-u)*B0,(k*Ts),((k+1)*Ts), 'ArrayValued', true);
x0=expm(A0*Ts)*x0+Q;
end

plot(t,x_value,'r-','linewidth',1);