```Date: Feb 22, 2013 4:10 AM
Author: Milos Milenkovic
Subject: Automating an iterative procedure

%Dear,%I need help with automating the following iterative procedure:%i=1,...,N; N is number of iterations, it depends on the convergence of%B0[i] and B1[i] but let assume for now  that it is given.    %Input parameters: K0=[2.337 48.551 0 2.337; 48.551 3072.487 0 48.551; 0 0 0 0; 2.337 48.551 0 2.337];K1=[-0.2582 -18.596 0 -0.2582; -34.737 -1300.188 0 -34.737; 0 0 0 0; -0.2582 -18.596 0 -0.2582]; %for i=1: B0[i=1]=[-1 0 0 0;  0 -1 0 0 ; 0 0 -1 0; 0 0 0 -1];    B1[i=1]=[0.01 0.02 0.07 0.03; 0.02 0.05 0.06 0.02; 0.07 0.09 0.01 0.03; 0.01 0.03 0.02 0.04];%Unknown parameters:%for i=1:N =>Ke0[i]; Ke1[i]; Ke[i];B0[i+1]=B0[i]+Delta0[i];B1[i+1]=B1[i]+Delta1[i];%Iteration i: condB1[i] = cond(B1[i]);condB0[i] = cond(B0[i]);condB1[i] = 1.4662e+017;condB0[i] = 1;Ke1[i] = pinv(B1[i])*K1/B0[i]';Ke0[i]=dlyap(B0[i],K0,[],B1[i]);Ke[i]=0.5*[Ke1[i]+Ke0[i]);A1[i]=B1[i]*Ke[i]*B1[i]'; A2[i]=Ke[i]*B1[i]'; A3[i]=B1[i]*Ke[i];A4[i]=K1-A2[i];n=length(A2[i]);P1[i] = zeros(n);P1[i](1:end)=1:n^2;Q1[i] = kron(speye(n),A3[i]);Q1[i] = Q1[i](:,P1[i]');M1[i] = kron(A2[i]',speye(n)) + Q1[i];Delta1[i]= M1[i] \ A4[i](:);Delta1[i] = reshape(Delta1[i],[n n]);norm( Delta1[i]*A2[i]+A3[i]*Delta1[i]'-A4[i]) / norm(A4[i]);A5[i]=B0[i]*Ke[i]*B0[i]';A6[i]=B0[i]*Ke[i];A7[i]=Ke[i]*B0[i]';A8[i]=K0-A5[i]-A1[i]-A3[i]*Delta1[i]'-Delta1[i]*A2[i];n=length(A7[i]);P0[i] = zeros(n);P0[i](1:end)=1:n^2;Q0[i]= kron(speye(n),A6[i]);Q0[i] = Q0[i](:,P0[i]');M0[i] = kron(A7[i]',speye(n)) + Q0[i];Delta0[i] = M0[i] \ A8[i](:);Delta0[i] = reshape(Delta0[i],[n n]);norm( Delta0[i]*A7[i]+A6[i]*Delta0[i]'-A8[i]) / norm(A8[i]);%for the next iteration     B0[i+1]=B0[i]+Delta0[i];    B1[i+1]=B1[i]+Delta1[i];    and the iteration again.     %Thanks!!
```