Calcul des zéros réels et / ou complexes d’un polynôme

rem ============================================================================
rem       Calcul des zéros réels et / ou complexes d’un polynôme
rem                de degré N à coefficients réels
rem                  par la méthode de Bairstow
rem                Auteur : Papydall - Avril 2020
rem ============================================================================
rem Pour des plus amples informations sur la méthode de Bairstow,
rem voici un lien pour un fichier PDF qui explique la méthode et que vous pouvez
rem télécharger.
rem http://www.pcastelan.ups-tlse.fr/EEA_CASTELAN/l2-calcul-sc/3-bairstow.pdf
rem ============================================================================
label saisie, resoudre, quitter, infos, Choix_Point_Fixe
dim a(99),b(99),c(99),p(99) : ' 99 est le degré maximal du polynôme, mais soyons
                            : ' raisonable ! Le programme tel quel permet de
                            : ' résoudre un polynôme jusqu'au degré 20
dim a1,b1,d,e,f,p1,q,r,t,u,v,x,y,z,i,j,k,m,n,epsilon,nb,nbc,tex$

epsilon = 1E-16
nb = 6 : ' nombre de chiffres (par défaut) après la virgule pour l'affichage des solutions
nbc = power(10,nb)

GUI()
end
rem ============================================================================
SUB GUI()
    create_hide
    width 0,800 : height 0,800 : color 0,100,100,50
    shape 5 : shape_rectangle_rounded 5 : top 5,55 : left 5,10 : width 5,760 : height 5,640
      color 5,50,100,100
    shape 6 : shape_rectangle_rounded 6 : top 6,700 : left 6,60 : width 6,660: height 6,50
      color 6,200,150,20  
 '   shape 7 : shape_rectangle_rounded 7 : top 7,05 : left 7,40 : width 7,740: height 7,40
    alpha 8 : top 8,650 : left 8,250 : font_bold 8  : color 8,0,255,255
    font_color 8,255,0,0 : font_size 8,14
    caption 8," Auteur : Papydall - Avril 2020 "
    alpha 10 : top 10,10 : left 10,50  : font_size 10,14
      font_bold 10 : font_color 10,255,0,0 :  color 10,0,255,255  
    tex$ = " Résolution du polynôme : a(0) + a(1)*x + a(2)*x² + ... + a(n)*x^n = 0 "
    caption 0,tex$ : caption 10,tex$ : application_title tex$

    alpha 20 : top 20,85 : left 20,150 : font_name 20,"arial black" : font_size 20,10
      color 20,255,255,255 : font_color 20,0,0,255 : caption 20," Degré du polynôme "
    spin 30 : top 30,80 : left 30,320 : width 30,50 : min 30,1 : max 30,20 : position 30,5
      font_color 30,255,0,0 : font_bold 30 : font_size 30,12 : on_change 30,Saisie
    alpha 35 : top 35,85 : left 35,400 : font_name 35,"arial black" : font_size 35,10
      color 35,255,255,255 : font_color 35,0,0,255 : caption 35," Chiffres après le point décimal "
    spin 36 : top 36,80 : left 36,650 : width 36,50 : min 36,3 : max 36,15 : position 36,6
      font_color 36,255,0,0 : font_bold 36 : font_size 36,12
    container 40 : top 40,140 : left 40,60 : height 40,500 : caption 40,"Coefficients"
      font_color 40,255,255,0 : font_bold 40
    
    for i = 0 to 20
        alpha 100+i : parent 100+i,40 : top 100+i,22*(i+1) : left 100+i,10
        font_bold 100+i : caption 100+i,"a(" + str$(i) + ")"
        edit 200+i : parent 200+i,40 : top 200+i,22*(i+1) : left 200+i,60
        width 200+i,100 : height 200+i,5 : font_color 200+i, 255,0,0
    next i
    
    container 50 : top 50,140 : left 50,250 : height 50,500 : caption 50,"Racines"
      font_color 50,255,255,0 : font_bold 50 : width 50,500
    alpha 300 : parent 300,50 : top 300,20 : left 300,100
      caption 300,"x réel" + string$(60," ") + "x imaginaire"
    list 301 : parent 301,50 : top 301,40 : left 301,10 : width 301,230 : height 301,450
      font_color 301,0,0,255 : font_size 301,10
    list 401 : parent 401,50 : top 401,40 : left 401,250 : width 401,230 :height 401,450
      font_color 401,0,0,255 : font_size 401,10
      
    show_all : set_focus 200
    for i = 0 to 20 : hide 100+i : hide 200+i : next i
    for i = 0 to position(30) : show 100+i : show 200+i : next i
    show 301 : show 401 : create_show
    
    button 60 : top 60,710 : left 60, 100 : width 60,150 : font_bold 60
      caption 60,"&Résoudre le polynôme" : on_click 60,Resoudre
    button 70 : top 70,710 : left 70,600 : font_bold 70
      caption 70,"&Quitter" : on_click 70,quitter
    button 80 : top 80,710 : left 80,380 : font_bold 80
      caption 80,"&Infos" : on_click 80, infos
      
