Solving a Linear Diophantine Equation(see descript

问题:

Let me start off by clarifying that(before you guys dismiss me), this is not a homework problem and I'm not a university student. :)

EDIT Thanks to @Klas and others, my question now boils down to a mathematical equation which needs to be solved programmatically.

I'm looking for an algorithm/code which solves Linear Diophantine Equation. For lesser mortals like me, here's how such an equation looks like:

Example 1: 3x + 4y + 5z = 25 (find all possible values of x,y,z)

Example 2: 10p + 5q + 6r + 11s = 224 (find all possible values of p,q,r,s)

Example 3: 8p + 9q + 10r + 11s + 12t = 1012 (find all possible values of p,q,r,s,t)

I tried googling to no avail. I would have thought that some code would already be written to solve this. Do let me know if you guys have come across some kind of library which has already implemented this. And if the solution is in Java, nothing can be cooler!. Algorithm/pseudo code will also do. Thanks much.

回答1:

Brute-force recursion is an option, depending on how large you will allow the value or number of values to become.

Assumptions: The user inputs (the multiplicands) are always distinct positive integers. The coefficients to be found must be non-negative integers.

Algorithm:

Of the multiplicands, let M be the largest.
Calculate C=floor(F/M).
If F=M*C, output solution of the form (0,0,...,C) and decrement C
If M is the only multiplicand, terminate processing
Loop from C down to 0 (call that value i)
  Let F' = F - i*M
  Recursively invoke this algorithm:
    The multiplicands will be the current set minus M
    The goal value will be F'
  For each solution output by the recursive call:
     append i to the coefficient list and output the solution

回答2:

This is a mathematical question rather than a programming one. Once you have a suitable algorithm, impelementing it shouldn't be too hard.

I suggest you google on Diophantine Equations.

I found an explanation for you.

回答3:

I happened to write Java code for this. Please help yourself. The solutions are not extensively tested, but it seems to work well so far.

package expt.qp;

import java.util.HashMap;
import java.util.LinkedHashMap;
import java.util.Map;

public class LinearDiophantine {

    private Map<Integer, Integer> sol = new LinkedHashMap<Integer, Integer>();
    private Map<Integer, Integer> coeff = new HashMap<Integer, Integer>();

    /**
     * @param args
     */
    public static void main(String[] args) {
        // Fill up the data
        // 3x + 4y + 5z + 3a = 25
        LinearDiophantine ld = new LinearDiophantine();
        ld.coeff.put(1, 1);ld.coeff.put(2, 2);ld.coeff.put(3, 3);ld.coeff.put(4, 4);
        Map<Integer, Integer> coeffCopy = new HashMap<Integer, Integer>(ld.coeff);
        int total=30;

        // Real algo begins here
        ld.findPossibleSolutions(total, coeffCopy);

    }

    private void findPossibleSolutions(int total, Map<Integer, Integer> coeff) {
        int index=returnLargestIndex(coeff);
        int range = (int) Math.floor(total/coeff.get(index));
        if(range*coeff.get(index) == total) {
            sol.put(index, range);
            displaySolution();
            //System.out.println();
            range--;
        }
        if(coeff.size() == 1) {
            return;
        }
        while(range>=0) {
            int remTotal = total - range*coeff.get(index);
            Map<Integer, Integer> coeffCopy = new HashMap<Integer, Integer>(coeff);
            coeffCopy.remove(index);
            sol.put(index, range);
            findPossibleSolutions(remTotal, coeffCopy);
            range--;
        }
    }

    private void displaySolution() {
        int total = 0;
        for(int i : sol.keySet()) {
            //System.out.print(coeff.get(i)+"("+sol.get(i)+"), ");
            total = total + (coeff.get(i)*sol.get(i));
        }
        if(total != 30)
            System.out.print(total+",");
    }

    /**
     * @param coeff
     */
    private int returnLargestIndex(Map<Integer, Integer> coeff) {
        int largestKey = coeff.keySet().iterator().next();
        for(int i : coeff.keySet()) {
            if(coeff.get(i)>coeff.get(largestKey)) {
                largestKey=i;
            }
        }
        return largestKey;
    }

}

回答4:

Adding on to Klas' very accurate answer:

Hilbert's 10th problem asked if an algorithm existed for determining whether an arbitrary Diophantine equation has a solution. Such an algorithm does exist for the solution of first-order Diophantine equations. However, the impossibility of obtaining a general solution was proven by Yuri Matiyasevich in 1970

taken from: Wolfram MathWorld

回答5:

A brute force algorithm is as follows (3 variable case):

int sum = 25;
int a1 = 3;
int a2 = 4;
int a3 = 5;
for (int i = 0; i * a1 <= sum; i++) {
    for (int j = 0; i * a1 + j * a2 <= sum; j++) {
        for (int k = 0; i * a1 + j * a2 + k * a3 <= sum; k++) {
            if (i * a1 + j * a2 + k * a3 == sum) {
                System.out.println(i + "," + j + "," + k);
            }
        }
    }
}

To generalize this for the N variable case, you need to convert into a recursive form.

This algorithm is O(f(size, a)^N) for some function f.

• We can place bounds on f as follows: size / MaxValue(a) <= f(size, a) <= size / MinValue(a).
• In the worst case (where all of the a[i]s are 1) f(size, a) is size.

Either way, this is pretty horrendous for large values of N. So while the recursive N variable algorithm would be more elegant, it is probably not very practical.

---

If you are willing to fork out 34 Euro's to Springer Verlag, here's a link to a paper which (according to the abstract) includes an algorithm for solving the N variable case.

回答6:

There are either no, or infinitely many solutions. It is often the case that you have an extra constraint that the solution must match. Is this the case in your problem?

Let's start with the most simple situation where there are two unkowns a*x + b*y = c:

The first step is using the Euclidean algorithm to find the GCD of a and b, let's call itd. As a bonus, the algorithm provides x' and y' such that a*x' + b*y' = d. If d doesn't divide c, then there is no solution. Otherwise, a solution is:

x = x' * (c/d)
y = y' * (c/d)

The second step is to find all solutions. This means we must find all p and q such that a*p + b*q = 0. For if both (x,y) and (X, Y) are solutions, then

a * (X-x) + b * (Y-y) = 0

The answer to this is p = b/d and q = -a/d where d = GCD(a,b) and is already calculated in step 1. The general solution now is:

x = x' * (c/d) + n * (b/d)
y = y' * (c/d) - n * (a/d)

where n is an integer.

The first step is easy to extend to multiple variables. I am not sure about generalizing the second step. My first guess would be to find a solution for all pairs of coefficients and combine these solutions.

回答7:

This is a solution in Perl. rather a hack by using Regex.

Following this blog post to solve algebraic equations using regex.

we can use the following script for 3x + 2y + 5z = 40

#!/usr/bin/perl
$_ = 'o' x 40;
$a = 'o' x 3;
$b = 'o' x 2;
$c = 'o' x 5;
$_ =~ /^((?:$a)+)((?:$b)+)((?:$c)+)$/;
print "x = ", length($1)/length($a), "\n";
print "y = ", length($2)/length($b), "\n";
print "z = ", length($3)/length($c), "\n";

output: x=11, y = 1, z = 1

the famous Oldest plays the piano puzzle ends up as a 3 variable equation

This method applies for a condition that the variables are actually positive and the constant is positive.

回答8:

The reason why there is no library that does this is similar to why you won't find a library to do permutations - you generate so much data that this is probably the wrong thing to do.

More precisely, if you have n variables whose sum totals X, then you will have O(Xn-1) answers. X and n don't have to be very large for this to become an issue.

That said, here is some Python that fairly efficiently figures out all of the necessary information to encode the answer:

def solve_linear_diophantine (*coeff_tuple):
    coeff = list(coeff_tuple)
    constant = coeff.pop()

    cached = []
    for i in range(len(coeff)):
        # Add another cache.
        cached.append({})

    def solve_final (i, remaining_constant):
        if remaining_constant in cached[i]:
            return cached[i][remaining_constant]
        answer = []
        if i+1 == len(coeff):
            if 0 == remaining_constant%coeff[i]:
                answer = [(remaining_constant/coeff[i], [])]
        else:
            for j in range(remaining_constant/coeff[i] + 1):
                finish = solve_final(i+1, remaining_constant - j*coeff[i])
                if finish is not None:
                    answer.append((j, finish))
        if 0 == len(answer):
            cached[i][remaining_constant] = None
            return None
        else:
            cached[i][remaining_constant] = answer
            return answer

    return solve_final(0, constant)

When I say "encode", the data structure looks like this. For each possible coefficient, we'll get an array of (coefficient, [subanswers]). Whenever possible it reuses subanswers to make the data structure as small as possible. If you cant you can recursively expand this back out into an array of answers, and in so doing you'll be pretty efficient. (In fact it is O(nX).) You'll do very little repeating of logic to discover the same facts over and over again. (However the returned list may get very large simply because there is a large list of answers to be returned.)