%Simple Matlab interface for Moransys;
%ypliu99 at mails.tsinghua.edu.cn  or ypliu26 at gmail.com;
%Electronics Engineering Department, Tsinghua University

%once computation use only;
[A,rows,cols,entries]=MOR2MATLAB('mor.A');
[B,rows,cols,entries]=MOR2MATLAB('mor.B');
[E,rows,cols,entries]=MOR2MATLAB('mor.E');
[C,rows,cols,entries]=MOR2MATLAB('mor.C');
[fix,rows,cols,entries]=MOR2MATLAB('fixstep.txt');

%global matrix transform;
Ein = inv(E);
AA = Ein*A;
BB = Ein*B;
global Am;
global Bm;
Am = AA;
Bm = BB;

%transient analysis;
y0=zeros(30,1); %initial temperature; 
atol=ones(30,1)*1e-5;
options = odeset('RelTol',1e-3,'AbsTol',atol);
%[t,Y] = ode45(@thrusys,[0 20],y0,options); %very slow;
[t,Y] = ode15s(@thrusys,[0 20],y0,options); %very fast;

%plot figure;
T=C*Y'; %turn down the output dimension by ode15 call; 
To = fix;
to = To(1,:);
figure(1);
%plot the results from reduce models
plot(t,T(1,:)','-g',t,T(2,:),'-.g',t,T(3,:),'.g',t,T(4,:),'--g');
hold on;
%plot the fix fixstep results from Ansys.
plot(to,To(2,:)','-r',to,To(3,:),'-.r',to,To(4,:),'.r',to,To(5,:),'--r');