4

我正在尝试为给定的 H 生成以下方程的所有解。

随着 H=4 ：

1) ALL solutions for x_1 + x_2 + x_3 + x_4 =4
2) ALL solutions for x_1 + x_2 + x_3 = 4
3) ALL solutions for x_1 + x_2 = 4
4) ALL solutions for x_1 =4

对于我的问题，总是有 4 个方程需要求解（独立于其他方程）。总共有 2^(H-1) 个解。对于上一个，以下是解决方案：

1) 1 1 1 1
2) 1 1 2 and 1 2 1 and 2 1 1
3) 1 3 and 3 1 and 2 2
4) 4

这是解决问题的R算法。

library(gtools)
H<-4
solutions<-NULL

for(i in seq(H))
{
    res<-permutations(H-i+1,i,repeats.allowed=T)
    resum<-apply(res,1,sum)
    id<-which(resum==H)

    print(paste("solutions with ",i," variables",sep=""))
    print(res[id,])
}

但是，该算法进行的计算超出了需要。我相信有可能走得更快。我的意思是不生成总和为 > H 的排列

对于给定的 H 有更好的算法的想法吗？

4

3 回答 3

5

与许多问题一样，当某些术语已知时，解决方案变得更容易找到/研究。

这些问题的解决方案被称为整数组合，它是整数分区的泛化（其中顺序无关紧要，即只考虑在排列下唯一的答案）。

例如，4的整数分区为：1+1+1+1、1+1+2、1+3、2+2、4，而4的整数组成为：1+1+1+1， 1+1+2、1+2+1、2+1+1、1+3、3+1、2+2、4。

有一些现成的实现（对与语言无关的算法的引用如下）：

  • 由于您在R中工作，因此该partitions包可以为您生成分区。您需要找到每个分区的唯一排列才能获得组合（请参阅这个 SO question）。
  • 如果您能够使用另一种语言（通过与 R 接口，或通过预先计算答案），那么Mathematica具有计算合成的功能：Compositions.
  • Sage是免费的（与 Mathematica 不同），并且还具有生成内置合成的功能：Compositions. 值得注意的是，这是使用生成器实现的，它可以更有效地使用内存。
  • Python：请参阅此Stack Overflow 问题（生成分区，然后您可以对其进行置换）。我在这里做了类似的事情（尽管它用于itertools查找排列，然后我需要过滤独特的排列，因此可以通过使用专门用于multisets的排列算法更有效地完成此操作）。

为了更好地理解算法（或自己实现），您可以查看这本不完整但有用的电子书：Frank Ruskey的组合生成，它展示了如何在恒定摊销时间(CAT) 中生成分区。由于您想要组合，您还可以使用 CAT 算法生成排列（也在书中）来生成每个整数分区的排列。

Ruskey还解释了如何对它们进行排名和取消排名，这对于存储/散列结果非常方便。

我相信 Knuth 的The Art of Computer Programming 第 4A 卷也很好地介绍了这些内容，如果你碰巧有的话。

ElKamina 建议递归解决它是一个很好的建议，但我不会对大 H 使用这种方法；由于R（以及 Python）没有优化 tail-calls，您最终可能会出现堆栈溢出。

于 2013-02-15T04:32:06.877 回答
2

这是C++ 中的一个实现

废话.cpp：

#include <stdlib.h>
#include <iostream>
#include <vector>

using namespace std;

vector<int> ilist;

void diophantine(int n)
{
    size_t i;
    if (n==0) 
    {
        for (i=0; i < ilist.size(); i++) cout << " " << ilist[i];
        cout << endl;
    }
    else
    {
        for (i=n; i > 0; i--)
        {
            ilist.push_back(i);
            diophantine(n-i);
            ilist.pop_back();
        }
    }          
}


int main(int argc, char** argv)
{
    int n;    

    if (argc == 2 && (n=strtol(argv[1], NULL, 10)))
    {
        diophantine(n);
    }
    else cout << "usage: " << argv[0] << " <Z+>" << endl;

    return 0;
}


命令行的东西：

$ g++ -oblah blah.cpp
$ ./blah 4
 4
 3 1
 2 2
 2 1 1
 1 3
 1 2 1
 1 1 2
 1 1 1 1
$


这是一个实现bash：

废话.sh：

#!/bin/bash

diophantine()
{
    local i
    local n=$1
    [[ ${n} -eq 0 ]] && echo "${ilist[@]}" ||
    {
        for ((i = n; i > 0; i--))
        do
            ilist[${#ilist[@]}]=${i}
            diophantine $((n-i))
            unset ilist[${#ilist[@]}-1]
        done               
    }    
}

RE_POS_INTEGER="^[1-9]+$"
[[ $# -ne 1 || ! $1 =~ $RE_POS_INTEGER ]] && echo "usage: $(basename $0) <Z+>" ||
{
    declare -a ilist=
    diophantine $1
}
exit 0


这是Python中的一个实现

废话.py：

#!/usr/bin/python

import time
import sys


def output(l):
    if isinstance(l,tuple): map(output,l) 
    else: print l,


#more boring faster way -----------------------
def diophantine_f(ilist,n):
    if n == 0:
        output(ilist)
        print
    else: 
        for i in xrange(n,0,-1):
            diophantine_f((ilist,i), n-i)


#crazy fully recursive way --------------------
def diophantine(ilist,n,i):
    if n == 0:
        output(ilist)
        print
    elif i > 0:
        diophantine(ilist, n, diophantine((ilist,i), n-i, n-i))
    return 0 if len(ilist) == 0 else ilist[-1]-1 


##########################
#main
##########################
try:

    if    len(sys.argv) == 1:  x=int(raw_input())
    elif  len(sys.argv) == 2:  x=int(sys.argv[1])
    else: raise ValueError 

    if x < 1: raise ValueError

    print "\n"
    #diophantine((),x,x)
    diophantine_f((),x)    
    print "\nelapsed: ", time.clock()

except ValueError:
    print "usage: ", sys.argv[0], " <Z+>"
    exit(1)
于 2012-06-06T01:42:13.443 回答
0

我假设您没有尝试同时求解方程。

您可以使用递归或动态编程来解决这个问题。

如果您使用递归，只需为第一个变量分配一个有效值并递归求解其余部分。

这里n是变量的数量，sum是总和。cursol 是部分解决方案（最初设置为 [] ）

def recSolve(n,sum, cursol):
    if n==1:
        print cursol + [sum]
        return
    if n == sum:
         print cursol + [1 for i in range(n)]
         return
    else:
        for i in range(1, sum-n+2):
             recsolve(n-1, sum-i, cursol+[i])

如果要使用动态规划，则必须记住 n 和 sum 的每个组合的解集。

于 2012-05-16T15:13:29.930 回答