END_SUB
rem ============================================================================
Saisie:
    if position(30) < 1 then position 30,1
    if position(30) > 20 then position 30,20
    for i = 0 to 20 : hide 100+i : hide 200+i : next i
    for i = 0 to position(30) : show 100+i : show 200+i : next i    
return
rem ============================================================================
Resoudre:  
    off_click 60  : inactive 60
    Init()  
    on_click 60,Resoudre
return
rem ============================================================================
Infos:
    tex$ = "Ce programme permet de résoudre un polynôme de degré N (c-à-d calculer toutes ses racines réelles et / ou complexes)" + chr$(13)
    tex$ = tex$ + "Pour utilser ce programme, choisissez d'abord le degré du polynôme à résoudre (le degré 5 est proposé par défaut)." + chr$(13)
    tex$ = tex$ + "Puis saisissez les coefficients a(n) du polynôme en commençant par la constante a(0)." + chr$(13)
    tex$ = tex$ + "Tous les coefficients du polynôme doivent être saisis et le dernier coefficient a(n) ne doit pas être nul !"+chr$(13)
    tex$ = tex$ + "Vous pouvez également choisir le nombres des chiffres après le point décimal pour l'affichage des solutions"+chr$(13)
    tex$ = tex$ + "Ce nombre est compris entre 3 et 15 (valeur 6 par défaut)." + chr$(13)
    tex$ = tex$ + "En demandant la résolution (click sur le bouton 'Résoudre le polynôme', Le programme se déroule sans aucune intervention de votre part,"+chr$(13)
    tex$ = tex$ + "et (si tout va bien) affiche finalement les racines présentées en liste partie réelle, partie imaginaire." + chr$(13)+chr$(13)
    tex$ = tex$ + "L'algorithme utilisé dans ce programme est basé sur la méthode de Bairstow."+chr$(13)
    tex$ = tex$ + "Brièvement, le principe de recherche d'une racine est le suivant :" + chr$(13)
    tex$ = tex$ + "A partir d'un couple (p1,q) défini fixe au départ (0,0) ou (0.1,0.1) ou autre couple, on construit une suite de points (x,y) qui"+chr$(13)
    tex$ = tex$ + "évolue à chaque itération en (x+dp, y+dq) : dp et dq sont les accroissements relatifs à chaque variable." + chr$(13)
    tex$ = tex$ + "Ces accroissements sont calculés à l'aide des coefficients du polynôme sélectionné dans le tableau C()." + chr$(13)
    tex$ = tex$ + "Lorsque les accroissements sont suffisamment petits (en valeur absolue), la suite de points(x,y) converge vers un point donné"+chr$(13)
    tex$ = tex$ + "qui est alors une racine dont x est la partie réelle, y la partie imaginaire."+chr$(13) + chr$(13)
    tex$ = tex$ + "Sachez que le programme effectue un maximum de 500 (constante k) itérations et si, au bout de ces 500 itérations, on n'atteint pas une convergence"+chr$(13)
    tex$ = tex$ + "convénable, on lancera une nouvelle tentative à partir d'un point fixe (p1,q) différent." + chr$(13)
    tex$ = tex$ + "Ce nombre de tentatives est limité à 10."+chr$(13)+chr$(13)
    tex$ = tex$ + "Il peut arriver, dans certains cas, que l'algorithme n'est pas en mesure de résoudre le polynôme (pas de convergence après 10*500 itérations."+chr$(13)
    tex$ = tex$ + "Dans ce cas, le programme rénonce au calcul et le signale par un message adéquat."+chr$(13)+ chr$(13)
    tex$ = tex$ + "Pour des plus amples information sur la méthode de Bairstow, voici un lien : "+chr$(13)
    tex$ = tex$ + "http://www.pcastelan.ups-tlse.fr/EEA_CASTELAN/l2-calcul-sc/3-bairstow.pdf" + chr$(13)
    tex$ = tex$ + "C'est un fichier .PDF qui explique la méthode de Baistow (c'est en français)"
    message tex$
