function [t,y,y_theta,x,x_last] = mtt_stime(system_name,theta,free);
## usage: [t,y,y_theta,x,x_last] = mtt_stime(system_name,theta);
##
## Simulate system with name system_name and parameter vector theta
## The order of components in theta is determined in system_numpar.txt:
## y_theta contains the corresponding sensitivity functions
## Assumes system generated by the sBG approach
## Copyright (C) 1999 by Peter J. Gawthrop
## $Id$
## Simulate using mtt-generated function
N = length(theta);
eval(sprintf("[nx,ny,nu,nz,nyz] = %s_def;", system_name));
if nargin<3
free = 1;
endif
y_theta = [];
if length(free)==0
free=1; # Make the loop happen once to get y and X
endif
[n,m] = size(free);
if m==1
free = free';
endif
for i=free
args=sprintf("%i",i);
for j=1:length(theta)
args = sprintf("%s %g",args, theta(j));
endfor
## Run system and replace NaN by 1e30 -- easier to handle
file_name = sprintf("%s_input.dat", system_name);
if exist(file_name)==2 # Then use data from this file ...
command = sprintf("./%s_ode2odes.out %s < %s | sed \'s/NAN/Inf/g\' >mtt_out_data.dat\n", ...
system_name, args, file_name);
else
command = sprintf("./%s_ode2odes.out %s | sed \'s/NAN/Inf/g\' >mtt_out_data.dat\n", ...
system_name, args);
endif
system(command);
## Retrieve out_data
load -force mtt_out_data.dat
y_theta = [y_theta mtt_out_data(:,3:2:1+ny)];
endfor
## System out_data
[n,m]=size(mtt_out_data);
t = mtt_out_data(:,1);
y = mtt_out_data(:,2:2:ny);
x = mtt_out_data(:,3+ny:m);
x_last = mtt_out_data(n,3+ny:m);
endfunction