Upload new file

CRA_Updated_Version/CRA_Toolbox_bash/start_simulation.m

+%Function for simulating the ODE model using the perturbed parameter
+%vector. It takes in input:
+%1. model which is the struct containing all the characheristics of the
+%model
+%2. Nr is the number of independent realizations to perform;
+%3. LBpi and UBpi are the lower and upper boundaries of the Latin
+%Hypercube;
+%4. Ns is the number of samples to generate for creating the Latin
+%Hypercube;
+%It returns in output the array AllResults containing all the output
+%variables simulated using the parameter vectors in AllPerturbations
+
+function [AllResults, AllPerturbations]=start_simulation(model,Nr,LBpi,UBpi,Ns)
+
+nominal_parameters=model.nominal_parameters;
+num_parameters=length(nominal_parameters);
+
+AllResults=cell(Nr,model.num_observables);
+
+AllPerturbations=cell(Nr,1);
+
+h = waitbar(0,'Please wait...');
+
+for k=1:Nr
+    
+  waitbar(k/Nr,h,['Realization number ',num2str(k)]);
+  
+  hypercube=LatinHypercube(LBpi,UBpi,Ns);
+  perturbation=hypercube.generate_sample(num_parameters);
+  perturbed_p=repmat(nominal_parameters,length(perturbation),1).*perturbation;
+  AllPerturbations{k,1}=perturbation;
+  disp(strcat('Generated a Latin Hypercube of size',{' '},sprintf('%.0f',size(perturbation,1)),'x',sprintf('%.0f',size(perturbation,2))));
+    %psim_all=cell(Ns,model.num_observables);
+  psim_all=cell(Ns,1);  
+  parfor i=1:Ns
+     pi=perturbed_p(i,:);
+     modelsim=@(t,x)feval(model.name,t,x,pi,model.input,model.total_proteins);
+     [tsimi,psimi]=feval(model.odesolver,modelsim,model.time,model.initial_conditions);
+    % for t=1:model.num_observables
+     %  psim_all{i,t}=psimi(:,t);
+     %end
+     msgstr=lastwarn;
+     if isequal(msgstr,'')        
+        psim_all{i,1}=psimi;
+     else 
+        psim_all{i,1}=NaN(length(model.time),model.num_observables);
+     end
+      lastwarn('')
+   end 
+  
+ for j=1:model.num_observables
+    for w=1:Ns
+     AllResults{k,j}{w,1}=psim_all{w,1}(:,j);
+    end
+  end
+  %AllResults{k,1:end}=psim_all;
+  disp(strcat('Completed one realization. The array of results has size',{' '},sprintf('%.0f',size(psim_all,1)),'x',sprintf('%.0f',size(psim_all,2))));
+ end
+
+delete(gcp) %close parallel pool
+close (h)  %close the waitbar
+end
\ No newline at end of file