return
rem ============================================================================
SUB Init()
    dim_local tentative
    tentative = 0
    clear 301 : clear 401 : item_add 301,""  : item_add 401,""
    n = position(30)
    m = n
    nb = val(text$(36)) : ' nombre de chiffres après le point décimal pour l'affiche des solutions
    nbc = power(10,nb)
    for i = 0 to n
        if numeric(text$(200+i)) = 0
           Message "!!! ERREUR !!! ... Vérifiez vos coefficients ... Ils doivent être tous numériques !!!"
           active 60 : exit_sub
        end_if
        a(n-i) = val(text$(200+i))
        p(i) = a(i)
    next i
    if a(0) = 0          
       beep_error
       message "!!! ERREUR !!! Le coefficient a(" + str$(n) + ") ne doit pas être nul !!!"
       active 60 : exit_sub
    end_if  
    k = 500    : ' nombre maxi d'itérations avant de passer à un autre choix du point fixe
    e = 1E-10  : ' précision des calculs
    Choix_Point_Fixe:  
    if tentative > 10 : ' Si après 10 tentatives, le système ne converge pas, on arrête
       message "!!! Impossible à calculer !!! Pas de convergence !!!"
       active 60 : exit_sub : ' terminate
    end_if
    p1 = rnd(1) : q = rnd(1) :  ' initialisation du point fixe
    ' boucle principale de factorisation
    while  n > 0    
        if n <= 2 then Dernier_Facteur(): exit_sub
        j = 0
        repeat    
            if j > k
               tentative = tentative + 1 : goto Choix_Point_Fixe : ' une autre tentative
            end_if
            j = j + 1
            ' calcul de b(i) et c(i)
            b(0) = 0 : b(1) = 0 : c(0) = 0 : c(1) = 0
            for i = 2 to n+2
                b(i) = a(i-2) - p1*b(i-1) - q*b(i-2)
                c(i) = 0-b(i-1) - p1*c(i-1) - q*c(i-2)
            next i
            ' calcul de dp = a1 et dq = b1
            x = b(n+1) : y = b(n+2) : z = c(n)
            t = c(n+1) : u = c(n+2)
            d = t*t - z*(u+x)
            if d = 0
               tentative = tentative + 1 : goto Choix_Point_Fixe : ' une autre tentative
            end_if
            a1 = (z*y - x*t) / d
            b1 = (0-x*(q*z + p1*t) - y*t) / d
            ' nouveaux p1 et q
            p1 = p1 + a1 : q = q + b1
            f = (abs(a1) + abs(b1)) / (abs(p1) + abs(q))
        until f <= e
        ' un facteur est trouvé
        Trinome()
        ' mise à jour  
        n = n-2
        for i = 0 to n : a(i) = b(i+2) : next i
    end_while
    Dernier_Facteur()  
END_SUB
rem ============================================================================
' 1er ou 2ème degré
SUB Dernier_Facteur()
    if n = 2  : ' deuxième degré
       p1 = a(1) / a(0) : q = a(2) / a(0)
       Trinome() : active 60 : exit_sub
    end_if
    x = 0-a(1) / a(0) : y = 0 : ' 1er degré
    Afficher_Solutions() : active 60
  
END_SUB
rem ============================================================================
' résoudre x² + p1.x + q = 0
SUB Trinome()
    dim_local d
    d = p1*p1 - 4*q
    if d >= 0
       d = sqr(d) : y = 0
       x = (0-p1+d)/2 : Afficher_Solutions()
       x = (0-p1-d)/2 : Afficher_Solutions()
    else
       d = sqr(0-d)/2 : x = 0-p1/2
       y = d   : Afficher_Solutions()
       y = 0-d : Afficher_Solutions()
    end_if
    
END_SUB
rem ============================================================================
SUB Afficher_Solutions()  
    dim_local a1,reel$,imag$  
    x = x + sgn(x)*epsilon : x = int(x*nbc+0.5)/nbc
    y = y + sgn(y)*epsilon : y = int(y*nbc+0.5)/nbc  
    if abs(x) < 1E-4 then x = 0
    if abs(y) < 1E-4 then y = 0
    if x >= 0 then reel$ = " " + str$(x) : else : reel$ = str$(x)  
    if y >= 0 then imag$ = " " + str$(y) : else : imag$ = str$(y)
    item_add 301,reel$
    if y <> 0
       item_add 401,imag$
    else
       item_add 401,"---"
    end_if
    
END_SUB
rem ============================================================================
Quitter:
    if message_information_yes_no("Vous voulez vraiment quitter ?") = 1 then terminate
return
rem ============================================================================