% Test script for ltpdamodel/inv
% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% HISTORY:     02-10-2008 L Ferraioli
%                 Creation
% 
% $Id: test_ltpdamodel_inv.m,v 1.2 2009/05/25 15:11:33 luigi Exp $
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Build numeric test matrix

a=1;b=2;c=3;
d=5;e=2;f=6;
g=1;h=3;i=5;
M=[a,b,c;d,e,f;g,h,i];
iM = inv(M);

%% Equivalent ltpdamodel matrix

A = [ltpdamodel(plist('expression','a','params',{'a'},'values',{1},'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m'))),...
  ltpdamodel(plist('expression','b','params','b','values',2,'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m'))),...
  ltpdamodel(plist('expression','c','params',{'c'},'values',{3},'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m')));...
 ltpdamodel(plist('expression','d','params',{'d'},'values',{5},'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m'))),...
  ltpdamodel(plist('expression','e','params','e','values',2,'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m'))),...
  ltpdamodel(plist('expression','f','params',{'f'},'values',{6},'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m')));...
  ltpdamodel(plist('expression','g','params',{'g'},'values',{1},'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m'))),...
  ltpdamodel(plist('expression','h','params','h','values',3,'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m'))),...
  ltpdamodel(plist('expression','i','params',{'i'},'values',{5},'xvar','x','xvals',[1:10].','xunits','Hz','yunits',unit('m')))];

%% Do inverse and eval

iA = inv(A);

%% Numeric eval

imat = zeros(3,3);
for ii=1:3
  for jj=1:3
    imat(ii,jj)=eval(iA(ii,jj).expr.s);
  end
end

imat-iM

%% Symbolic test

B = A*iA;

%% Numeric eval

imat2 = zeros(3,3);
for ii=1:3
  for jj=1:3
    imat2(ii,jj)=eval(B(ii,jj).expr.s);
  end
end