2014-01-08 11 views
5

Juliaでは、2つの有理数行列の行列除算が浮動小数点行列を返します。代わりに有理数の行列を得る方法はありますか?Juliaの合理的な行列分割

convert(Array{Rational}, A \ b)を使用することはできません。浮動小数点数に関連する精度が低下するためです。

答えて

5

有理行列のためにピボットされたQR factorizationアルゴリズムを実装する必要がありますが、これはかなり単純ではありませんが、不可能ではありません。 Juliaは、これを行うためにLAPACK DGEQP3 functionを64ビット浮動小数点行列に使用しています。あなたが合理的なQR分解を働かせても、それはLAPACKアルゴリズムのどこよりも速くはありません。だから私はあなたが正確な合理的な行列の分割が必要なのか尋ねるだろうと思う。

脇に:julia-users listのように質問をするほうが実り多いかもしれません。あなたはそれについての会話をすることができるからです。これは本当に「尋ねられ、答えられた」問題ではありません。

更新:あなたはRational{Int128}を使用する必要がありますので、

julia> A = [rand(1:10)//rand(1:10) for i=1:4, j=1:4] 
4x4 Array{Rational{Int64},2}: 
5//3 1//2 10//1 8//9 
1//1 3//2 6//7 2//3 
4//1 8//9 6//7 10//3 
7//2 5//2 1//2 5//1 

julia> b = collect(1:4) 
4-element Array{Int64,1}: 
1 
2 
3 
4 

julia> c = A\b 
4-element Array{Rational{Int64},1}: 
    42055//62497 
    61344//62497 
    -2954//62497 
-19635//124994 

julia> A*c 
4-element Array{Rational{Int64},1}: 
1//1 
2//1 
3//1 
4//1 

は、Rational{Int}がオーバーフローすることはかなり傾向があること、しかし、注意してください:ジェネリックに旋回QRは今、ジュリアに存在するためこれは今「だけの作品」オーバーフローを避けるにはRational{BigInt}を使用してください。このアルゴリズムは完全に汎用的であると同様に、よりエキゾチックな数値型の作品:これはかなり興味深い四元数の乗算は可換ではないことを

julia> using Quaternions # install with Pkg.add("Quaternions") 

julia> A = [Quaternion(rand(1:10), rand(1:10), rand(1:10), rand(1:10)) for i=1:4, j=1:4] 
4x4 Array{Quaternions.Quaternion{Int64},2}: 
    4 + 3im + 5jm + 8km 9 + 7im + 10jm + 3km 9 + 3im + 1jm + 7km 8 + 4im + 8jm + 5km 
    1 + 4im + 2jm + 1km 5 + 4im + 1jm + 4km 1 + 5im + 8jm + 2km 7 + 2im + 5jm + 3km 
10 + 1im + 4jm + 4km 2 + 4im + 9jm + 2km 2 + 10im + 4jm + 10km 2 + 3im + 4jm + 8km 
    7 + 4im + 3jm + 7km 8 + 3im + 5jm + 9km 6 + 5im + 1jm + 3km 10 + 10im + 3jm + 1km 

julia> b = collect(1:4) 
4-element Array{Int64,1}: 
1 
2 
3 
4 

julia> c = A\b 
4-element Array{Quaternions.Quaternion{Float64},1}: 
    0.18112 + 0.019288355350921868im - 0.002638049486498678jm + 0.11047233503816825km 
0.000578119 + 0.08073035854610844im - 0.023278758601757744jm - 0.16761078955242836km 
    -0.0501257 - 0.02428891715971441im - 0.11059096793480247jm - 0.059017235311989824km 
    0.0394953 - 0.16771397199827004im - 0.019823106776170954jm + 0.05251791790026253km 

julia> A*c 
4-element Array{Quaternions.Quaternion{Float64},1}: 
            1.0 + 1.1102230246251565e-16im + 0.0jm + 0.0km 
            2.0 + 2.220446049250313e-16im + 0.0jm + 0.0km 
3.0 + 4.440892098500626e-16im - 2.220446049250313e-16jm + 8.881784197001252e-16km 
4.0 + 2.220446049250313e-16im - 2.220446049250313e-16jm + 6.661338147750939e-16km 

julia> norm(b - A*c) 
1.680072297996111e-15 

注意、。

+1

これは単なる問題です。 Juliaを使用して、高次積分法を導き出します。係数は、線形システムの解であるAx = bです。正確な有理係数は、浮動小数点係数よりもはるかに数学的解析に適しています。 –

+0

BLAS/LAPACK機能のかなりの部分集合の良い合理的な類似点を持つことは確かに興味深いプロジェクトですが、私が着手している人は誰もわかっていません。 Mathematicaを除いて、このようなことをしている既存のシステムはわかりませんが、私はそれを試していませんが、私はそれを証明することはできません。合理的なオーバーフローは、この種の事業にとって大きな問題となる可能性が高い。 – StefanKarpinski

+0

julia-usersにこの質問を投稿しました:https://groups.google.com/forum/#!topic/julia-users/RLqRiHm-MpA –