ADDED mttroot/mtt/bin/trans/m/mtt_identify.m Index: mttroot/mtt/bin/trans/m/mtt_identify.m ================================================================== --- /dev/null +++ mttroot/mtt/bin/trans/m/mtt_identify.m @@ -0,0 +1,38 @@ +function [theta,Theta,Error,Y] = mtt_identify (system_name,y_s,theta_0,criterion,max_iterations) + + ## usage: [theta,Theta,Error] = mtt_identify (system_name,theta_0,criterion,max_iterations) + ## y_s actual system output + ## theta_0 initial parameter estimate + ## criterion convergence criterion + ## max_iterations limit to number of iterations + + if nargin<4 + criterion = 1e-5; + endif + + if nargin <5 + max_iterations = 20; + endif + + alpha = 1.0; + e_last = 1e20; + error=1e10; + theta = theta_0; + Theta = []; + Error = []; + Y = []; + iterations = 0; + while abs(e_last-error)>criterion + theta + iterations = iterations + 1; + e_last = error; + [t,y,y_theta] = mtt_ssimulate(system_name,theta); # Simulate system + Theta = [Theta theta]; # Save parameters + Y = [Y y]; # Save output + E = (y - y_s); # Error(t) + error = (E'*E); # Sum the error + Error = [Error error]; + theta = theta - alpha*((y_theta'*y_theta)\(y_theta'*E)); + endwhile + +endfunction ADDED mttroot/mtt/bin/trans/m/mtt_ssimulate.m Index: mttroot/mtt/bin/trans/m/mtt_ssimulate.m ================================================================== --- /dev/null +++ mttroot/mtt/bin/trans/m/mtt_ssimulate.m @@ -0,0 +1,34 @@ +function [t,y,y_theta] = mtt_ssimulate(system_name,theta); + ## usage: [t,y,y_theta] = mtt_ssimulate(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 + y_theta = []; + for i=1:length(theta) + args=""; + for j=1:length(theta) + i_sensitivity=(j==i); + args = sprintf("%s%i %g ",args, i_sensitivity, theta(j)); + endfor + + command = sprintf("./%s_ode2odes.out %s > mtt_data.dat\n", system_name, args); + system(command); + + ## Retrieve data + load -force mtt_data.dat + t = mtt_data(:,1); + y = mtt_data(:,2); + y_theta = [y_theta mtt_data(:,3)]; + endfor + +endfunction + + +