CVodeMonitor is the default CVODES monitoring function.
   To use it, set the Monitor property in CVodeSetOptions to
   'CVodeMonitor' or to @CVodeMonitor and 'MonitorData' to mondata
   (defined as a structure).
  
   With default settings, this function plots the evolution of the step 
   size, method order, and various counters.
   
   Various properties can be changed from their default values by passing
   to CVodeSetOptions, through the property 'MonitorData', a structure
   MONDATA with any of the following fields. If a field is not defined, 
   the corresponding default value is used.

   Fields in MONDATA structure:
     o stats [ {true} | false ]
         If true, report the evolution of the step size and method order.
     o cntr [ {true} | false ]
         If true, report the evolution of the following counters:
         nst, nfe, nni, netf, ncfn (see CVodeGetStats)
     o mode [ {'graphical'} | 'text' | 'both' ] 
         In graphical mode, plot the evolutions of the above quantities.
         In text mode, print a table.
     o sol  [ true | {false} ]
         If true, plot solution components.
     o sensi [ true | {false} ]
         If true and if FSA is enabled, plot sensitivity components.
     o select [ array of integers ]
         To plot only particular solution components, specify their indeces in
         the field select. If not defined, but sol=true, all components are plotted.
     o updt [ integer | {50} ]
         Update frequency. Data is posted in blocks of dimension n.
     o skip [ integer | {0} ]
         Number of integrations steps to skip in collecting data to post.
     o dir [ {1} | -1 ]
         Specifies forward or backward integration.
     o post [ {true} | false ]
         If false, disable all posting. This option is necessary to disable
         monitoring on some processors when running in parallel.

   See also CVodeSetOptions, CVMonitorFn

   NOTES:
     1. The argument mondata is REQUIRED. Even if only the default options
        are desired, set mondata=struct; and pass it to CVodeSetOptions.
     2. The yQ argument is currently ignored.

