![]() |
|
![]() |
|
|
Thread Tools | Display Modes |
|
|
#1 |
|
Hobbyist Programmer
|
Solving Diophantine Equations
I have a project where I am to write a programme that illustrates how computers can be used for solving iterative tasks, using Euclid's algorithm as an example, and presenting and explaining it to a class in which most students have little to no knowledge about programming.
Thus, I need to keep it as simple as possible, and I chose to use a high-level language like Python. Euclid's algorithm went smoothly, but I was additionally requested to add a function for solving linear diophantine equations with two variables and one equation, and I'm not even sure how to begin tackling this. I'm not asking for a solution at all, since this is my task; I only need a way to get started with this. For example, given the equation 75x + 42y = 3, the greatest common divisor is 3. 33 = 75 - 1 * 42 9 = 42 - 1 * 33 6 = 33 - 3 * 9 3 = 9 - 1 * 6 Having stored all the numbers in the process (in the form of [[a,b,c,d],[a,b,c,d],[a,b,c,d],...]), how can I work my way 'back up'? Do I have to treat the equation as a string and keep replacing on my way up, or is there a simpler way? Here's the function I wrote to find gdc(a,b) (in Norwegian): # Finn største felles divisor av tallene a og b.
def Finn_sfd(a, b):
global stop, eualg
if a % b == 0: # Sjekk om 'b' går opp i 'a'
print "%d modulo %d = 0\n" % (a, b)
if stop == True:
raw_input("Trykk Enter...\n")
return b # Hvis ja, ferdig; send tilbake resultat
k = a / b # Finn k (antall ganger 'b' går opp i 'a')
r = a - b * k # Finn rest
# a = k * b + r
print "\t%d = %d * %d + %d" % (a, k, b, r)
#a - k * b = r
print "\t%d - %d * %d = %d" % (a, k, b, r)
# sfd(a, b) = sfd(b, r)
print "sfd(%d, %d) = sfd(%d, %d)\n" % (a, b, b, r)
if stop == True:
raw_input("Trykk Enter...\n")
eualg.append([r, a, k, b]) # Lagre tallene for senere bruk
return Finn_sfd(b, r)Sample run of the gcd function (also Norwegian): Finner sfd(108, 63)...
Trykk Enter...
108 = 1 * 63 + 45
108 - 1 * 63 = 45
sfd(108, 63) = sfd(63, 45)
63 = 1 * 45 + 18
63 - 1 * 45 = 18
sfd(63, 45) = sfd(45, 18)
45 = 2 * 18 + 9
45 - 2 * 18 = 9
sfd(45, 18) = sfd(18, 9)
18 modulo 9 = 0 |
|
|
|
|
|
#2 | |
|
Programming Guru
![]() Join Date: Aug 2005
Location: England
Posts: 1,499
Rep Power: 4
![]() |
Quote:
|
|
|
|
|
|
|
#3 |
|
Hobbyist Programmer
|
Thanks.
I had the extended Euclidean algorithm in mind for solving this, though I still have to find the two variables before that. I'll check the article and links from there, though, thanks. |
|
|
|
|
|
#4 |
|
Programming Guru
![]() Join Date: Aug 2005
Location: England
Posts: 1,499
Rep Power: 4
![]() |
Can't you just ask the user for three numbers, e.g.
python Syntax (Toggle Plain Text)
|
|
|
|
|
|
#5 |
|
Hobbyist Programmer
|
That's what the program does already, but it's ax + by = c, and I have to find (i.e. calculate) a possible solution for x and y, then generalise it and find all possible integral solutions.
|
|
|
|
|
|
#6 |
|
Hobbyist Programmer
|
I took a quick stab at it. I'm not done, but heres a strategy I was working on.
First, find the maximum value for x (by minimizing y) and then a max value for y (by minimizing x) to find the max value of x, write a loop that breaks when a solution is found, set y and x to 1, and increase x until the result is greater then c, then increment y and repeat (hope thats understandable) To find all solutions, then, you could just bruteforce it and try every combination of numbers less then max x and max y. Theres certainly cleverer ways to do it, thats just as far as i've gotten so far and I have to go to class now ![]() good luck, hope that hint helps. |
|
|
|
|
|
#7 | |
|
Programming Guru
![]() Join Date: Aug 2005
Location: England
Posts: 1,499
Rep Power: 4
![]() |
Quote:
![]() Let's go over what we have so far. We have the Euclidean function gcd, and the Extended Euclidean function ex_gcd. Given the equation 75x + 42y = 6, we can deduce the following: python Syntax (Toggle Plain Text)
And the output we want is a complete range of solutions like: x = -10 + 14n y = 18 - 25n You basically have all the information you need. You just need to use it to generate the right set of output ranges. |
|
|
|
|
|
|
#8 |
|
Hobbyist Programmer
|
I've mostly figured it out. Thanks for the help!
![]() The only things bothering me now are printing every step described in painstaking detail, and finding solutions for positive x and y--I'm on the right track, but it keeps messing up for some numbers, and working for others! I'll figure it out, though... |
|
|
|
|
|
#9 |
|
Hobbyist Programmer
|
Okay, it's pretty much done now. I've posted a sample run at the bottom of this post (translated to English).
I might post the code later, if anyone's interested, but the code needs a lot more documentation and spiffing up (or down, in some cases) before being presentable. Also, I realise that my usage of the implication, equivalence and inequality signs near the end is a bit clumsy and hard to read. If you have a better suggestion, do post it. Input three integers a, b, c in the equation ax + bx = c
a = 108
b = 63
c = 765
------------------------------------------------------------
Finding gcd(108, 63)...
108 = 1 * 63 + 45
108 - 1 * 63 = 45
gcd(108, 63) = gcd(63, 45)
63 = 1 * 45 + 18
63 - 1 * 45 = 18
gcd(63, 45) = gcd(45, 18)
45 = 2 * 18 + 9
45 - 2 * 18 = 9
gcd(45, 18) = gcd(18, 9)
18 modulo 9 = 0
gcd(108, 63) = 9
108 / 9 = 12 + 0
63 / 9 = 7 + 0
------------------------------------------------------------
Finding lcm(108, 63)...
a * b = gcd(a, b) * lcm(a, b)
lcm(a, b) = a * b / gcd(a, b)
lcm = 108 * 63 / 9
lcm = 756
lcm(108, 63) = 756
------------------------------------------------------------
Checking if 108x + 63y = 765 has a solution...
765 modulo 9 = 0
Has solution! Solving 108x + 63y = 765...
108 = 1 * 63 + 45
63 = 1 * 45 + 18
45 = 2 * 18 + 9
18 = 2 * 9 + 0
Returning 1 and -2...
Returning -2 and 3...
Returning 3 and -5...
9 = (3)*108 + (-5)*63
Finding a general solution for 108x + 63y = 765...
x[0] = 3 * (765 / 9) = 3 * 85 = 255
y[0] = -5 * (765 / 9) = -5 * 85 = -425
x = x[0] + b / d * z
y = y[0] - a / d * z
x = 255 + 63 / 9 * z = 255 + 7z
y = -425 - 108 / 9 * z = -425 - 12z
108(255 + 7z) + 63(-425 - 12z) = 765
Finding positive integral solutions...
... for x...
x >= 0 =>
255 + 7z >= 0 <=>
7z >= -255 =>
z >= -36.429
z >= -36
... for y...
y >= 0 =>
-425 - 12z >= 0 <=>
-12z <= 425 <=>
12z <= -425 =>
z <= -35.417
z <= -35
There are 2 positive integral solutions for z between -35 and -36.
z = -36
108(255 + 7 * (-36)) + 63(-425 - 12 * (-36))
= 108 * (3) + 63 * (7)
= 765
z = -35
108(255 + 7 * (-35)) + 63(-425 - 12 * (-35))
= 108 * (10) + 63 * (-5)
= 765
------------------------------------------------------------
Done! |
|
|
|
![]() |
| Bookmarks |
| Currently Active Users Viewing This Thread: 1 (0 members and 1 guests) | |
| Thread Tools | |
| Display Modes | |
|
|
Similar Threads
|
||||
| Thread | Thread Starter | Forum | Replies | Last Post |
| Win PCs, xboxes, zunes and more by solving puzzles! | bulio | Coder's Corner Lounge | 0 | Jan 9th, 2007 8:10 PM |
| solving warning messages, in a program | kalulu | C | 4 | Nov 25th, 2005 8:41 PM |