Équation diophantienne ax + by = c

rem ============================================================================
rem              Équation diophantienne ax + by = c
rem                 Papydall le 01 / 11 / 2019
rem ============================================================================
rem L’équation ax + by = c, où les coefficients a, b et c sont trois entiers relatifs
rem (a et b non tous deux nuls) et où les inconnues x et y sont des entiers relatifs,
rem est une des équations diophantiennes les plus simples à résoudre.
rem Sa résolution s’appuie sur l’algorithme d’Euclide, le théorème de Bachet-Bézout
rem (qui correspond au cas, appelé aussi identité de Bézout, où c est égal au PGCD
rem de a et b) et le théorème de Gauss.
rem Dans l’ensemble des entiers relatifs, une telle équation possède, ou bien aucune
rem solution, ou bien une infinité de solutions.
rem ============================================================================
dim true, false
dim a,b,c,x0,y0,p,q,t$,plus$
true = 1 : false = 0
caption 0,"Résolution dans Z de l'équation diophantienne ax + by = c"
memo 10 : full_space 10 : color 10,200,150,50 : font_name 10,"arial" : font_bold 10
font_size 10,12 : font_color 10,0,0,255
item_add 10,""
item_add 10,string$(10," ") + "Résolution dans Z de l'équation diophantienne ax + by = c"
rem ============================================================================
' Exemples d'équations diophantiennes
 a = 3   : b = -2  : c = 7   : ' <---- X = 3 + 2 * k    ; Y = 1 + 3 * k
' a = 3   : b = 7   : c = 4   : ' <---- X = -1 + 7 * k   ; Y = 1 - 3 * k
' a = 63  : b = 85  : c = 1   : ' <---- X = 27 + 85 * k  ; Y = -20 - 63 * k
' a = 65  : b = 104 : c = 26  : ' <---- X = 2 + 8 * k    ; Y = -1 - 5 * k
' a = 217 : b = 34  : c = 2   : ' <---- X = -26 + 34 * k ; Y = 166 - 217 * k
' a = 544 : b = 944 : c = 160 : ' <---- X = 35 + 59 * k  ; Y = -20 - 34 * k
' a = 0   : b = 1   : c = 3   : ' <-----------------------  Pas de solutions !!!
rem ============================================================================
item_add 10,""
if b > 0 then plus$ = "+ " : else : plus$ = ""
t$ = " Pour l'équation : " + str$(a) + "x " + plus$ + str$(b) + "y = " + str$(c)
item_add 10, t$
if Diophancian(a,b,c) = true  
   item_add 10," Les solutions sont :" : item_add 10,"================="
   item_add 10,"   X = " + str$(x0) + " + " + str$(abs(q)) + " * k"
   t$ = "   Y = " + str$(y0)  
   if p*q > 0 then t$ = t$ + " -" : else : t$ = t$ + " + "
   item_add 10,t$ + str$(abs(p)) + " * k"
   item_add 10, "avec k dans l'ensemble Z des entiers rélatifs"    
else
   item_add 10,"Désolé ! Pas de solutions"
end_if
end
rem ============================================================================
' Retourne le Plus Grand Commun Diviseur de deux nombres entiers
FNC PGCD(a,b)
    dim_local r,x
    a = int(abs(a)) : b = int(abs(b))
    if (a > 1E10) or (b > 1E10) then result 1 : exit_fnc
    if (a = 0) or (b = 0) then result 1 : exit_fnc
    if a < b then x = a : a = b : b = x
    repeat
        r = a - b * int(a/b) : a = b : b = r
    until abs(r) < 1E-10
    result a  
END_FNC
rem ============================================================================
' Résolution de l'équation ax + by = c
' Paramètres d'entrée : a,b,c coefficients de l'équation (nombres entiers)
' Retourne 1 (true)
'        si des solutions existent (si le PGCD de a et b est aussi diviseur de c
' Retourne 0 (false) dans le cas contraire
' Cette fonction modifient les variables globales x0, y0, p et q qui contiendront
' les solutions X et Y si elles existent sous forme :
' X = x0 + k*p et Y = y0 - k*q avec k = 0,1,2,3, ... ou k = -1,-2,-3 ...
FNC Diophancian(a,b,c)
    dim_local pg,x1,x2,y1,y2,trouve,ret
    ret = false
    if (a = 0) or (b = 0) then result ret : exit_fnc
    pg = PGCD(a,b)
    a = a/pg : b = b/pg : c = c/pg
    if c <> int(c) then result ret : exit_fnc : ' pg doit être aussi diviseur de c
    x1 = 0 : y2 = 0 : trouve = false
    repeat
        y1 = (c-a*x1) / b
        if y1 = int(y1)
           x0 = x1 : y0 = y1 : trouve = true
        else
           x1 = 0-x1 : if x1 >= 0 then x1 = x1+1
           x2 = (c-b*y2) / a
           if x2 = int(x2)
              x0 = x2 : y0 = y2 : trouve = true
           else
              y2 = 0-y2 : if y2 >= 0 then y2 = y2+1
           end_if
        end_if
    until trouve = true
    p = a : q = b
    result true    
END_FNC
rem ============================================================================