Programmation en C

dp · January 2023

Oui, oui. Excuse moi d'être en weekend et fatigué.

parisse · January 2023

pas de soucis, l'essentiel c'est de ne pas laisser de code inefficace qui risquerait d'entrainer chatgpt.

dp · January 2023

Euh… comment dire que j'en ai rien à fiche de chatGPT ? Je ne bosse pas pour OpenAI.

dp · February 2023

Bon, c’est un peu après la bataille, mais je me dis que @gebrane ne sera peut-être pas le seul à ne pas savoir comment utiliser le code que j’ai mis en première page de cette discussion. J’ai donc rapidement pris le temps de rajouter le support de getopt_long afin d’avoir des arguments nommés. J’en ai profité pour rajouter une fonction usage, de sorte qu’à chaque fois que le programme ne reçoit pas les arguments voulus, il affiche son fonctionnement. En espérant que ça soit suffisant pour que plus personne ne se retrouve dépourvu devant celui-ci.

gebrane · February 2023

Voici Un TP avec un problème similaire
Question 1 Suivez l'algorithme ci-dessous pour programmer en C l'interpolation de Lagrange avec n points entrés par l'utilisateur en un point x donné également entré par l'utilisateur.
Algorithme

    Entrez les n+1 points (x0, y0), (x1, y1),.... (xn, yn)
    Entrez la valeur x pour laquelle vous souhaitez évaluer la fonction interpolée
    Initialisation P(x)=0.
    Pour chaque i de 0 à n :
    a. Initialisation L_i(x) = 1.
    b. Pour chaque j de 1 à n, avec j ≠ i :
    c. L_i(x) = L_i(x)+(x - xj) / (xi - xj);
    d. P(x)=P(x)+yi*L_i(x);

    Sortie : Le polynôme interpolant P en x est :

Voici mon programme en C ( Il n'est pas rapide, mais ce n'est pas le but pour le moment)

    #include <stdio.h>
    #include <math.h>

    double Interpolation_Lagrange(double x, double x_values[], double y_values[], int n)
    {
         double P_x = 0;
         for (int i = 0; i < n; i++)
        {
           double Li = 1;
           for (int j = 0; j < n; j++)
        {
            if (j != i)
            {
                Li *= (x - x_values[j])/(x_values[i] - x_values[j]);
            }
        }
          P_x += y_values[i] * Li;
        }
        return P_x;
    }

    int main()
    {
        int n;
        printf("Entrer le nombre de points n : ");
        scanf("%d", &n);

        double x_values[n];
        double y_values[n];

        printf("Entrez (l'abscisse suivie de la touche espace puis l'ordonnée suivie de la touche Entrée ) les n points (xi, yi) pour i=0 à n-1 : \n");
        for (int i = 0; i < n; i++)
        {
            scanf("%lf%lf", &x_values[i], &y_values[i]);
        }

        double x;
        printf("Entrez la valeur x pour laquelle vous souhaitez évaluer la fonction interpolée: ");
        scanf("%lf", &x);

        printf("L'interpolation en x est : %lf", Interpolation_Lagrange(x, x_values, y_values, n));

        return 0;
    }

Question 2 Modifier le programme pour afficher le polynôme de Lagrange
Pour la question 2, il est nécessaire d'utiliser une bibliothèque, j'opte pour tinyexprs, travail en cours.

dp · February 2023

EDIT: oups, je ne sais pas pourquoi j'ai mis un pointeur à la ligne 33, sûrement le fait qu'il est tard. Comme tu connais le nombre de points, il faut plutôt mettre

point_t points[n] ;

Considère l'utilisation de structures pour te faciliter la vie et rendre ton code un peu plus naturel.

#include <stdio.h>
#include <math.h>

typedef struct point_t {
   double x;
   double y;
} point_t;

double Interpolation_Lagrange(double x, point_t points[], int n)
{
     double P_x = 0;
     for (int i = 0; i < n; i++)
    {
       double Li = 1;
       for (int j = 0; j < n; j++)
    {
        if (j != i)
        {
            Li *= (x - points[j].x)/(points[i].x - points[j].x);
        }
    }
      P_x += points[i].y * Li;
    }
    return P_x;
}

