% FILENAME: SIMULT.M
% DATE: NOV 1993 
% VERSION: 2.0
% REVISED : January 1995
% COPYRIGHT  1994-5 by T. Marlin and R. Cipolla
% This M-file enables the user to define time-related simulation 
% parameters, number of cases, call the simulation, and display results
% *********************************************************************

Y0 = Y0i;


while 1
  clc
  disp('*************************************************************')
  disp('*                 DYNAMIC SIMULATION                        *')
  disp('*                                                           *')
  disp('*************************************************************')
  disp('')
  disp('SELECT THE APPROPRIATE MENU ITEM')
  disp('')
  disp('MODIFY...')
  disp('                                                   PRESENT VALUES')
  fprintf('1) Time to run simulation                          %4.2f\n',time)
  fprintf('2) Time step                                       %4.3f\n',deltat)
  fprintf('3) Specify New Set Point                           %4.2f\n',SP)
  fprintf('4) Specify Disturbance Size                        %4.2f\n',D)
  disp('5) Execute dynamic simulation for Disturbance')
  disp('6) Execute dynamic simulation for Set Point change')
  disp('7) Run Model Mismatch Cases for a Disturbance')
  disp('8) Run Model Mismatch Cases for a Set Point change')
  disp('9) Plot Results' )  
  disp('10) Save graph(s) to META-file')
  disp('11) Return to Main Menu')      
  disp('')
  Sim = input('Enter the desired selection:   ');
  
  if isempty(Sim), break, end
      
  Sim = round(Sim);
  if (Sim<1) | (Sim>11)
      disp('')
      disp('Not a valid selection')
      disp('')
      disp('Press ENTER to continue')
      pause;
  end    
   
  if Sim == 1
      junk = input('Enter value for end time of integration:  ') ;
      
      if isempty(junk)
          junk=time;
      end

      if junk < 0
       disp('')
       disp(' You cannot have a negative total time')
       disp('')
       disp(' Press any key to continue')
       pause
       junk=time;
      end

      time=junk;
      if (time/deltat) > lim_it
          disp('')
          disp('Exceeded matrix limitations, total time has been recalculated') 
          disp('')
          disp('Press Enter to continue') 
          time=(lim_it*deltat);
          pause;
      end
      mvtest=2;  cn_t=1;  scn_t=1; maxit=ceil(time/deltat);
  elseif Sim == 2
      junk = input('Enter value for time step:  ') ;
      
      if isempty(junk)
          junk=deltat;
      end
      
      if junk <= 0
       disp('')
       disp(' You cannot have a negative or zero time step')
       disp('')
       disp(' Press any key to continue')
       pause
       junk=deltat;
      end

      deltat=junk;
      if (time/deltat) > lim_it
          disp('')
          disp('Exceeded matrix limitations, time step has been recalculated')
          disp('')
          disp('Press Enter to continue')
          deltat=(time/lim_it);
          pause;
      end
      mvtest=2;  cn_t=1; scn_t=1; maxit=ceil(time/deltat);

  elseif Sim == 3
      junk = input('Enter a value for the SET POINT :  ');
      
      if isempty(junk)
          junk=SP;
      end
      
      SP=junk;
      disp ('Enter the time at which to implement the set point'); 
      junk= input('Default (blank) is 5% of total simulation time:  ');
      
       if isempty(junk)
           junk=0.05*time;
       end
      
      if junk < 0
       disp('')
       disp(' Must be positive')
       disp('')
       disp(' Press any key to continue')
       pause
       junk=tset;
      end
      if junk < theta(1)
       clc
       disp('')
       disp(' Due to program limitations, T set-point must be > dead time')
       disp('')
       disp(' Press any key to continue ...')
       pause
       junk=tset;
      end

       tset=junk;
       mvtest=2;  cn_t=1; scn_t=1;

   elseif Sim == 4      
        
        
        junk = input('Enter the size (+ or -) of the disturbance:  ');
        
        if isempty(junk)
           junk=D;
        end
        D=junk;
       
        disp('Enter the time of the disturbance:   '); 
        junk = input('Default (blank) is 5% of total simulation time:  ');
        
        if isempty(junk)
          junk=0.05*time;
        end
        
        if junk < 0
          disp('')
          disp(' Must be positive')
          disp('')
          disp(' Press any key to continue')
          pause
          junk=tdist;
        end
        if junk < theta(1)
          clc
          disp('')
          disp(' Due to program limitations, Time of disturbance must be > dead time')
          disp('')
          disp(' Press any key to continue ...')
          pause
          junk=tset;
        end

       tdist=junk;
       mvtest=2;  cn_t=1; scn_t=1;

   elseif Sim == 5 | Sim == 7   
 
 %  ensure that the current values of kc ti, td are used
       junk2 = 2;
       Y0 = [ kc ti];
 
       pp=0; op_ts=0; o_func=0;
       change = 1;
       betaa=tdist;
       usegrh = 1;
       if mvtest == 0 & cn_t > 1
         useopt=1;
         mm=1;
         if Sim == 7, useopt=0; usegrh=0; mm=0; end
         integ4A(Y0);
       else    
          p=1;
          mvbound1
          useopt=1;
          mm=1;
          if Sim == 7, useopt=0; usegrh=0; mm=0; end
          integ4A(Y0);
       end
   
   elseif Sim== 6 | Sim==8 
 %  ensure that the current values of kc ti, td are used
       junk2 = 2;
       Y0 = [ kc ti];

       pp=0; op_ts=0; o_func=0;
       change = 2; 
       betaa=tset;
       usegrh = 1;
       if mvtest == 0 & scn_t > 1
          useopt=2;
          mm=1;
          if Sim == 8, useopt=0; usegrh=0; mm=0; end
          integ4A(Y0);     
       else
          p=1;
          mvbound1
          useopt=2;
          mm=1;
          if Sim == 8, useopt=0; usegrh=0; mm=0; end
          integ4A(Y0);
       end
   elseif Sim == 9     
       if usegrh == 2 
         disp('')
         disp('Nothing to plot, run a simulation before using this option') 
         disp('')
         disp('Press Enter to continue')
         pause;
       else 
         save_graph = 0;
         grph2
       end
   
   elseif Sim == 10   % Saves current graphs to META-file
%      save_graph=1;  % This facility not in MATLAB 4.2c
%      grph2

       disp ('  ')
       disp (' This option is not available in this version of MACTUNE')
       disp (' press enter to continue')
       pause

   else
       break         % quit menu
  end
end   



