mode(-1);

getf('Ordre.sci');
getf('aren.sci');
%ODEOPTIONS=[2,0,0,%inf,0,2,500,12,5,0,-1,-1];

ok=1;
while ok do
  str=x_mdialog(['Résolution d''un problème à trois corps'; ...
	  'Comparaison de méthodes à pas fixes et d''une méthode avec contrôle du pas'], ...
      ['Ordre (0 pour variable)'; 'Nombre de pas de temps'; 'Temps final'], ...
      ['4'; '5000'; '17.1']);

  ok = (str~=[]);
  if ok then
    ordre=evstr(str(1)); 
    n=evstr(str(2));
    tf=evstr(str(3));
    t0=0;
    y0=[0.994; 0; 0; -2.0015851063791];
    h=tf/n;
  
    if ordre ~=0 then
      [tvals, yvals] = FixedRK(aren, t0, y0, h, ordre, n);
      yvref= ode(y0, t0, tf, aren);
      xbasc();
      plot2d(yvals(1,:), yvals(2,:), style=1);
      plot2d(yvref(2,:), yvref(3,:), style=2, leg='Pas variable');
      xtitle('Solutions du problème à trois corps');
    else
      yvref= ode(y0, t0, tf, aren);      
      nts= size(yvref,2);
      x_message(['évolution du pas de temps et de la solution'; string(nts)+' pas de temps']);    
      xbasc();
      dts=yvref(1,2:nts)-yvref(1,1:nts-1);
      xsetech([0. 0. 1. 0.33333333]);
      plot2d(dts);
      xsetech([0. 0.3333333 1. 0.33333333]);
      plot2d(yvref(2,:));
      xsetech([0. 0.6666666 1. 0.33333333]);
      plot2d(yvref(3,:));
    end
  end
end