int main()
{
    int n;
    printf("Entrer le nombre de points n : ");
    scanf("%d", &n);

    point_t points[n] ;

    printf("Entrez (l'abscisse suivie de la touche espace puis  l'ordonnée suivie de la touche Entrée ) les n points  (xi, yi)  pour i=0 à n-1 : \n");
    for (int i = 0; i < n; i++)
    {
        scanf("%lf%lf", &points[i].x, &points[i].y);
    }

    double x;
    printf("Entrez la valeur x pour laquelle vous souhaitez évaluer la fonction interpolée: ");
    scanf("%lf", &x);

    printf("L'interpolation en x est : %lf", Interpolation_Lagrange(x, points, n));

    return 0;
}

gebrane · February 2023

Merci pour l'astuce @dp

dp · February 2023

Remarque, je ne sais pas si tu comptes faire écrire ces programmes à tes étudiants. Dans le doute, si jamais c’est le cas, tu peux leur cacher certaines choses, comme la création des structures dans des fichiers header (.h). Ils n’auront qu’à importer ces headers et utiliser tranquillement la structure sans jamais avoir à se soucier de comment elle est créée.

gebrane · February 2023

Bonne remarque, je vais voir pour le comment. Merci

Alain24 · February 2023

Rescassol

&nbsp est le résultat de deux OS (operating system en anglo-améraicain et en français système d'exploitation) incompatibles.

Si on utilise un bon vieux PC (personal computer en anglo-américain et en français ordinateur personnel) ces symboles apparaissent lors de tentatives de traduction LinuxIsNotUniX versus Windows.

parisse · February 2023

Il n'est pas nécessaire d'avoir une bibliothèque de polynômes pour faire de l'interpolation de Lagrange. On peut stocker les polynômes sous la forme de la liste de leurs coefficients (ce sera plus facile avec la STL en C++ avec un vector<double>).

Mais en fait ce n'est pas la bonne méthode. La bonne méthode c'est de calculer une fois pour toutes et de stocker la liste des différences divisées et de faire une évaluation "à la Horner" pour calculer la valeur de P(x) à partir des différences divisées et des x_i et de x.
https://www-fourier.univ-grenoble-alpes.fr/~parisse/doc/fr/algo.html#sec302

Sato · February 2023

Il serait peut-être bon de séparer les entrées sorties au clavier de calcul proprement dit dans le programme. Tant qu’à faire. Il faut qu’on puisse l’utiliser avec un pipe ou un fichier de données de test. C’est un minimum.

gebrane · February 2023

Avec l’entête du fichier de @dp, un collègue à réussi la bonne syntaxe avec tinyexprs pour trouver le polynôme de Lagrange.

dp · February 2023

Ça fait mal au cœur de voir ça, mais tant que c'est fonctionnel j'imagine que c'est suffisant.

Toutefois je me demande bien pourquoi le C est utilisé dans ce contexte. Tu as l'explication @gebrane ?

gebrane · February 2023

Bonjour parisse

parisse a dit : Il n'est pas nécessaire d'avoir une bibliothèque de polynômes pour faire de l'interpolation de Lagrange.

Je n'ai pas dit le contraire. Mais le C ne sait pas donner le polynôme sous forme mathématiques. On a besoin d'une bibliothèque (non pas des polynômes) qui permet de donner des résultats sous forme de fonctions

parisse a dit : Mais en fait ce n'est pas la bonne méthode

Oui on le sait, l'interpolation de Lagrange est mauvaise si le nombre de points est grand.
Tous cela viendra par la suite step by step

gebrane · February 2023

@dp c'est pour montrer les limites de C peut-être. Une vraie torture.

Je ne vais plus vous embêter et je vous souhaite une bonne journée à tous et à toutes.

dp · February 2023

Au fait dans le fichier que tu as envoyé ici il n’y a aucune utilisation de tinyexpr. Pourquoi en avoir inclus les près de 1000 lignes ?

#include <stdio.h>
#include <string.h>


double x_values[100];
double y_values[100];
char expression[1000];
int n;

void calculate_expression(int i)
{
sprintf(expression, "y%d = ", i);
for (int j = 0; j < n; j++)
{
if (j == i)
{
continue;
}
char tmp[100];
sprintf(tmp, "(x - %g)/(%g - %g)", x_values[j], x_values[i], x_values[j]);
strcat(expression, tmp);
if (j != n-1)
{
strcat(expression, " * ");
}
}
strcat(expression, " * ");
char tmp[100];
sprintf(tmp, "%g", y_values[i]);
strcat(expression, tmp);
}

int main()
{
printf("Entrer le nombre de points n : ");
scanf("%d", &n);



printf("Entrer les n points  (xi, yi)  pour i=0 Ã  n-1 : \n");
for (int i = 0; i < n; i++)
{
    scanf("%lf%lf", &x_values[i], &y_values[i]);
}

strcpy(expression, "");
for (int i = 0; i < n; i++)
{
    calculate_expression(i);
    if (i != n-1)
    {
        strcat(expression, " + ");
    }
}
printf("Le polynome P(x) est : %s\n", expression);

return 0;

}

gebrane · February 2023

Ah! Surprise. Puisque ça ne venait pas de moi, je vais voir.

merci pour cette remarque pertinente.

dp · February 2023

En montrer les limites ça me parait tordu. Le langage C n’a aucune réelle limite hormis notre imagination.

À la limite il serait peut-être pertinent de montrer un langage de programmation ultra rapide afin de comparer à Python. Et encore, même dans ce cas-là je trouve ça idiot d’utiliser le langage C, qui demande beaucoup (trop) d’efforts, alors qu’il existe go ou même, dans une autre mesure, rust qui seraient beaucoup plus adaptés à ton public. De plus, même si c’est pour faire plonger un orteil dans le bain de la programmation système à tes étudiants, le langage rust reste toujours beaucoup plus pertinent que le langage C dans le contexte d’un cours universitaire. M’enfin, le choix est de toute façon acté pour cette année. Toutefois pensez-y pour l’année prochaine !

gebrane · February 2023

Mais ce sont les premiers TP, après c'est un autre monde je ne dirais pas plus . Je vais me plonger sur rust car je ne le connais pas

gebrane · February 2023

@dp ou autres ( ce n'est pas une obligation ) , peux-tu me donner un programme en rust qui donne le polynôme de Lagrange pour voir que c'est plus intéressant que le C. D’après mes premières lectures, ce n'est pas évident la syntaxe. Merci

dp · February 2023

Tu as un exemple ici.

gebrane · February 2023

Bonjour, ce programme affiche un résultat bizarre, il est mal conçu

Dans ce cas précis, je préfère le C, la syntaxe de rust est aussi compliqué

Image: https://les-mathematiques.net/vanilla/uploads/editor/yl/6hn5xyou5nig.jpg

Sato · February 2023

Il serait peut-être bon de séparer les entrées sorties au clavier de calcul proprement dit dans le programme. Tant qu’à faire. Il faut qu’on puisse l’utiliser avec un pipe ou un fichier de données de test. C’est un minimum.

dp · February 2023

Ce n'est pas la syntaxe de rust qui est compliquée dans le cas présent, mais bien ton programme qui l'est inutilement. Nuance.

dp · February 2023

Tu sais, le but quand tu changes de langage de programmation n’est pas de refaire la même chose que dans les précédents, mais de s’adapter à ce nouveau langage et ses paradigmes. Sinon autant rester sur l’ancien.

gebrane · February 2023

@dp Mon programme C donne le polynôme d'une manière correcte et le programme rust du site que tu cites donne un résultat bizarre.

Pourquoi mon programme est inutilement compliqué, peux-tu faire mieux ou c'est la question que tu n'aimes pas ( Donner le polynôme en question et non l’évaluation du polynôme en un point) ) En maths, j'ai vu des complications plus pires et cela n'a décourager personne.

@sato Les données seront appelés d'une base de données mais la , on donne le programme basique

dp · February 2023

Je n’ai jamais dit que celui que je te donnais était parfait. En réalité j’ai juste tapé « Lagrange interpolation rust » sur Google est c’est le premier lien qui est tombé. Rien ne t’empêche de trouver mieux.

Ton programme est inutilement compliqué, car tu veux faire du C en rust. c’est débile. Le langage rust a été conçu, entre autres, pour s’affranchir d’un bon nombre de contraintes du C. Par exemple, quel est l’intérêt (en rust) de la ligne

let mut expression = String::with_capacity(MAX_EXPRESSION_LEN) ;

pour un programme aussi simple ?

PS. Je parle bien dans le cas présent que tu cites être « le programme basique ».

Enfin, arrête de croire que la programmation c’est des mathématiques appliquées ou que sais-je, car justement

En maths, j'ai vu des complications plus pires et cela n'a décourager personne.

ceci est la raison pour laquelle il existe des ingénieurs pour faire des choses propres et efficaces

parisse · February 2023

@gebrane. Écrire un polynôme d'interpolation dans la base (1,x,...,x^n) et l'évaluer ensuite comme une expression symbolique est mauvais non seulement d'un point de vue efficacité, mais à cause de la perte de précision numérique que cela peut entrainer. Par exemple le polynôme d'interpolation P:=lagrange([1000,1001,1002],[1,0,2]) en 1000, 1001 et 1002 valant 1 puis 0 puis 2 est dans la base de Newton (x-1000)*(3*(x-1001)/2-1)+1 et dans la base canonique 3/2*x^2-6005/2*x+1502501. Si on l'évalue en x=1000.6, dans le premier cas on perd environ 3 chiffres significatif, dans le 2ème cas on en perd environ 6 (en raison de différences de nombre proches).
Sur l'interpolation mauvaise, si on utilise trop de points d'interpolation, il faut bien sûr interpoler par intervalles en utilisant des points de Tchebycheff. Mais pour avoir une précision raisonnable en précision double et pas trop d'intervalles il faut quand même une vingtaine de points. Donc l'algorithme des différences divisées et le calcul du polynôme dans la base de Newton s'impose. L'écriture dans la base canonique est je pense un biais dans notre utilisation usuelle de polynômes...

gebrane · February 2023

@parisse. Bonjour et merci pour ton retour, je crois qu'on ne se comprend pas sur le but (gradué). On commence l’interpolation avec les polynômes de Lagrange. Une première étape c'est de maitriser l'algorithme et la programmation, puis on montre son inefficacité lorsque le nombre de points est grand, après on atteint l'algorithme de Neville-Aitken qui est plus efficace.

@dp. Je ne comprends pas cette phrase "Ton programme est inutilement compliqué car tu veux faire du C en rust. c'est débile."

Qui te dit que je voulais faire du C en rust. Si tu relis bien les messages j'ai proposé un programme en C, tu es intervenu pour me dire que le langage le plus adapté est rust ; pour comprendre : je t'ai demandé (et ce n'est pas une obligation) comment résoudre la question avec rust pour voir la différence.

dp · February 2023

Le problème n’est pas de vouloir adapter ton programme c en rust mais plutôt de vouloir programmer en rust comme on programme en c (donc de manière inadéquate) et venir dire que la syntaxe (ici c’est cette dernière mais ça pourrait être autre chose) est compliquée. Essaie réellement de comprendre comment rust se programme et tu verras que tu peux dans 95% des cas produire du code à la fois court, compréhensif et efficace.

gebrane · February 2023

C'est ce que je vais faire avec rust ( je commence à maîtriser C)

Merci et bonne journnée.

dp · February 2023

je commence à maitriser C

Cette phrase est la preuve que non.

gebrane · February 2023

Je mets la négation de ma proposition comme tu le réclames, cela donne : je ne commence pas à maîtriser C

dp · February 2023

Voici quelles que ressources si tu souhaites réellement t’intéresser à rust :

la liste de tous les manuels : https://lborb.github.io/book/
la liste de toutes les bibliothèques (tierces) : https://lib.rs/
la liste de celles qui sont "bénies" : https://blessed.rs/crates
Un ancien étudiant du MIT (travaillant maintenant chez Amazon) qui fait beaucoup de contenu orienté rust : https://www.youtube.com/@JonGjengset
Le meilleur livre pour commencer : https://www.oreilly.com/library/view/programming-rust-2nd/9781492052586/
Le meilleur livre pour appliquer : https://www.manning.com/books/rust-in-action
Le meilleur livre pour s'aventurer dans les arcanes : https://nostarch.com/rust-rustaceans

gebrane · February 2023

Merci, en plus je suis fan de jeu rust ( je ne sais pas s'il y a un lien ) . Je te remercie au nom de tous ceux et celles qui vont utiliser ces références. .

dp · February 2023

https://www.youtube.com/shorts/C9he5SmNkyE

parisse · February 2023

@gebrane : montrer l'existence et l'unicité du polynôme d'interpolation en ajoutant les points un par un (et donc la base de Newton) est je trouve la méthode la plus naturelle. Et son évaluation sous cette forme est tout aussi naturelle.

La seule chose qui est plus avancée c'est qu'on peut faire le calcul des coefficients avec les différences divisées, plutôt qu'avec la formule qui vient naturellement dans la construction par récurrence.

gebrane · February 2023

Bonsoir @parisse Merci pour ton retour. Tu me donnes une idée, je vais combiner les deux méthodes dans un seul programme.

Je vais considérer pour exemple le segment [0,10] que je vais subdiviser en $x_i$ avec $x_0=0$ et $x_n=10$ avec un pas h=0.1 (donc 101 points) on va prendre pour $y_i$ les images des $x_i$ par une fonction $f$, Je vais demander d’implémenter dans un même programme l'algorithme de Lagrange et de Newton (avec différences divisées). Le programme doit tracer la fonction f et le polynôme d'interpolation des deux méthodes.

gebrane · February 2023

Pour f(x) =x², j'obtiens des résultats bizarres : Newton est meilleur avant $x_{20}$, mais après c'est Lagrange qui devient meilleur.

parisse · February 2023

C'est probablement une conséquence de la précision des flottants en précision double et de la taille des coefficients qui apparaissent. Si on fait les calculs avec disons 70 décimales, les différences divisées sont nettement plus précises que l'évaluation du polynôme développé (via Horner). Cf. https://xcas.univ-grenoble-alpes.fr/forum/viewtopic.php?f=30&t=2854

gebrane · February 2023

Bonjour @parisse, je te mets en lien mon programme avec implémentation des deux méthodes Lagrange et Newton (avec différences divisées) et j'ai ajouté la valeur exacte et la différence entre la méthode et la valeur exacte. Je vais regarder ton lien.

Merci edit j'avais écrasé le bon fichier et laisser le mauvais par inattention

Après je vais essayer d’implémenter les 4 méthodes d'interpolation polynomiale dans un même programme.

dp · February 2023

Je m’ennuyais du coup je me suis demandé comment faire tes deux TP en rust. J’ai voulu me mettre dans la peau de tes étudiants, du coup je suis allé au plus simple : ce n’est pas du grand rust – juste du rust fonctionnel et simple à comprendre – et la seule notion un peu avancée que je vais utiliser ici étant la programmation orientée objets.

Commençons donc ici par le premier TP. Pour ce TP, tu as besoin de deux bibliothèques tierces : clap et meval. Tu peux très facilement les installer en tapant

cargo add clap
cargo add meval

dans ton terminal.

Pour ce qui est du code en lui-même, je l’ai réparti en deux fichiers : integral.rs et main.rs.

integral.rs ressemble à ceci :

#[derive(Debug)]
pub(crate) struct Integral {
    pub(crate) left_rect: f64,
    pub(crate) right_rect: f64,
    pub(crate) middle_rect: f64,
    pub(crate) trapeze: f64,
    pub(crate) simpson: f64
}

impl Integral {
    pub fn new() -> Self {
        Self { left_rect:   0.0,
               right_rect:  0.0,
               middle_rect: 0.0,
               trapeze:     0.0,
               simpson:     0.0
             }
    }

    pub fn integrate(&mut self, f: &str, inf: f64, sup: f64, n: usize) {
        let expression: meval::Expr = f.parse().unwrap();
        let function = expression.bind("x").unwrap();

        self.left_rect   = Self::integrate_left_rect(&function, inf, sup, n);
        self.right_rect  = Self::integrate_right_rect(&function, inf, sup, n);
        self.middle_rect = Self::integrate_middle_rect(&function, inf, sup, n);
        self.trapeze     = Self::integrate_trapeze(self);
        self.simpson     = Self::integrate_simpson(self);
    }

    fn integrate_left_rect(f: impl Fn(f64) -> f64, a: f64, b: f64, n: usize) -> f64 {
        let mut result = 0.0 ;
        let h = (b-a)/(n as f64) ;
        for i in 0..n { result += f(a+h*(i as f64)) ; }

        result*h
    }

    fn integrate_right_rect(f: impl Fn(f64) -> f64, a: f64, b: f64, n: usize) -> f64 {
        let mut result = 0.0 ;
        let h = (b-a)/(n as f64) ;
        for i in 1..n+1 { result += f(a+h*(i as f64)) ; }

        result*h
    }

    fn integrate_middle_rect(f: impl Fn(f64) -> f64, a: f64, b: f64, n: usize) -> f64 {
        let mut result = 0.0 ;
        let h = (b-a)/(n as f64) ;
        for i in 0..n { result += f(a+h*(0.5 + i as f64)) ; }

        result*h
    }

    fn integrate_trapeze(&self) -> f64 {
        (self.left_rect+self.right_rect)/2.0
    }

    fn integrate_simpson(&self) -> f64 {
        (self.left_rect+self.right_rect+4.0*self.middle_rect)/6.0
    }
}

main.rs quant-à lui ressemble à ceci :

mod integral;
use integral::Integral;

use clap::Parser ;

#[derive(Parser,Default)]
struct Args {
    #[arg(short, long)]
    function: String,
    #[arg(short, long)]
    inf: f64,
    #[arg(short, long)]
    sup: f64,
    #[arg(short)]
    n: usize
}

fn main() {
    let args = Args::parse() ;
    let mut i = Integral::new();
    i.integrate(args.function.as_str(), args.inf, args.sup, args.n);
    
    println!("gauche: {}",  i.left_rect) ;
    println!("droite: {}",  i.right_rect) ;
    println!("milieu: {}",  i.middle_rect) ;
    println!("trapeze: {}", i.trapeze) ;
    println!("simpson: {}", i.simpson) ;
}

Petite subtilité, dans cargo.toml, tu auras besoin de changer la ligne correspondante à la dépendance de clap par

clap = {version = "~4", features = ["derive"]}

Finalement, tu peux l'utiliser presque comme celui que j’ai écrit en C à la première page :

cargo run -- --function "e^(cos(x)+sin(x*x))-23*x^2+17*x-194" --inf 2 --sup 4 -n 10000

Je m’en vais faire le deuxième un peu plus tard dans la soirée.

parisse · February 2023

@gebrane: il y a un truc que je ne comprends pas dans votre code C, pour moi les différences divisées correspondent à un premier calcul (en O(n^2)) qui dépend des x_i et y_i mais ne dépend pas de x. Ensuite on peut évaluer le polynôme d'interpolation pour toutes les valeurs de x en O(n) opérations. C'est ce qui est fait dans le lien que j'ai donné (avec des flottants multiprécision pour conserver un minimum de précision sur le résultat).

dp · February 2023

Deuxième TP. Ce TP est amusant, car tu demandes d’obtenir un visuel. On peut donc s’amuser à voir ce que cela donne pour différents points et différents jeux de donnés.

Bref, pour ce TP je vais y aller avec deux messages. Dans le premier, celui-ci, on entrera les données en dur dans le programme; dans le second on cherchera une méthode simple pour les lires depuis un fichier. C’est tipar !

Pour ce TP on aura besoin de gnuplot et sscanf ! Pour les installer, c’est toujours aussi simple :

cargo add gnuplot
cargo add sscanf

Pour ce qui est du code en lui-même, je me suis basé sur ton fichier par manque de temps (je modifierais sans doute ce soir si 'y'a des problèmes). Tout a été fait en trois fichiers : point.rs, interpolation.rs et main.rs.

point.rs est très simple : il s’agit simplement d’avoir une structure Point

use crate::point::Point ;

use gnuplot::{Figure, Caption, Color} ;


pub(crate) struct Interpolation {
        data:   String,
    pub points: Vec<Point>
}


impl Interpolation {
    pub fn new(data: String) -> Self {
        Self {
            data,
            points: vec![]
        }
    }

    pub fn generate_points(&mut self) {
        for point in self.data.lines() {
            self.points.push(Point { x: sscanf::scanf!(point, "({},{})", f64, f64).unwrap().0,
                                     y: sscanf::scanf!(point, "({},{})", f64, f64).unwrap().1
                                   }
                            );
        }
    }

    pub fn interpolate(&self, p: f64) -> (f64, f64) {
        let lagrange = self.lagrange_interpolation(p);
        let newton   = self.newton_interpolation(p);

        (lagrange, newton)
    }

    pub fn plot(&mut self, xmin: i32, xmax: i32) {
        /* Generates arrays containing interpolated points */
        let mut xs: Vec<f64> = vec![];
        let mut yl: Vec<f64> = vec![];
        let mut yn: Vec<f64> = vec![];
        
        for i in (xmin*1000)..(xmax*1000) {
            let x = f64::from(i) * 0.001 ;
            let interpolated = self.interpolate(x) ;

            xs.push(x) ;
            yl.push(interpolated.0); // lagrange
            yn.push(interpolated.1); // newton
        }

        /* Generates the plot with gnuplot */
        let mut fig = Figure::new();
        fig.axes2d().lines(&xs, &yl, &[Caption("Interpolation de Lagrange"), Color("red")])
                    .lines(&xs, &yn, &[Caption("Interpolation de Newton"), Color("blue")]);

        /* Show the plot */
        fig.show();
    }

    fn lagrange_interpolation(&self, p: f64) -> f64 {
        let n: usize = self.points.len() ;
        let mut result: f64 = 0.0 ;

        for i in 0..n {
            let mut term: f64 = self.points[i].y ;

            for j in 0..n {
                if i != j {
                    let numerator   = p - self.points[j].x ;           // avoid magic
                    let denominator = self.points[i].x - self.points[j].x ; // avoid magic

                    term *= numerator / denominator ;
                }
            }
            result += term ;
        } 

        result
    }

    fn divided_differences(&self) -> Vec<f64> {
        let n: usize = self.points.len() ;
        let mut diffs: Vec<f64> = vec![self.points[0].y] ;

        for i in 1..n {
            let mut term = self.points[i].y ;

            for j in 0..i {
                let numerator   = term - diffs[j] ;
                let denominator = self.points[i].x - self.points[j].x ;
                term = numerator / denominator ;
            }

            diffs.push(term) ;
        }

        diffs
    }

    fn newton_interpolation(&self, p: f64) -> f64 {
        let d: Vec<f64> = self.divided_differences();
        let n = self.points.len() ;

        let mut result: f64 = 0.0 ;

        for i in 0..n {
            let mut term = d[i] ;

            for j in 0..i { 
                term *= p - self.points[j].x ;
            }

            result += term ;
        }

        result
    }
}

On remarquera que interpolate réaliser les interpolations en un seul point tandis que plot réalise le graphique à l’aide de gnuplot.

Finalement mais.rs est classique

mod interpolation;
mod point;
use crate::interpolation::Interpolation;


fn main() {
    let data="(0.0,0.0)
              (0.5,0.15)
              (1.0,0.30)
              (1.5,0.45)
              (2.0,0.60)";

    let mut interpol = Interpolation::new(String::from(data));
    interpol.generate_points() ;

    interpol.plot(0, 4);
}

Lorsque tu lances le programme, tu obtiens ce joli petit graphique :

dp · February 2023

Deuxième partie.

Cette fois-ci, on cherche à lire les données depuis un fichier.

On aura besoin de clap que l'on peu installer comme suit

cargo add clap

Ensuite, tout va se passer dans le fichier main.rs que l'on va modifier de sorte à obtenir le contenu suivant

mod interpolation;
mod point ;
use crate::interpolation::Interpolation;

use std::fs::File;
use std::io::{BufReader, Read};

use clap::Parser ;

#[derive(Parser,Default)]
struct Args {
    #[arg(short, long)]
    file: String,
    #[arg(short, long)]
    inf: i32,
    #[arg(short, long)]
    sup: i32
}

fn main(){
    let args = Args::parse() ;
    let mut input: BufReader = BufReader::new(File::open(args.file).expect("Could not open the file..."));
    let mut data = String::new();
    input.read_to_string(&mut data).expect("Could not read the data...");

    let mut interpol = Interpolation::new(data.to_string());
    interpol.generate_points() ;
    interpol.plot(args.inf, args.sup);
}

Et voilà !

Maintenant on peu tester. Si on entre le jeu de données suivant dans un fichier points.txt

(0.0,0.0)
(0.5,0.15)
(1.0,0.30)
(1.5,0.45)
(2.0,0.60)

et qu'on appelle notre programme avec la commande suivante

cargo r -- --file points.txt --inf 0 -- sup 2

on se retrouve avec ce graphique

Image: https://les-mathematiques.net/vanilla/uploads/editor/id/mdiwowdmfuj0.png

dp · February 2023

Okay. Donc comme je disais, ça fait un bail que je n’ai pas touché aux méthodes d’interpolations et je suis donc rouillé. De fait, j’ai cherché un peu sur le net et la méthode de Newton semble donner quelque chose comme ça :

fn divided_differences(&self) -> Vec<f64> {
    let n: usize = self.points.len() ;
    let mut diffs: Vec<f64> = vec![self.points[0].y] ;

    for i in 1..n {
        let mut term = self.points[i].y ;

        for j in 0..i {
            term = (term - diffs[j]) / (self.points[i].x - self.points[j].x) ;
        }

        diffs.push(term) ;
    }

    diffs
}

fn newton_interpolation(&self, p: f64) -> f64 {
    let d: Vec<f64> = self.divided_differences();
    let n = self.points.len() ;

    let mut result: f64 = 0.0 ;

    for i in 0..n {
        let mut term = d[i] ;

        for j in 0..i { 
            term *= p - self.points[j].x ;
        }

        result += term ;
    }

    result
}

Ce qui fonctionne enfin ! En effet, si on teste avec ce jeu de données

(1,1)
(2,0)
(4,2)

on récupère le graphique suivant

raoul.S · February 2023

dp a dit :

Ce TP est amusant...

dp a dit :

Si une personne arrive à voir où je me suis foiré, je lui en serais très reconnaissant.

Mais la résolution des bugs fait partie intégrante des joies du développeur cher dp, on ne va pas te priver de cette joie...

marco · February 2023

diffs[i]=numerator/denominator
Est-ce qu'alors diffs[i] ne prend pas, à la fin de la boucle sur $j$, la valeur relative au dernier $j$, c'est-à-dire $j=i-1$, ou alors $j=i$, auquel cas diffs[i]=numerator/denominator=0/0 ?

dp · February 2023

@raoul.S Certes, héhé

. J’ai pu à l’aide du message de @marco changer deux trois petites choses (voir le message précédent que j'ai édité) pour me rapprocher du résultat voulu. Néanmoins, il y a encore quelque chose qui cloche et je soupçonne que ce soit du fait de quelques subtilités que j’aurais mal comprises !

Programmation en C

Réponses

Bonjour!

Catégories

In this Discussion

Qui est en ligne 6