Линейное диофантово уравнение с двумя переменными - это уравнение вида ax + by = c , где a , b и c - постоянные целые числа, а x и y - целочисленные переменные.

Для многих естественных диофантовых уравнений x и y представляют величины, которые не могут быть отрицательными.

задача

Напишите программу или функцию, которая принимает коэффициенты a , b и c в качестве входных данных и возвращает произвольную пару натуральных чисел (0, 1, 2,…) x и y, которые проверяют уравнение ax + by = c , если такая пара существуют.

Дополнительные правила

• Вы можете выбрать любой формат для ввода и вывода, который включает только нужные целые числа и, необязательно, массив / список / матрицу / кортеж / векторную нотацию вашего языка, если вы не вставляете какой-либо код во входные данные.

• Вы можете предположить, что коэффициенты a и b не равны нулю.

• Ваш код должен работать для любой тройки целых чисел от -2 ⁶⁰ до 2 ⁶⁰ ; он должен завершиться менее чем за минуту на моей машине (Intel i7-3770, 16 ГБ ОЗУ).

• Вы не можете использовать любые встроенные модули, которые решают диофантовы уравнения и, таким образом, упрощают эту задачу, такие как Mathematica FindInstanceили FrobeniusSolve.

• Ваш код может вести себя так, как вы хотите, если решение не может быть найдено, если оно соответствует ограничению по времени и его вывод не может быть перепутан с допустимым решением.

• Применяются стандартные правила игры в гольф .

Примеры

1. Приведенные ниже примеры иллюстрируют действительный ввод / вывод для уравнения 2x + 3y = 11 , которое имеет ровно два действительных решения ( (x, y) = (4,1) и (x, y) = (1,3) ).

   Input:  2 3 11
   Output: [4 1]
   
   Input:  (11 (2,3))
   Output: [3],(1)
   

2. Единственным верным решением 2x + 3y = 2 является пара (x, y) = (1,0) .

3. Приведенные ниже примеры иллюстрируют действительный ввод / вывод для уравнения 2x + 3y = 1 , которое не имеет действительных решений .

   Input:  (2 3 1)
   Output: []
   
   Input:  1 2 3
   Output: -1
   
   Input:  [[2], [3], [1]]
   Output: (2, -1)
   

4. Для (a, b, c) = (1152921504606846883, -576460752303423433, 1) все правильные решения (x, y) удовлетворяют тому, что (x, y) = (135637824071393749 - bn, 271275648142787502 + an) для некоторого неотрицательного целого числа n ,

Pyth, 92 байта

I!%vzhK%2u?sm,ed-hd*ed/F<G2cG2@G1G+~Q,hQ_eQj9 2)J*L/vzhKtKeoSNm-VJ/RhK_*LdQsm+LdtM3/V*LhK_JQ

Это довольно монстр.

Попробуйте онлайн: Демонстрация . Формат ввода c\n[a,b]и выходной формат [x,y].

В случае, если целочисленного решения не существует, я ничего не буду печатать, а в случае, если естественного целочисленного решения не существует, я просто выведу случайное целочисленное решение.

Пояснение (приблизительный обзор)

1. Сначала я найду целочисленное решение уравнения ax + by = gcd(a,b)с помощью расширенного евклидова алгоритма.

2. Затем я изменю решение (мое умножение aи bс c/gcd(a,b)), чтобы получить целочисленное решение ax + by = c. Это работает, если c/gcd(a,b)является целым числом. В противном случае не существует решения.

3. Все остальные целочисленные решения имеют вид a(x+n*b/d) + b(y-n*a/d) = c с d = gcd(a,b)целым числом n. Используя два неравенства x+n*b/d >= 0и y-n*a/d >= 0я могу определить , 6 возможных значений n. Я попробую все 6 из них и распечатаю решение с самым высоким наименьшим коэффициентом.

Пояснение (подробно)

Первый шаг - найти целочисленное решение уравнения ax' + by' = gcd(a,b). Это можно сделать с помощью расширенного евклидова алгоритма. Вы можете получить представление о том, как это работает в Википедии . Разница лишь в том, что вместо использования 3 columns ( r_i s_i t_i) я буду использовать 6 columns ( r_i-1 r_i s_i-1 s_i t_i-1 t_i). Таким образом, мне не нужно хранить последние две строки в памяти, только последнюю.