0001 function [new_data] = CVodeMonitor(call, T, Y, YQ, YS, data)
0002 %CVodeMonitor is the default CVODES monitoring function.
0003 %   To use it, set the Monitor property in CVodeSetOptions to
0004 %   'CVodeMonitor' or to @CVodeMonitor and 'MonitorData' to mondata
0005 %   (defined as a structure).
0006 %
0007 %   With default settings, this function plots the evolution of the step
0008 %   size, method order, and various counters.
0009 %
0010 %   Various properties can be changed from their default values by passing
0011 %   to CVodeSetOptions, through the property 'MonitorData', a structure
0012 %   MONDATA with any of the following fields. If a field is not defined,
0013 %   the corresponding default value is used.
0014 %
0015 %   Fields in MONDATA structure:
0016 %     o stats [ {true} | false ]
0017 %         If true, report the evolution of the step size and method order.
0018 %     o cntr [ {true} | false ]
0019 %         If true, report the evolution of the following counters:
0020 %         nst, nfe, nni, netf, ncfn (see CVodeGetStats)
0021 %     o mode [ {'graphical'} | 'text' | 'both' ]
0022 %         In graphical mode, plot the evolutions of the above quantities.
0023 %         In text mode, print a table.
0024 %     o sol  [ true | {false} ]
0025 %         If true, plot solution components.
0026 %     o sensi [ true | {false} ]
0027 %         If true and if FSA is enabled, plot sensitivity components.
0028 %     o select [ array of integers ]
0029 %         To plot only particular solution components, specify their indeces in
0030 %         the field select. If not defined, but sol=true, all components are plotted.
0031 %     o updt [ integer | {50} ]
0032 %         Update frequency. Data is posted in blocks of dimension n.
0033 %     o skip [ integer | {0} ]
0034 %         Number of integrations steps to skip in collecting data to post.
0035 %     o dir [ {1} | -1 ]
0036 %         Specifies forward or backward integration.
0037 %     o post [ {true} | false ]
0038 %         If false, disable all posting. This option is necessary to disable
0039 %         monitoring on some processors when running in parallel.
0040 %
0041 %   See also CVodeSetOptions, CVMonitorFn
0042 %
0043 %   NOTES:
0044 %     1. The argument mondata is REQUIRED. Even if only the default options
0045 %        are desired, set mondata=struct; and pass it to CVodeSetOptions.
0046 %     2. The yQ argument is currently ignored.
0047 
0048 % Radu Serban <radu@llnl.gov>
0049 % Copyright (c) 2005, The Regents of the University of California.
0050 % $Revision: 1.4 $Date: 2006/10/05 22:12:20 $
0051 
0052 
0053 new_data = [];
0054 
0055 if call == 0
0056 
0057 % Initialize unspecified fields to default values.
0058   data = initialize_data(data);  
0059 
0060 % Open figure windows
0061   if data.post
0062 
0063     if data.grph
0064       if data.stats | data.cntr
0065         data.hfg = figure;
0066       end
0067 %     Number of subplots in figure hfg
0068       if data.stats
0069         data.npg = data.npg + 2;
0070       end
0071       if data.cntr
0072         data.npg = data.npg + 1;
0073       end
0074     end
0075     
0076     if data.text
0077       if data.cntr | data.stats
0078         data.hft = figure;
0079       end
0080     end
0081 
0082     if data.sol | data.sensi
0083       data.hfs = figure; 
0084     end
0085   
0086   end
0087   
0088 % Initialize other private data
0089   data.i = 0;
0090   data.n = 1;  
0091   data.t = zeros(1,data.updt);
0092   if data.stats
0093     data.h = zeros(1,data.updt);
0094     data.q = zeros(1,data.updt);
0095   end
0096   if data.cntr
0097     data.nst = zeros(1,data.updt);
0098     data.nfe = zeros(1,data.updt);
0099     data.nni = zeros(1,data.updt);
0100     data.netf = zeros(1,data.updt);
0101     data.ncfn = zeros(1,data.updt);
0102   end
0103 
0104   data.first = true;        % the next one will be the first call = 1
0105   data.initialized = false; % the graphical windows were not initalized
0106   
0107   new_data = data;
0108   
0109   return;
0110 
0111 else
0112 
0113 % If this is the first call ~= 0,
0114 % use Y and YS for additional initializations
0115   
0116   if data.first
0117 
0118     if isempty(YS)
0119       data.sensi = false;
0120     end
0121     
0122     if data.sol | data.sensi
0123       
0124       if isempty(data.select)
0125       
0126         data.N = length(Y);
0127         data.select = [1:data.N];
0128         
0129       else
0130         
0131         data.N = length(data.select);
0132         
0133       end
0134 
0135       if data.sol
0136         data.y = zeros(data.N,data.updt);
0137         data.nps = data.nps + 1;
0138       end
0139         
0140       if data.sensi
0141         data.Ns = size(YS,2);
0142         data.ys = zeros(data.N, data.Ns, data.updt);
0143         data.nps = data.nps + data.Ns;
0144       end
0145       
0146     end
0147     
0148     data.first = false;
0149   
0150   end
0151   
0152 % Extract variables from data
0153 
0154   hfg  = data.hfg;
0155   hft  = data.hft;
0156   hfs  = data.hfs;
0157   npg  = data.npg;
0158   nps  = data.nps;
0159   i    = data.i;
0160   n    = data.n;
0161   t    = data.t;
0162   N    = data.N;
0163   Ns   = data.Ns;
0164   y    = data.y;
0165   ys   = data.ys;
0166   h    = data.h;
0167   q    = data.q;
0168   nst  = data.nst;
0169   nfe  = data.nfe;
0170   nni  = data.nni;
0171   netf = data.netf;
0172   ncfn = data.ncfn;
0173   
0174 end
0175 
0176 
0177 % Load current statistics?
0178 
0179 if call == 1
0180 
0181   if i ~= 0
0182     i = i-1;
0183     data.i = i;
0184     new_data = data;
0185     return;
0186   end
0187 
0188   if data.dir == 1
0189     si = CVodeGetStats;
0190   else
0191     si = CVodeGetStatsB;
0192   end
0193 
0194   t(n) = si.tcur;
0195   
0196   if data.stats
0197     h(n) = si.hlast;
0198     q(n) = si.qlast;
0199   end
0200   
0201   if data.cntr
0202     nst(n) = si.nst;
0203     nfe(n) = si.nfe;
0204     nni(n) = si.nni;
0205     netf(n) = si.netf;
0206     ncfn(n) = si.ncfn;
0207   end
0208 
0209   if data.sol
0210     for j = 1:N
0211       y(j,n) = Y(data.select(j));
0212     end
0213   end
0214 
0215   if data.sensi
0216     for k = 1:Ns
0217       for j = 1:N
0218         ys(j,k,n) = YS(data.select(j),k);
0219       end
0220     end
0221   end
0222   
0223 end
0224 
0225 % Is it time to post?
0226 
0227 if data.post & (n == data.updt | call==2)
0228 
0229   if call == 2
0230     n = n-1;
0231   end
0232   
0233   if ~data.initialized
0234 
0235     if (data.stats | data.cntr) & data.grph
0236       graphical_init(n, hfg, npg, data.stats, data.cntr, data.dir, ...
0237                      t, h, q, nst, nfe, nni, netf, ncfn);
0238     end
0239     
0240     if (data.stats | data.cntr) & data.text
0241       text_init(n, hft, data.stats, data.cntr, ...
0242                 t, h, q, nst, nfe, nni, netf, ncfn);
0243     end
0244 
0245     if data.sol | data.sensi
0246       sol_init(n, hfs, nps, data.sol, data.sensi, data.dir, ...
0247                N, Ns, t, y, ys);
0248     end
0249     
0250     data.initialized = true;
0251   
0252   else
0253 
0254     if (data.stats | data.cntr) & data.grph
0255       graphical_update(n, hfg, npg, data.stats, data.cntr, ...
0256                        t, h, q, nst, nfe, nni, netf, ncfn);
0257     end
0258 
0259     if (data.stats | data.cntr) & data.text
0260       text_update(n, hft, data.stats, data.cntr, ...
0261                   t, h, q, nst, nfe, nni, netf, ncfn);
0262     end
0263     
0264     if data.sol
0265       sol_update(n, hfs, nps, data.sol, data.sensi, N, Ns, t, y, ys);
0266     end
0267       
0268   end
0269 
0270   if call == 2
0271     
0272     if (data.stats | data.cntr) & data.grph
0273       graphical_final(hfg, npg, data.cntr, data.stats);
0274     end
0275     
0276     if data.sol | data.sensi
0277       sol_final(hfs, nps, data.sol, data.sensi, N, Ns);
0278     end
0279 
0280     return;
0281   
0282   end
0283   
0284   n = 1;
0285 
0286 else
0287 
0288   n = n + 1;
0289 
0290 end
0291 
0292 
0293 % Save updated values in data
0294 
0295 data.i    = data.skip;
0296 data.n    = n;
0297 data.npg  = npg;
0298 data.t    = t;
0299 data.y    = y;
0300 data.ys   = ys;
0301 data.h    = h;
0302 data.q    = q;
0303 data.nst  = nst;
0304 data.nfe  = nfe;
0305 data.nni  = nni;
0306 data.netf = netf;
0307 data.ncfn = ncfn;
0308   
0309 new_data = data;
0310 
0311 return;
0312 
0313 %-------------------------------------------------------------------------
0314 
0315 function data = initialize_data(data)
0316 
0317 if ~isfield(data,'mode')
0318   data.mode = 'graphical';
0319 end
0320 if ~isfield(data,'updt')
0321   data.updt = 50;
0322 end
0323 if ~isfield(data,'skip')
0324   data.skip = 0;
0325 end
0326 if ~isfield(data,'stats')
0327   data.stats = true;
0328 end
0329 if ~isfield(data,'cntr')
0330   data.cntr = true;
0331 end
0332 if ~isfield(data,'sol')
0333   data.sol = false;
0334 end
0335 if ~isfield(data,'sensi')
0336   data.sensi = false;
0337 end
0338 if ~isfield(data,'select')
0339   data.select = [];
0340 end
0341 if ~isfield(data,'dir')
0342   data.dir = 1;
0343 end
0344 if ~isfield(data,'post')
0345   data.post = true;
0346 end
0347 
0348 data.grph = true;
0349 data.text = true;
0350 if strcmp(data.mode,'graphical')
0351   data.text = false;
0352 end
0353 if strcmp(data.mode,'text')
0354   data.grph = false;
0355 end
0356 
0357 if ~data.sol & ~data.sensi
0358   data.select = [];
0359 end
0360   
0361 % Other initializations
0362 data.npg = 0;
0363 data.nps = 0;
0364 data.hfg = 0;
0365 data.hft = 0;
0366 data.hfs = 0;
0367 data.h = 0;
0368 data.q = 0;
0369 data.nst = 0;
0370 data.nfe = 0;
0371 data.nni = 0;
0372 data.netf = 0;
0373 data.ncfn = 0;
0374 data.N = 0;
0375 data.Ns = 0;
0376 data.y = 0;
0377 data.ys = 0;
0378 
0379 %-------------------------------------------------------------------------
0380 
0381 function [] = graphical_init(n, hfg, npg, stats, cntr, dir, ...
0382                              t, h, q, nst, nfe, nni, netf, ncfn)
0383 
0384 fig_name = 'CVODES run statistics';
0385 
0386 % If this is a parallel job, look for the MPI rank in the global
0387 % workspace and append it to the figure name
0388 
0389 global sundials_MPI_rank
0390 
0391 if ~isempty(sundials_MPI_rank)
0392   fig_name = sprintf('%s (PE %d)',fig_name,sundials_MPI_rank);
0393 end
0394 
0395 figure(hfg);
0396 set(hfg,'Name',fig_name);
0397 set(hfg,'color',[1 1 1]);
0398 pl = 0;
0399 
0400 % Time label and figure title
0401 
0402 if dir==1
0403   tlab = '\rightarrow   t   \rightarrow';
0404 else
0405   tlab = '\leftarrow   t   \leftarrow';
0406 end
0407 
0408 % Step size and order
0409 if stats
0410   pl = pl+1;
0411   subplot(npg,1,pl)
0412   semilogy(t(1:n),abs(h(1:n)),'-');
0413   hold on;
0414   box on;
0415   grid on;
0416   xlabel(tlab);
0417   ylabel('|Step size|');
0418   
0419   pl = pl+1;
0420   subplot(npg,1,pl)
0421   plot(t(1:n),q(1:n),'-');
0422   hold on;
0423   box on;
0424   grid on;
0425   xlabel(tlab);
0426   ylabel('Order');
0427 end
0428   
0429 % Counters
0430 if cntr
0431   pl = pl+1;
0432   subplot(npg,1,pl)
0433   plot(t(1:n),nst(1:n),'k-');
0434   hold on;
0435   plot(t(1:n),nfe(1:n),'b-');
0436   plot(t(1:n),nni(1:n),'r-');
0437   plot(t(1:n),netf(1:n),'g-');
0438   plot(t(1:n),ncfn(1:n),'c-');
0439   box on;
0440   grid on;
0441   xlabel(tlab);
0442   ylabel('Counters');
0443 end
0444 
0445 drawnow;
0446 
0447 %-------------------------------------------------------------------------
0448 
0449 function [] = graphical_update(n, hfg, npg, stats, cntr, ...
0450                                t, h, q, nst, nfe, nni, netf, ncfn)
0451 
0452 figure(hfg);
0453 pl = 0;
0454   
0455 % Step size and order
0456 if stats
0457   pl = pl+1;
0458   subplot(npg,1,pl)
0459   hc = get(gca,'Children');
0460   xd = [get(hc,'XData') t(1:n)];
0461   yd = [get(hc,'YData') abs(h(1:n))];
0462   set(hc, 'XData', xd, 'YData', yd);
0463   
0464   pl = pl+1;
0465   subplot(npg,1,pl)
0466   hc = get(gca,'Children');
0467   xd = [get(hc,'XData') t(1:n)];
0468   yd = [get(hc,'YData') q(1:n)];
0469   set(hc, 'XData', xd, 'YData', yd);
0470 end
0471 
0472 % Counters
0473 if cntr
0474   pl = pl+1;
0475   subplot(npg,1,pl)
0476   hc = get(gca,'Children');
0477 % Attention: Children are loaded in reverse order!
0478   xd = [get(hc(1),'XData') t(1:n)];
0479   yd = [get(hc(1),'YData') ncfn(1:n)];
0480   set(hc(1), 'XData', xd, 'YData', yd);
0481   yd = [get(hc(2),'YData') netf(1:n)];
0482   set(hc(2), 'XData', xd, 'YData', yd);
0483   yd = [get(hc(3),'YData') nni(1:n)];
0484   set(hc(3), 'XData', xd, 'YData', yd);
0485   yd = [get(hc(4),'YData') nfe(1:n)];
0486   set(hc(4), 'XData', xd, 'YData', yd);
0487   yd = [get(hc(5),'YData') nst(1:n)];
0488   set(hc(5), 'XData', xd, 'YData', yd);
0489 end
0490 
0491 drawnow;
0492 
0493 %-------------------------------------------------------------------------
0494 
0495 function [] = graphical_final(hfg,npg,stats,cntr)
0496 
0497 figure(hfg);
0498 pl = 0;
0499 
0500 if stats
0501   pl = pl+1;
0502   subplot(npg,1,pl)
0503   hc = get(gca,'Children');
0504   xd = get(hc,'XData');
0505   set(gca,'XLim',sort([xd(1) xd(end)]));
0506   
0507   pl = pl+1;
0508   subplot(npg,1,pl)
0509   ylim = get(gca,'YLim');
0510   ylim(1) = ylim(1) - 1;
0511   ylim(2) = ylim(2) + 1;
0512   set(gca,'YLim',ylim);
0513   set(gca,'XLim',sort([xd(1) xd(end)]));
0514 end
0515 
0516 if cntr
0517   pl = pl+1;
0518   subplot(npg,1,pl)
0519   hc = get(gca,'Children');
0520   xd = get(hc(1),'XData');
0521   set(gca,'XLim',sort([xd(1) xd(end)]));
0522   legend('nst','nfe','nni','netf','ncfn',2);
0523 end
0524 
0525 %-------------------------------------------------------------------------
0526 
0527 function [] = text_init(n,hft,stats,cntr,t,h,q,nst,nfe,nni,netf,ncfn)
0528 
0529 fig_name = 'CVODES run statistics';
0530 
0531 % If this is a parallel job, look for the MPI rank in the global
0532 % workspace and append it to the figure name
0533 
0534 global sundials_MPI_rank
0535 
0536 if ~isempty(sundials_MPI_rank)
0537   fig_name = sprintf('%s (PE %d)',fig_name,sundials_MPI_rank);
0538 end
0539 
0540 figure(hft);
0541 set(hft,'Name',fig_name);
0542 set(hft,'color',[1 1 1]);
0543 set(hft,'MenuBar','none');
0544 set(hft,'Resize','off');
0545 
0546 % Create text box
0547 
0548 margins=[10 10 50 50]; % left, right, top, bottom
0549 pos=get(hft,'position');
0550 tbpos=[margins(1) margins(4) pos(3)-margins(1)-margins(2) ...
0551        pos(4)-margins(3)-margins(4)];
0552 tbpos(tbpos<1)=1;
0553 
0554 htb=uicontrol(hft,'style','listbox','position',tbpos,'tag','textbox');
0555 set(htb,'BackgroundColor',[1 1 1]);
0556 set(htb,'SelectionHighlight','off');
0557 set(htb,'FontName','courier');
0558 
0559 % Create table head
0560 
0561 tpos = [tbpos(1) tbpos(2)+tbpos(4)+10 tbpos(3) 20];
0562 ht=uicontrol(hft,'style','text','position',tpos,'tag','text');
0563 set(ht,'BackgroundColor',[1 1 1]);
0564 set(ht,'HorizontalAlignment','left');
0565 set(ht,'FontName','courier');
0566 newline = '   time         step      order  |    nst   nfe   nni  netf  ncfn';
0567 set(ht,'String',newline);
0568 
0569 % Create OK button
0570   
0571 bsize=[60,28];
0572 badjustpos=[0,25];
0573 bpos=[pos(3)/2-bsize(1)/2+badjustpos(1) -bsize(2)/2+badjustpos(2)...
0574       bsize(1) bsize(2)];
0575 bpos=round(bpos);
0576 bpos(bpos<1)=1;
0577 hb=uicontrol(hft,'style','pushbutton','position',bpos,...
0578              'string','Close','tag','okaybutton');
0579 set(hb,'callback','close');
0580 
0581 % Save handles
0582 
0583 handles=guihandles(hft);
0584 guidata(hft,handles);
0585 
0586 for i = 1:n
0587   newline = '';
0588   if stats
0589     newline = sprintf('%10.3e   %10.3e     %1d    |',t(i),h(i),q(i));
0590   end
0591   if cntr
0592     newline = sprintf('%s %5d %5d %5d %5d %5d',...
0593                       newline,nst(i),nfe(i),nni(i),netf(i),ncfn(i));
0594   end
0595   string = get(handles.textbox,'String');
0596   string{end+1}=newline;
0597   set(handles.textbox,'String',string);
0598 end
0599 
0600 drawnow
0601 
0602 %-------------------------------------------------------------------------
0603 
0604 function [] = text_update(n,hft,stats,cntr,t,h,q,nst,nfe,nni,netf,ncfn)
0605 
0606 figure(hft);
0607 
0608 handles=guidata(hft);
0609 
0610 for i = 1:n
0611    if stats
0612     newline = sprintf('%10.3e   %10.3e     %1d    |',t(i),h(i),q(i));
0613   end
0614   if cntr
0615     newline = sprintf('%s %5d %5d %5d %5d %5d',...
0616                       newline,nst(i),nfe(i),nni(i),netf(i),ncfn(i));
0617   end
0618   string = get(handles.textbox,'String');
0619   string{end+1}=newline;
0620   set(handles.textbox,'String',string);
0621 end
0622 
0623 drawnow
0624 
0625 %-------------------------------------------------------------------------
0626 
0627 function [] = sol_init(n, hfs, nps, sol, sensi, dir, N, Ns, t, y, ys)
0628 
0629 fig_name = 'CVODES solution';
0630 
0631 % If this is a parallel job, look for the MPI rank in the global
0632 % workspace and append it to the figure name
0633 
0634 global sundials_MPI_rank
0635 
0636 if ~isempty(sundials_MPI_rank)
0637   fig_name = sprintf('%s (PE %d)',fig_name,sundials_MPI_rank);
0638 end
0639 
0640 
0641 figure(hfs);
0642 set(hfs,'Name',fig_name);
0643 set(hfs,'color',[1 1 1]);
0644 
0645 % Time label
0646 
0647 if dir==1
0648   tlab = '\rightarrow   t   \rightarrow';
0649 else
0650   tlab = '\leftarrow   t   \leftarrow';
0651 end
0652 
0653 % Get number of colors in colormap
0654 map = colormap;
0655 ncols = size(map,1);
0656 
0657 % Initialize current subplot counter
0658 pl = 0;
0659 
0660 if sol
0661 
0662   pl = pl+1;
0663   subplot(nps,1,pl);
0664   hold on;
0665 
0666   for i = 1:N
0667     hp = plot(t(1:n),y(i,1:n),'-');
0668     ic = 1+(i-1)*floor(ncols/N);
0669     set(hp,'Color',map(ic,:));
0670   end
0671   box on;
0672   grid on;
0673   xlabel(tlab);
0674   ylabel('y');
0675   title('Solution');
0676 
0677 end
0678 
0679 if sensi
0680   
0681   for is = 1:Ns
0682     
0683     pl = pl+1;
0684     subplot(nps,1,pl);
0685     hold on;
0686       
0687     ys_crt = ys(:,is,1:n);
0688     for i = 1:N
0689       hp = plot(t(1:n),ys_crt(i,1:n),'-');
0690       ic = 1+(i-1)*floor(ncols/N);
0691       set(hp,'Color',map(ic,:));
0692     end
0693     box on;
0694     grid on;
0695     xlabel(tlab);
0696     str = sprintf('s_{%d}',is); ylabel(str);
0697     str = sprintf('Sensitivity %d',is); title(str);
0698     
0699   end
0700   
0701 end
0702 
0703 
0704 drawnow;
0705 
0706 %-------------------------------------------------------------------------
0707 
0708 function [] = sol_update(n, hfs, nps, sol, sensi, N, Ns, t, y, ys)
0709 
0710 figure(hfs);
0711 
0712 pl = 0;
0713 
0714 if sol
0715   
0716   pl = pl+1;
0717   subplot(nps,1,pl);
0718   
0719   hc = get(gca,'Children');
0720   xd = [get(hc(1),'XData') t(1:n)];
0721 % Attention: Children are loaded in reverse order!
0722   for i = 1:N
0723     yd = [get(hc(i),'YData') y(N-i+1,1:n)];
0724     set(hc(i), 'XData', xd, 'YData', yd);
0725   end
0726 
0727 end
0728   
0729 if sensi
0730   
0731   for is = 1:Ns
0732     
0733     pl = pl+1;
0734     subplot(nps,1,pl);
0735 
0736     ys_crt = ys(:,is,:);
0737     
0738     hc = get(gca,'Children');
0739     xd = [get(hc(1),'XData') t(1:n)];
0740 %   Attention: Children are loaded in reverse order!
0741     for i = 1:N
0742       yd = [get(hc(i),'YData') ys_crt(N-i+1,1:n)];
0743       set(hc(i), 'XData', xd, 'YData', yd);
0744     end
0745     
0746   end
0747   
0748 end
0749 
0750 
0751 drawnow;
0752 
0753 
0754 %-------------------------------------------------------------------------
0755 
0756 function [] = sol_final(hfs, nps, sol, sensi, N, Ns)
0757 
0758 figure(hfs);
0759 
0760 pl = 0;
0761 
0762 if sol
0763 
0764   pl = pl +1;
0765   subplot(nps,1,pl);
0766   
0767   hc = get(gca,'Children');
0768   xd = get(hc(1),'XData');
0769   set(gca,'XLim',sort([xd(1) xd(end)]));
0770 
0771   ylim = get(gca,'YLim');
0772   addon = 0.1*abs(ylim(2)-ylim(1));
0773   ylim(1) = ylim(1) + sign(ylim(1))*addon;
0774   ylim(2) = ylim(2) + sign(ylim(2))*addon;
0775   set(gca,'YLim',ylim);
0776   
0777   for i = 1:N
0778     cstring{i} = sprintf('y_{%d}',i);
0779   end
0780   legend(cstring);
0781   
0782 end
0783 
0784 if sensi
0785   
0786   for is = 1:Ns
0787     
0788     pl = pl+1;
0789     subplot(nps,1,pl);
0790 
0791     hc = get(gca,'Children');
0792     xd = get(hc(1),'XData');
0793     set(gca,'XLim',sort([xd(1) xd(end)]));
0794 
0795     ylim = get(gca,'YLim');
0796     addon = 0.1*abs(ylim(2)-ylim(1));
0797     ylim(1) = ylim(1) + sign(ylim(1))*addon;
0798     ylim(2) = ylim(2) + sign(ylim(2))*addon;
0799     set(gca,'YLim',ylim);
0800   
0801     for i = 1:N
0802       cstring{i} = sprintf('s%d_{%d}',is,i);
0803     end
0804     legend(cstring);
0805     
0806   end
0807   
0808 end
0809 
0810 drawnow