! un programme d'illustration de l'interpolation de Lagrange
! on donne les points par lesquels l'interpolant doit passer: x et f
! et on evalue l'interpolant en tous points (une boucle sur x)

program Lagrange
implicit none
integer, parameter :: n = 3
integer :: i,j
real :: xinterp, p
real, dimension(0:n) :: x, f, l

! les points par lesquels doit passer l'interpolant
x=(/-2,0,4,6/) ! abcisse
f=(/3,5,8,5/)  ! ordonnée

! on ouvre le fichier dans lequel on sauve les resultats
open(unit=100,file="interp.dat")

! boucle sur xinterp, pour evaluer la valeur de l'interpolant et en faire un graph
! on va de -3 a 7 pour couvrir les abcisses interessantes
xinterp=-3
do while (xinterp<7)
   p=0
   do j=0,n
      l(j)=1 ! initialisation de l pour ce point
      do i=0,n

         if (i/=j) then ! verification que i est different de j
            l(j)=l(j)*(xinterp-x(i))/(x(j)-x(i)) ! on ajoute un terme au produit
         end if

      end do
      p=p+f(j)*l(j) ! on ajoute la contribution de ce polynome a l'interpolant
   end do

   ! on ecrit le point x et la valeur de l'interpolant dans le fichier
   write(100,*) xinterp, p, l

   xinterp=xinterp+0.01 ! on incremente xinterp de 0.01
end do

! fermeture du fichier dans lequel on a ecrit les resultats
close(100)

end program