K%2u?sm,ed-hd*ed/F<G2cG2@G1G+~Q,hQ_eQj9 2)   implicit: Q = [a,b] (from input)
                                     j9 2    convert 9 to base 2: [1,0,0,1]
                            + Q              add to Q => [a,b,1,0,0,1]
                                             this is the initial row
   u                                     )   start with G = ^ and update G repeatedly
                                             by the following expression, until
                                             the value of G doesn't change anymore
    ?                   @G1                    if G[1] != 0:
                     cG2                         split G into parts of 2
      m                                          map the parts d to:
       ,                                           the pair 
        ed                                           d[1]
          -hd*ed/F<G2                                d[0]-d[1]*G[0]/G[1]
     s                                           unfold
                                               else:
                           G                     G (don't change it, stop criterion for u)
 %2                                          take every second element
                                             we get the list [gcd(a,b),x',y']
K                                            store this list in K
                             ~Q,hQ_eQ        afterwards change Q to [Q[0],-Q[1]] = [a,-b]
                                             This will be important for the other parts. 

Теперь я хочу найти решение ax + by = c. Это возможно только тогда, когда c mod gcd(a,b) == 0. Если это уравнение выполнено, я просто умножаюсь x',y'на c/gcd(a,b).

I!%vzhK...J*L/vzhKtK   implicit: z = c in string format (from input)
  %vzhK                evaluated(z) mod K[0] (=gcd(a,b))
I!                     if not ^ than: 
             /vzhK        c/K[0]
           *L     tK      multipy ^ to each element in K[1:] (=[x',y'])
          J               and store the result in J, this is now [x,y]

У нас есть целочисленное решение для ax + by = c. Обратите внимание, что x, yили оба могут быть отрицательными. Поэтому наша цель - преобразовать их в неотрицательные.

Хорошая вещь о диофантовых уравнениях состоит в том, что мы можем описать все решения, используя только одно начальное решение. Если (x,y)есть решение, что все другие решения имеют вид (x-n*b/gcd(a,b),y+n*a/gcd(a,b))на nцелое число.

Поэтому мы хотим найти n, где x-n*b/gcd(a,b) >= 0и y+n*a/gcd(a,b >= 0. После некоторого преобразования мы получаем два неравенства n >= -x*gcd(a,b)/bи n >= y*gcd(a,b)/a. Обратите внимание, что символ неравенства может выглядеть в другом направлении из-за деления с потенциальным отрицательным aили b. Меня это не сильно волнует, я просто говорю, что одно число -x*gcd(a,b)/b - 1, -x*gcd(a,b)/b, -x*gcd(a,b)/b + 1определенно удовлетворяет неравенству 1, а одно число y*gcd(a,b)/a - 1, y*gcd(a,b)/a, y*gcd(a,b)/a + 1удовлетворяет неравенству 2. Если есть a n, которое удовлетворяет обоим неравенствам, одно из 6 чисел также делает.

Затем я рассчитываю новые решения (x-n*b/gcd(a,b),y+n*a/gcd(a,b))для всех 6 возможных значений n. И я печатаю решение с самым высоким наименьшим значением.

eoSNm-VJ/RhK_*LdQsm+LdtM3/V*LhK_JQ
                               _J    reverse J => [y,x]
                           *LhK      multiply each value with K[0] => [y*gcd,x*gcd]
                         /V      Q   vectorized division => [y*gcd/a,-x*gcd/b]
                  m                  map each d of ^ to:
                      tM3              [-1,0,1]
                   +Ld                 add d to each ^
                 s                   unfold
                                     these are the possible values for n
    m                                map each d (actually n) of ^ to:
             *LdQ                      multiply d to Q => [a*n,-b*n]
            _                          reverse => [-b*n,a*n]
        /RhK                           divide by K[0] => [-b*n/gcd,a*n/gcd]
     -VJ                               vectorized subtraction with J
                                       => [x+b*n/gcd,y-a*n/gcd]
 oSN                                 order the solutions by their sorted order
e                                    print the last one

Сортировка по порядку сортировки работает следующим образом. Я использую пример2x + 3y = 11

Я сортирую каждое из 6 решений (это называется ключами) и сортирую исходные решения по ключам:

solutions: [1, 3], [4, 1], [7, -1], [-5, 7], [-2, 5], [1, 3]
keys:      [1, 3], [1, 4], [-1, 7], [-5, 7], [-2, 5], [1, 3]
sort by key:
solutions: [-5, 7], [-2, 5], [7, -1], [1, 3], [1, 3], [4, 1]
keys:      [-5, 7], [-2, 5], [-1, 7], [1, 3], [1, 3], [1, 4]

Это сортирует полное неотрицательное решение до конца (если есть).

• после замечаний Денниса, из-за которых моя предыдущая идея перевернулась с ног на голову, мне пришлось сменить код с его корней, и это заняло у меня длительную отладку и стоило мне вдвое больше n ° байтов: '(.

Матлаб (660)

a=input('');b=input('');c=input('');if((min(a*c,b*c)>c*c)&&a*c>0&&b*c>0)||(a*c<0&&b*c<0),-1,return,end,g=abs(gcd(a,b));c=c/g;a=a/g;b=b/g;if(c~=floor(c)),-1,return,end,if(c/a==floor(c/a)&&c/a>0),e=c/a-b;if(e>0),e,a,return,else,c/a,0,return,end,end,if(c/b==floor(c/b)&&c/b>0),e=c/b-a;if(e>0),b,e,return,else,0,c/b,return,end,end,f=max(abs(a),abs(b));if f==abs(a),f=b;b=a;a=f;g=0.5;end,e=(c-b)/a;f=(c-2*b)/a;if(e<0&&f<e),-1,elseif(e<0&&f>e),for(i=abs(c*a):abs((c+1)*a)),e=(c-i*b);if(mod(e,a)==0)if(g==0.5),i,e/a;else,e/a,i,end,return,end,end,else for(i=1:abs(a)),e=(c-i*b);if(e/a<0),-1,elseif(mod(e,a)==0),if(g==0.5),i,e/a,else,e/a,i,end,return,end,end,end,-1

• Ну, я знаю, что это не игра в гольф, так как этот тип языков не приспособлен для уменьшения длины кода, но я могу гарантировать, что сложность времени в лучшем виде.

Объяснение:

• код принимает в качестве входных данных три инварианта a, b, c, эти последние подчиняются паре условий перед тем, как продолжить вычисление:

  1 - если (a + b> c) и (a, b, c> 0) нет решения!

  2- если (a + b <c), (a, b, c <0) нет решения!

  3 - если (a, b) имеют общие противоположные признаки c: нет решения!

  4 - если GCD (a, b) не делится на c, то снова нет решения! иначе делим все варианты по GCD.

• после этого мы должны проверить другое условие, оно должно облегчить и сократить путь к желаемому решению.

  5- если c разделить a или b, решение s = (x или y) = (c- [ax, yb]) / [b, a] = C / [b, a] + [ax, yb] / [b , a] = S + [ax, yb] / [b, a], где S является естественным, поэтому у ax / b или / a будут впредь неотрицательные прямые решения, которые соответственно x = b или y = a. (обратите внимание, что решения могут быть просто нулевыми значениями в случае, если предыдущие произвольные решения выявлены как отрицательные)

• когда программа достигает этой стадии, более узкий диапазон решений для x = (c-yb) / a вместо этого, благодаря конгруэнтности, охватывается большими диапазонами чисел, которые периодически повторяются регулярными циклами. самое большое поле поиска - [xa, x + a], где a - делитель.

ПОПЫТАЙСЯ

Аксиома, 460 байт

w(a,b,x,u)==(a=0=>[b,x];w(b rem a,a,u,x-u*(b quo a)))
d(a,b,k)==(o:List List INT:=[];a=0 and b=0=>(k=0=>[1,1];[]);a=0=>(k=0=>[[1,0]];k rem b=0=>[1,k quo b];[]);b=0=>(k=0=>[[0,1]];k rem a=0=>[k quo a,1];[]);r:=w(a,b,0,1);q:=k quo r.1;(y,x,u,v):=(q*(r.1-r.2*a)quo b,q*r.2,b quo r.1,a quo r.1);m:=min(80,4+abs(k)quo min(abs(a),abs(b)));l:=y quo v;x:=x+l*u;y:=y-l*v;for n in -m..m repeat(t:=x+n*u;z:=y-n*v;t>=0 and z>=0 and t*a+z*b=k=>(o:=cons([t,z],o)));sort(o))

негольф и какой-то тест

-- input a b and k for equation a*x+b*y=k
-- result one List of List of elments [x,y] of solution of  
-- that equation with x and y NNI (not negative integers) 
-- or Void list [] for no solution
diopanto(a,b,k)==
  o:List List INT:=[]
  a=0 and b=0=>(k=0=>[1,1];[])
  a=0=>(k=0=>[[1,0]];k rem b=0=>[1,k quo b];[])
  b=0=>(k=0=>[[0,1]];k rem a=0=>[k quo a,1];[])
  r:=w(a,b,0,1)
  q:=k quo r.1
  (y,x,u,v):=(q*(r.1-r.2*a)quo b,q*r.2,b quo r.1,a quo r.1)
  m:=min(80,4+abs(k)quo min(abs(a),abs(b)))
  l:=y quo v           -- center the interval
  x:=x+l*u; y:=y-l*v
  for n in -m..m repeat
     t:=x+n*u;z:=y-n*v
     t>=0 and z>=0 and t*a+z*b=k=>(o:=cons([t,z],o))
  sort(o)

 ------------------------------------------------------
(4) -> d(0,-9,0)
   (4)  [[1,0]]
                                                  Type: List List Integer
(5) -> d(2,3,11)
   (5)  [[4,1],[1,3]]
                                                  Type: List List Integer
(6) -> d(2,3,2)
   (6)  [[1,0]]
                                                  Type: List List Integer
(7) -> d(2,3,1)
   (7)  []
                                                  Type: List List Integer
(8) -> d(1152921504606846883,-576460752303423433,1)
   (8)
   [[135637824071393749,271275648142787502],
    [712098576374817182,1424197152749634385],
    [1288559328678240615,2577118657356481268],
    [1865020080981664048,3730040161963328151],
    [2441480833285087481,4882961666570175034]]
                                                  Type: List List Integer

В других возможных «решениях» была ошибка, потому что она пыталась сохранить бесконечные решения в одном Списке; теперь наложено ограничение в 80 решений макс