/*---------------------------------------------------*/
/* Résolution numériques d'équations différentielles */
/* ----------- Méthode d'euler ----------------------*/
/* ----------- Laure Gonnord, Novembre 2004 ---------*/
/*---------------------------------------------------*/

/*compiler avec gcc -L /usr/X11R6/lib -lm -lX11 -lrt -L . -lgraphdsu euler.c -o euler*/
/*ou, si on est sur thales ou aspasi : gccgraph -o euler euler.c*/

#include <stdio.h>
#include "graphdsu.h" /*pour la visualisation graphique*/

#define MAX 1000 /*max du nb de points*/
/*taille de la fenêtre graphique*/
#define TAILLEX 600 
#define TAILLEY 600

/*la fonction f*/
double f(double x)
{
  return (x);
}

/*résolution d'une équa diff avec euler*/
/*on suppose que la fonction f existe et est écrite plus haut*/
/*avec ce qu'on sait, on ne peut pas la passer en paramètre à une fonction*/
/*y'=f(y), avec CI : y(t0)=y0, nb pas = N, taille pas p */
/* on stocke les points dans abs/ord*/
void euler(double t0, double y0, int N, double p, double abs[MAX], double ord[MAX])
{
  int k;
  abs[0] = t0;
  ord[0] = y0;

  for (k=1;k<N;k++)
    {
      abs[k] = t0 + ((double) k)*p;
      ord[k] = ord[k-1] + f(ord[k-1])*p;
    }

  /*impression des deux tableaux*/
    for (k=0;k<N;k++)
    {
      printf(" %lf ",abs[k]);
    }
    printf("\n");
    for (k=0;k<N;k++)
    {
      printf(" %lf ",ord[k]);
    }
    printf("\n");
    /*fin impression*/
}

/*dessin d'une fonction affine par morceaux*/
/*on donne le tableau d'ascisses, le tableau d'ordonnees*/
/*le nb de points (N<MAX), le range pour x (xmin/xmax), et pour y*/
/*on suppose que la fenetre graphique est ouverte*/
void joli_dessin(double abs[MAX], double ord[MAX], int N, int xmin, int xmax, int ymin, int ymax)
{
  int i,x1,x2,y1,y2;
  /*calcul du repere*/
  double ratiox = TAILLEX/(xmax-xmin);
  double ratioy = TAILLEY/(ymax-ymin);

  ChangeCouleur(3);

  /*dessin des segments*/
  for (i=0;i<N-1;i++)
    {
      /*calcul des coordonnées dans le dessin*/
      x1=(int) (ratiox * abs[i]); 
      y1=(int) (TAILLEX - (ratioy * ord[i]));
/*tableau ord --> ordonnee y1 --> ordonnee dans le dessin*/
/*attention axe ordonnee dans la mauvaise direction*/
      x2=(int) (ratiox * abs[i+1]);
      y2=(int) (TAILLEY -(ratioy * ord[i+1]));

      /*dessin du segment*/
      //      printf("ligne (%d,%d)---(%d,%d)\n",x1,y1,x2,y2);
      Ligne(x1,y1,x2,y2);
    }

}

int main (void)
{
  double abs[MAX];
  double ord[MAX];
  double d = 0.01;
  int N = 100;
  double t0 = 0;
  double y0 = 1;

  euler(t0,y0,N,d,abs,ord);
  Initialiser(TAILLEX,TAILLEY);

  /*pour f = id, penser à changer les bornes pour d'autres f*/
  joli_dessin(abs,ord,N,0,3,0,3);

  Cliquer();
  Clore();

  return 0;
}