2

多項係数を計算したい:

ここに画像の説明を入力

満たされているところn=n0+n1+n2

この演算子の Matlab 実装は、関数で簡単に実行できます。

function N = nchooseks(k1,k2,k3)
    N = factorial(k1+k2+k3)/(factorial(k1)*factorial(k2)*factorial(k3)); 
end

ただし、インデックスが 170 より大きい場合、階乗は無限になりNaN、場合によっては生成され180!/(175! 3! 2!) -> Inf/Inf-> NaNます。

他の投稿では、 CとPythonのこのオーバーフローの問題を解決しました。

  • Cの場合：「すべての階乗からリストを作成し、リスト内のすべての数字の素因数分解を見つけ、数字が完全に削減されるまで、一番上のすべての数字を一番下の数字でキャンセルできます" .
  • Python の場合: 「factorial(n) = gamma(n+1) という事実を利用し、ガンマ関数の対数を使用し、乗算の代わりに加算を使用し、除算の代わりに減算を使用します」 .

最初の解決策は非常に遅いように思われるため、2 番目のオプションを試しました。

function N = nchooseks(k1,k2,k3)
    N = 10^(log_gamma(k1+k2+k3)-(log_gamma(k1)+log_gamma(k2)+log_gamma(k3))); 
end
function y = log_gamma(x),  y = log10(gamma(x+1));  end

元の実装と log_gamma の実装を次のコードと比較します。

% Calculate
N=100; err = zeros(N,N,N);
for n1=1:N,
    for n2=1:N,
        for n3=1:N,
            N1 = factorial(n1+n2+n3)/(factorial(n1)*factorial(n2)*factorial(n3)); 
            N2 = 10^(log10(gamma(n1+n2+n3+1))-(log10(gamma(n1+1))+log10(gamma(n2+1))+log10(gamma(n3+1)))); 
            err(n1,n2,n3) = abs(N1-N2); 
        end
    end
end
% Plot histogram of errors
err_ = err(~isnan(err));
[nelements,centers] = hist(err_(:),1e2);
figure; bar(centers,nelements./numel(err_(:)));

ただし、次のヒストグラムに示すように、一部のケースでは結果がわずかに異なります。

ここに画像の説明を入力

したがって、私の実装が正しいと仮定する必要がありますか、それとも数値エラーが数の相違を正当化しないでしょうか?

4

4 に答える 4

2

古い投稿を復活させて申し訳ありませんが、将来の検索者のために、ほぼ確実に、多項係数を二項係数の積として記述し、組み込みメソッドを使用して二項係数を計算する必要があります (または、パスカルの三角形などを使用して独自の方法を記述します)方法）。関連する式は、ウィキペディアの多項係数に関するセクションの最初の段落に記載されています。（ここに書きたいのですが、LaTeXをレンダリングする方法がないようです。）

このアプローチのもう 1 つの利点は、係数がすべて整数であるため、オーバーフローを回避できることです。多項式の係数を計算する場合、本質的に除算する必要はありません。

于 2015-09-28T05:10:08.877 に答える
0

@jemidiah から提供されたヒントを使用して、

ここに画像の説明を入力

そしてここにコードがあります

function c = multicoeff (k), 
    c = 1; 
    for i=1:length(k), 
      c = c* bincoeff(sum(k(1:i)),k(i)); 
    end; 
end

およびいくつかの使用例:

octave:88> multicoeff([2 2 2])
ans =  90
octave:89> factorial(6)/(factorial(2)*factorial(2)*factorial(2))
ans =  90
octave:90> multicoeff([5 4 3])
ans =  27720
octave:91> factorial(12)/(factorial(5)*factorial(4)*factorial(3))
ans =  27720
于 2019-04-08T14:28:01.687 に答える