Consider the following S-system model: XaX. where a and b are constant parameter
ID: 3110907 • Letter: C
Question
Consider the following S-system model: XaX. where a and b are constant parameters and Xa is the independent variable. 2.1. To find the steady state solution (Xi rearrange the terms in the steady state equations to have one term on each side of the equality. Take logarithms on both sides in each equation. Introduce new variables y, In(xis) and y2 In(Xis) and solve the equations to find y, and y in terms of a, b, and X3. Finally, express Xts and X2 in terms of a, b, and X3 2.2. Linearize the model near the steady state, i.e. assume that X1(t XI AX1(t) and X2(t) Xi AX2(t) and find the equations for AX1(t) and AX20t) in the form disregarding all nonlinear terms of the Taylor series. J must be expressed in terms of a, b, X3, Xis, and X 2.3. Now, suppose a 2, b 0.5, and X, 1. Find X X and the eigenvalues of J 2.4 For these values of the parameters implement the model in PLAS on a time interval 0 to 20 with the step 0.1. For initial conditions use Xi 0.5 rand[1] X2 0.5 rand[1]Explanation / Answer
function result=main_function(X,S,beta,Bmet,lbH,ubH,lbB,ubB,iter)
eps=2;% Distance of the boundary
[Vt,vt]=feasible_beta_space(beta,X,S,eps,Bmet);
thres=0.099;
if abs(vt(2:end)) < 0.1111
%**************************************************************************
ct=1;
Bmet2=Bmet;
if vt(1) > ubB
ubB=vt(1)+5;
end
[result2(ct)]=opt_Beta_Atr_diff(X,S,Bmet2,vt,lbH,ubH,lbB,ubB,iter);
vet_error(ct)=result2(ct).error;
aux1=find(abs(result2.h) < thres);
if length(aux1) == 0
aux1=Bmet2;
end
while length(aux1) < length(Bmet2) & length(aux1) ~= 0
ct=ct+1;
Bmet2(aux1)=[];
[Vt1,vt]=feasible_beta_space(beta,X,S,eps,Bmet2);
if vt(1) > ubB
ubB=vt(1)+5;
end
[result2(ct)]=opt_Beta_Atr_diff(X,S,Bmet2,vt,lbH,ubH,lbB,ubB,iter);
aux1=find(abs(result2(ct).h(Bmet2)) < thres);
vet_error(ct)=result2(ct).error;
end
[aa,bb]=min(vet_error);
resultA=result2(bb);
%*************************************************************
if Vt(2:end) >= lbH & Vt(2:end) <= ubH
ct=1;clear vet_error;
if Vt(1) > ubB
ubB=Vt(1)+5;
end
[result1(ct)]=opt_Beta_Atr_diff(X,S,Bmet,Vt,lbH,ubH,lbB,ubB,iter);
vet_error(ct)=result1(ct).error;
aux1=find(abs(result1.h) < thres);
if length(aux1) == 0
aux1=Bmet;
end
while length(aux1) < length(Bmet) & length(aux1) ~= 0
ct=ct+1;
Bmet(aux1)=[];
[Vt,vt]=feasible_beta_space(beta,X,S,eps,Bmet);
if Vt(1) > ubB
ubB=Vt(1)+5;
end
[result1(ct)]=opt_Beta_Atr_diff(X,S,Bmet,Vt,lbH,ubH,lbB,ubB,iter);
aux1=find(abs(result1(ct).h(Bmet)) < thres);
vet_error(ct)=result1(ct).error;
end
[aa,bb]=min(vet_error);
resultB=result1(bb);
else
resultB.error=inf;
end
if resultA.error < resultB.error
result=resultA;
else
result=resultB;
end
else
if Vt(2:end) >= lbH & Vt(2:end) <= ubH
ct=1;clear vet_error;
if Vt(1) > ubB
ubB=Vt(1)+5;
end
[result1(ct)]=opt_Beta_Atr_diff(X,S,Bmet,Vt,lbH,ubH,lbB,ubB,iter);
vet_error(ct)=result1(ct).error;
aux1=find(abs(result1.h) < thres);
if length(aux1) == 0
aux1=Bmet;
end
while length(aux1) < length(Bmet) & length(aux1) ~= 0
ct=ct+1;
Bmet(aux1)=[];
[Vt,vt]=feasible_beta_space(beta,X,S,eps,Bmet);
if Vt(1) > ubB
ubB=Vt(1)+5;
end
[result1(ct)]=opt_Beta_Atr_diff(X,S,Bmet,Vt,lbH,ubH,lbB,ubB,iter);
aux1=find(abs(result1(ct).h(Bmet)) < thres);
vet_error(ct)=result1(ct).error;
end
[aa,bb]=min(vet_error);
result=result1(bb);
else
while sum(abs(vt(2:end)) < 0.1111) < (numel(vt)-1)
beta=vt(1)+1;
[Vt,vt]=feasible_beta_space(beta,X,S,eps,Bmet);
end
ct=1;
Bmet2=Bmet;
if vt(1) > ubB
ubB=vt(1)+5;
end
[result2(ct)]=opt_Beta_Atr_diff(X,S,Bmet2,vt,lbH,ubH,lbB,ubB,iter);
vet_error(ct)=result2(ct).error;
aux1=find(abs(result2.h) < thres);
if length(aux1) == 0
aux1=Bmet2;
end
while length(aux1) < length(Bmet2) & length(aux1) ~= 0
ct=ct+1;
Bmet2(aux1)=[];
[Vt1,vt]=feasible_beta_space(beta,X,S,eps,Bmet2);
if vt(1) > ubB
ubB=vt(1)+5;
end
[result2(ct)]=opt_Beta_Atr_diff(X,S,Bmet2,vt,lbH,ubH,lbB,ubB,iter);
aux1=find(abs(result2(ct).h(Bmet2)) < thres);
vet_error(ct)=result2(ct).error;
end
[aa,bb]=min(vet_error);
result=result2(bb);
end
end
function S=savePLAS(S,S1)
[nl1,nc]=size(S1);
[nl,nc]=size(S);
[snf,spf,typefile]=uiputfile({'*.xml','SBML';'*.txt','Text Files format PLAS'});
if isequal(snf,0) || isequal(spf,0)
disp('User pressed cancel')
elseif typefile == 2
disp(['User selected (Save)', fullfile(spf, snf)])
fid=fopen(strcat(spf,snf,'.txt'),'wt');
for i=1:nl
fprintf(fid,'%s ',S(i,:));
end
fclose(fid)
else
disp(['User selected (Save)', fullfile(spf, snf)])
fid=fopen('temp.txt','wt');
for i=1:nl1
fprintf(fid,'%s ',S1(i,:));
end
fclose(fid)
sbm=SBmodel('temp.txt');
SBexportSBMLmod(snf,sbm)
end
function S= SsystemPLAS(int,alfa,g,beta,h)
S='';
sss='';
nm=length(alfa);
for i=1:nm
ctp=0;
ctd=0;
stpa='';
stda='';
for j=1:nm
if g(i,j) ~= 0
stp1=strcat(' X',num2str(j),'^',num2str(g(i,j)));
stpa=strcat(stpa,stp1);
ctp=ctp+1;
end
if h(i,j) ~= 0
std1=strcat(' X',num2str(j),'^',num2str(h(i,j)));
stda=strcat(stda,std1);
ctd=ctd+1;
end
end
S1=strcat('X',num2str(i),'''','=', num2str(alfa(i)),stpa,' - ', num2str(beta(i)) , stda);
S=strvcat(S,S1);
ss=strcat('X',num2str(i),'=',num2str(int(i)));
sss=strvcat(sss,ss);
end
aux=strvcat('t0=0','hr=0.1','tf=5');
S=strvcat('S-system',S,' ','Initial Values',sss,' ','INITIAL TIME, FINAL TIME and REPORT INTERVAL:',aux);
function SBcreateODEfile(varargin)
% SBcreateODEfile: creates an m-file that can be simulated by
% using standard MATLAB integrators (ODE45, ODE23s, etc.)
%
% USAGE:
% ======
% [] = SBcreateODEfile(sbm)
% [] = SBcreateODEfile(sbm,filename)
% [] = SBcreateODEfile(sbm,filename,dataFileFlag)
% [] = SBcreateODEfile(sbm,filename,dataFileFlag,eventFlag)
%
% sbm: SBmodel to convert to an ODE file description
% filename: filename for the created ODE file (without extension)
% dataFileFlag: =1: creates an additional file allowing to determine
% variable and reaction rate values for given state values
% and time vector.
% =0: doe not create the additional file
% eventFlag: =1: creates 2 additional files for the handling of events in the model.
% =0: does not create these two files.
% THE EVENT FILES ARE ONLY CREATED IN CASE THAT EVENTS ARE
% PRESENT IN THE MODEL
%
% DEFAULT VALUES:
% ===============
% filename: constructed from the models name
% dataFileFlag: 0
% eventFlag: 0
% Information:
% ============
% Systems Biology Toolbox for MATLAB
% Copyright (C) 2005 Henning Schmidt, FCC, henning@fcc.chalmers.se
%
% This program is free software; you can redistribute it and/or
% modify it under the terms of the GNU General Public License
% as published by the Free Software Foundation; either version 2
% of the License, or (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program; if not, write to the Free Software
% Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PROCESS VARIABLE INPUT ARGUMENTS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
dataFileFlag = 0;
eventFlag = 0;
if nargin == 1,
sbm = varargin{1};
% if no filename provided then use the name of the SBmodel as filename
% remove white spaces from the name.
functionName = strrep(sbm.name,' ','');
functionName = strrep(functionName,'-','');
functionName = strrep(functionName,'+','');
functionName = strrep(functionName,'*','');
functionName = strrep(functionName,'/','');
filename = strcat(functionName,'.m');
elseif nargin == 2,
sbm = varargin{1};
% extract filename from input argument number 2 - needed since
% SBsimulate gives a whole path as filename
[PATHSTR,functionName,EXT,VERSN] = fileparts(varargin{2});
if length(PATHSTR) == 0,
filename = strcat(functionName,'.m');
else
filename = strcat(PATHSTR,'/',functionName,'.m');
end
elseif nargin == 3,
sbm = varargin{1};
% extract filename from input argument number 2 - needed since
% SBsimulate gives a whole path as filename
[PATHSTR,functionName,EXT,VERSN] = fileparts(varargin{2});
if length(PATHSTR) == 0,
filename = strcat(functionName,'.m');
else
filename = strcat(PATHSTR,'/',functionName,'.m');
end
dataFileFlag = varargin{3};
elseif nargin == 4,
sbm = varargin{1};
% extract filename from input argument number 2 - needed since
% SBsimulate gives a whole path as filename
[PATHSTR,functionName,EXT,VERSN] = fileparts(varargin{2});
if length(PATHSTR) == 0,
filename = strcat(functionName,'.m');
else
filename = strcat(PATHSTR,'/',functionName,'.m');
end
dataFileFlag = varargin{3};
eventFlag = varargin{4};
else
error('Incorrect number of input arguments.');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% CHECK IF MODEL CONTAINS STATES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if length(sbm.states) == 0,
error('The model does not contain any states');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% HANDLE THE EVENT FLAG
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% if the eventflag is set to a non-zero value AND at least on event is
% present in the model an event function file is created that has to be
% used when calling the integrator. additionally a function is created
% that determines the new values to which the states have to be set.
if eventFlag == 1 && length(sbm.events) > 0,
filenameEvent = strrep(filename,functionName,strcat('event_',functionName));
filenameEventAssignment = strrep(filename,functionName,strcat('event_assignment_',functionName));
SBcreateEventFunction(sbm,filenameEvent);
SBcreateEventAssignmentFunction(sbm,filenameEventAssignment);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% HANDLE THE DATAFILE FLAG
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if dataFileFlag == 1,
filenameDataFile = strrep(filename,functionName,strcat('datafile_',functionName));
SBcreateSimulationDataFunction(sbm,filenameDataFile);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% WRITE WARNING MESSAGE IN CASE EVENT FLAG IS NOT SET AND EVENTS ARE PRESENT
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if eventFlag == 0 && length(sbm.events) > 0,
if length(sbm.events) == 1,
disp(sprintf('Warning: 1 event is present in the model. For this analysis it is not taken into account.'));
else
disp(sprintf('Warning: %d events are present in the model. For this analysis they are not taken into account.',length(sbm.events)));
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% OPEN FILE FOR WRITING
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fid = fopen(filename,'w');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% WRITE FUNCTION DEFINITION
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fprintf(fid,'function [output] = %s(varargin) ',functionName);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% WRITE THE HEADER
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% %s ',sbm.name);
fprintf(fid,'%% Generated: %s ',datestr(now));
fprintf(fid,'%% ');
fprintf(fid,'%% [output] = %s() => output = initial conditions in column vector ',functionName);
fprintf(fid,'%% [output] = %s(''states'') => output = species names in cell-array ',functionName);
fprintf(fid,'%% [output] = %s(''parameters'') => output = parameter names in cell-array ',functionName);
fprintf(fid,'%% [output] = %s(''parametervalues'') => output = parameter values in column vector ',functionName);
fprintf(fid,'%% [output] = %s(time,statevector) => output = time derivatives in column vector ',functionName);
fprintf(fid,'%% ');
fprintf(fid,'%% State names and ordering: ');
fprintf(fid,'%% ');
for k = 1:length(sbm.states),
fprintf(fid,'%% statevector(%d): %s ',k,sbm.states(k).name);
end
fprintf(fid,'%% ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,' ');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% HANDLE VARARGINS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% HANDLE VARIABLE INPUT ARGUMENTS ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'if nargin == 0, ');
% Return initial conditions for state variables
fprintf(fid,' %% Return initial conditions of the state variables ');
if length(sbm.states) > 0,
fprintf(fid,' output = [');
else
fprintf(fid,' output = []; ');
end
count = 0;
for k = 1:length(sbm.states)-1,
fprintf(fid,'%s, ',num2str(sbm.states(k).initialCondition,20));
count = count + 1;
if count == 10, fprintf(fid,'... '); count = 0; end
end
if length(sbm.states) > 0,
fprintf(fid,'%f]; ',sbm.states(length(sbm.states)).initialCondition);
end
fprintf(fid,' output = output(:); ');
fprintf(fid,' return ');
fprintf(fid,'elseif nargin == 1, ');
fprintf(fid,' if strcmp(varargin{1},''states''), ');
fprintf(fid,' %% Return state names in cell-array ');
if length(sbm.states) > 0,
fprintf(fid,' output = {');
else
fprintf(fid,' output = {}; ');
end
count = 0;
for k = 1:length(sbm.states)-1,
fprintf(fid,'''%s'', ',sbm.states(k).name);
count = count + 1;
if count == 10, fprintf(fid,'... '); count = 0; end
end
if length(sbm.states) > 0,
fprintf(fid,'''%s''}; ',sbm.states(length(sbm.states)).name);
end
fprintf(fid,' elseif strcmp(varargin{1},''parameters''), ');
fprintf(fid,' %% Return parameter names in cell-array ');
if length(sbm.parameters) > 0,
fprintf(fid,' output = {');
else
fprintf(fid,' output = {}; ');
end
count = 0;
for k = 1:length(sbm.parameters)-1,
fprintf(fid,'''%s'', ',sbm.parameters(k).name);
count = count + 1;
if count == 10, fprintf(fid,'... '); count = 0; end
end
if length(sbm.parameters) > 0,
fprintf(fid,'''%s''}; ',sbm.parameters(length(sbm.parameters)).name);
end
fprintf(fid,' elseif strcmp(varargin{1},''parametervalues''), ');
fprintf(fid,' %% Return parameter values in column vector ');
if length(sbm.parameters) > 0,
fprintf(fid,' output = [');
else
fprintf(fid,' output = []; ');
end
count = 0;
for k = 1:length(sbm.parameters)-1,
fprintf(fid,'%s, ',num2str(sbm.parameters(k).value, 20));
count = count + 1;
if count == 10, fprintf(fid,'... '); count = 0; end
end
if length(sbm.parameters) > 0,
fprintf(fid,'%s]; ',num2str(sbm.parameters(length(sbm.parameters)).value,20));
end
fprintf(fid,' else ');
fprintf(fid,' error(''Wrong input arguments! Please read the help text to the ODE file.''); ');
fprintf(fid,' end ');
fprintf(fid,' output = output(:); ');
fprintf(fid,' return ');
fprintf(fid,'elseif nargin == 2, ');
fprintf(fid,' time = varargin{1}; ');
fprintf(fid,' statevector = varargin{2}; ');
fprintf(fid,'elseif nargin == 3, ');
fprintf(fid,'%% For the reason of initialization and finishing, NARGIN is sometimes 3. ');
fprintf(fid,'%% The third input argument id of no interest for this toolbox! ');
fprintf(fid,' time = varargin{1}; ');
fprintf(fid,' statevector = varargin{2}; ');
fprintf(fid,'else ');
fprintf(fid,' error(''Wrong input arguments! Please read the help text to the ODE file.''); ');
fprintf(fid,'end ');
fprintf(fid,' ');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% INITIALIZE THE STATES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% STATES ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
for k = 1:length(sbm.states)
fprintf(fid,'%s = statevector(%d); ',sbm.states(k).name,k);
end
fprintf(fid,' ');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% INITIALIZE THE PARAMETES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if length(sbm.parameters),
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% PARAMETERS ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
for k = 1:length(sbm.parameters),
fprintf(fid,'%s = %s; ',sbm.parameters(k).name,num2str(sbm.parameters(k).value,20));
end
fprintf(fid,' ');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% INITIALIZE THE VARIABLES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if length(sbm.variables),
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% VARIABLES ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
for k = 1:length(sbm.variables),
fprintf(fid,'%s = %s; ',sbm.variables(k).name,sbm.variables(k).formula);
end
fprintf(fid,' ');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% INITIALIZE THE REACTIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if length(sbm.reactions),
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% REACTION KINETICS ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
for k = 1:length(sbm.reactions),
fprintf(fid,'%s = %s; ',sbm.reactions(k).name,sbm.reactions(k).formula);
end
fprintf(fid,' ');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% WRITE THE ODES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% DIFFERENTIAL EQUATIONS ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
for k = 1:length(sbm.states),
fprintf(fid,'%s_dot = %s; ',sbm.states(k).name,sbm.states(k).ODE);
end
fprintf(fid,' ');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% WRITE RETURN VALUES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% RETURN VALUES ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
for k = 1:length(sbm.states),
fprintf(fid,'output(%d) = %s_dot; ',k,sbm.states(k).name);
end
fprintf(fid,'%% return a column vector ');
% Return a column vector
fprintf(fid,'output = output(:); ');
fprintf(fid,'return ');
fprintf(fid,' ');
fprintf(fid,' ');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% WRITE THE FUNCTIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if length(sbm.functions),
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% FUNCTIONS ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
for k = 1:length(sbm.functions),
fprintf(fid,'function [result] = %s(%s) ',sbm.functions(k).name,sbm.functions(k).arguments);
fprintf(fid,'result = %s; ',sbm.functions(k).formula);
fprintf(fid,'return ');
fprintf(fid,' ');
end
fprintf(fid,' ');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% WRITE THE MATLAB FUNCTIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if ~isempty(sbm.functionsMATLAB),
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%% MATLAB FUNCTIONS ');
fprintf(fid,'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ');
fprintf(fid,'%s',sbm.functionsMATLAB);
fprintf(fid,' ');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% CLOSE THE FILE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
status = fclose(fid);
% Return
return
Related Questions
drjack9650@gmail.com
Navigate
Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.