diff --git a/MatlabFiles/MSV/fwz_msv_msre.m b/MatlabFiles/MSV/fwz_msv_msre.m index 79ffb6b25c74f535404e0b92889cab710a58c2ac..5b97593340fcbd8b010922b4b22947f4172d3817 100755 --- a/MatlabFiles/MSV/fwz_msv_msre.m +++ b/MatlabFiles/MSV/fwz_msv_msre.m @@ -17,8 +17,10 @@ function [F1, F2, G1, G2, V, err] = fwz_msv_msre(P, A, B, Psi, s, x, max_count, % % eta(t) = F2{s(t)}*x(t-1) + G2{s(t)}*epsilont(t) % -% err is set to 0 upon sucess. A positive value is the number of -% iterations before the method was terminated. +% A positive value of err is the number of iterations needed to obtain +% convergence and indicates success. A negitive value of err is the number of +% iterations before the method terminated without convergence and indicates +% failure. % h=size(P,1); @@ -38,20 +40,9 @@ if nargin <= 5 end f=zeros(h*s*(n-s),1); - Is=eye(s); I=eye(n-s); -% The following code would be used if Pi were passed instead of the -% assumption that Pi' = [zeros(s,n-s) eye(s)]. -% for i=1:h -% [Q,R] = qr(Pi{i}); -% W=[zeros(n-s,s) I; inv(R(1:s,1:s)); zeros(s,n-s)]*Q'; -% A{i}=W*A{i}; -% B{i}=W*B{i}; -% Psi{i}=W*Psi{i}; -% end - U=cell(h,1); for i=1:h U{i}=inv(A{i}); @@ -109,9 +100,9 @@ while cont end if (norm(f) < tol) - err=0; -else err=count; +else + err=-count; end F1=cell(h,1); diff --git a/MatlabFiles/MSV/fwz_pi2eye.m b/MatlabFiles/MSV/fwz_pi2eye.m new file mode 100644 index 0000000000000000000000000000000000000000..eb58610077f5080b39ea16f8ed8f19352b3fff68 --- /dev/null +++ b/MatlabFiles/MSV/fwz_pi2eye.m @@ -0,0 +1,20 @@ +function [A_new, B_new, Psi_new] = fwz_pi2eye(A, B, Psi, Pi) +% [A, B, Psi] = fwz_pi2eye(A, B, Psi, Pi) +% Converts the representation +% +% A(s(t) x(t) = B(s(t) x(t-1) + Psi(s(t)) epsilon(t) + Pi{s(t)} eta(t) +% +% to a form compatible with fwz_msv_msre(). This is the form +% +% A(s(t) x(t) = B(s(t) x(t-1) + Psi(s(t)) epsilon(t) + Pi eta(t) +% +% where Pi' = [zeros(s,n-s) eye(s)] +% + +for i=1:h + [Q,R] = qr(Pi{i}); + W=[zeros(n-s,s) I; inv(R(1:s,1:s)); zeros(s,n-s)]*Q'; + A_new{i}=W*A{i}; + B_new{i}=W*B{i}; + Psi_new{i}=W*Psi{i}; +end \ No newline at end of file diff --git a/MatlabFiles/MSV/verify_solution.m b/MatlabFiles/MSV/verify_solution.m new file mode 100644 index 0000000000000000000000000000000000000000..dc0737f24cd88bcdcbeeda0406e6386423308f4c --- /dev/null +++ b/MatlabFiles/MSV/verify_solution.m @@ -0,0 +1,52 @@ +function diff = verify_solution(P,A,B,Psi,s,F1,F2,G1,G2,V) +% +% Verifies that +% +% x(t) = V(s(t))(F1(s(t))x(t-1) + G1(s(t))epsilon(t) +% eta(t) = F2(s(t))x(t-1) + G2(s(t))epsilon(t) +% +% is a solution of +% +% A(s(t)) x(t) = B(s(t)) x(t-1) + Psi(s(t) epsilon(t) + Pi eta(t) +% +% where Pi' = [zeros(s,n-s) eye(s)] by verifying that +% +% [A(j)*V(j) Pi] [ F1(j) +% F2(j) ] = B(j) +% +% [A(j)*V(j) Pi] [ G1(j) +% G2(j) ] = Psi(j) +% +% (P(i,1)*F2(1) + ... + P(i,h)*F2(h))*V(i) = 0 +% + +h=size(P,1); +n=size(A{1},1); +r=size(Psi{1},2); + +diff=0; +Pi=[zeros(n-s,s); eye(s)]; +for j=1:h + X=inv([A{j}*V{j} Pi]); + + tmp=norm(X*B{j} - [F1{j}; F2{j}]); + if tmp > diff + diff=tmp; + end + + tmp=norm(X*Psi{j} - [G1{j}; G2{j}]); + if (tmp > diff) + diff=tmp; + end +end + +for i=1:h + X=zeros(s,n); + for j=1:h + X=X+P(i,j)*F2{j}; + end + tmp=norm(X*V{i}); + if tmp > diff + diff=tmp; + end; +end