احسب محدد هاسكل لمصفوفة

قررت نشر رمز حساب المحددات. يعمل الكود ، على الرغم من أنه لا يتظاهر بأنه مبدع. كان حل هذه المشكلة في هاسكل أمرًا مثيرًا للاهتمام. يتم النظر في طريقتين لحل المشكلة: العودية البسيطة والطريقة الغاوسية.

قليلا من النظرية

كما تعلم ، محدد المصفوفة المربعة n * n هو مجموع n! المصطلحات ، كل منها عبارة عن منتج يحتوي على عنصر مصفوفة واحد بالضبط من كل عمود وواحد بالضبط من كل صف. علامة القطعة التالية:

$$

$${a}_{1,i1}*{a}_{2,i2}*...{a}_{n,in}$$

يتم تحديده من خلال تعادل الاستبدال:

$$\begin{pmatrix}1 & 2 &… & n \\ {i}_{1} & {i}_{2} &… & {i}_{n} \end{pmatrix}$$



تتمثل الطريقة المباشرة لحساب المحدد في تحليله بواسطة عناصر صف أو عمود إلى مجموع منتجات عناصر صف أو عمود بواسطة مكملاتهم الجبرية. بدوره ، المكمل الجبري لعنصر المصفوفة

$$

$${a}_{i,j}$$

يوجد

$$

$${(-1)}^{i+j}*{M}_{i,j}$$

حيث

$$

$${M}_{i,j}$$

- يوجد عنصر ثانوي (i، j) أي المحدد الذي تم الحصول عليه من المحدد الأصلي عن طريق حذف الصف الأول والعمود ي.



تولد هذه الطريقة عملية تكرارية تسمح بحساب أي محدد. لكن أداء هذه الخوارزمية يترك الكثير مما هو مرغوب فيه - O (n!). لذلك ، يتم استخدام مثل هذا الحساب المباشر ، باستثناء الحسابات الرمزية (ومع محددات ليست عالية الترتيب).



تبين أن طريقة غاوس أكثر إنتاجية. يرتكز جوهرها على الأحكام التالية:



1. محدد المصفوفة المثلثية العليا $$\begin{pmatrix}{a}_{1,1} & {a}_{1,2} &… & {a}_{1,n} \\ 0 & {a}_{2,2} &… & {a}_{2,n} \\ 0 & 0 &… & ...\\ 0 & 0 &… & {a}_{n,n} \\\end{pmatrix}$$ يساوي حاصل ضرب عناصره القطرية. تأتي هذه الحقيقة مباشرة من تحلل المحدد من حيث عناصر الصف الأول أو العمود الأول.



2. إذا أضفت عناصر صف آخر في المصفوفة لعناصر صف واحد ، مضروبة في نفس الرقم ، فلن تتغير قيمة المحدد.



3. إذا تم تبادل صفين (أو عمودين) في المصفوفة ، فإن قيمة المحدد ستغير علامتها إلى العكس.



يمكننا ، باختيار المعاملات ، إضافة الصف الأول من المصفوفة مع جميع الصفوف الأخرى والحصول على الأصفار في العمود الأول في جميع المواضع باستثناء الأول. للحصول على صفر في السطر الثاني ، عليك إضافة الأول مضروبًا في السطر الثاني

$$

$$-{a}_{2,1}/{a}_{1,1}$$

للحصول على صفر في السطر الثالث ، أضف السطر الأول مضروبًا في

$$

$$-{a}_{3,1}/{a}_{1,1}$$

إلخ في النهاية ، سيتم تقليل المصفوفة إلى شكل تكون فيه جميع العناصر

$$

$${a}_{n,1}$$

لـ n> 1 يساوي صفرًا.



إذا كان العنصر في المصفوفة

$$

$${a}_{1,1}$$

تبين أنه يساوي صفرًا ، ثم يمكنك العثور على عنصر غير صفري في العمود الأول (لنفترض أنه في المكان k) وقم بتبديل الصفين الأول والصف k. مع هذا التحول ، يغير المحدد ببساطة علامته ، والتي يمكن أخذها في الاعتبار. إذا لم تكن هناك عناصر غير صفرية في العمود الأول ، فإن المحدد يساوي صفرًا.



علاوة على ذلك ، بالتصرف بطريقة مماثلة ، يمكنك الحصول على أصفار في العمود الثاني ، ثم في العمود الثالث ، إلخ. من المهم ألا تتغير الأصفار التي تم الحصول عليها مسبقًا عند إضافة السلاسل. إذا لم يكن من الممكن العثور على عنصر غير صفري للمقام لأي سطر ، فإن المحدد يساوي الصفر ويمكن إيقاف العملية. ينتج عن الإكمال الطبيعي لعملية Gaussian مصفوفة تكون فيها جميع العناصر الموجودة أسفل القطر الرئيسي مساوية للصفر. كما هو مذكور أعلاه ، فإن محدد مثل هذه المصفوفة يساوي حاصل ضرب العناصر القطرية.

دعنا ننتقل إلى البرمجة.

نحن نعمل مع بيانات الفاصلة العائمة. نحن نمثل المصفوفات كقوائم من السلاسل. أولاً ، دعنا نحدد نوعين:

type Row    = [Double]
type Matrix = [Row]

العودية البسيطة

لا مزيد من التردد ، سنقوم بتوسيع المحدد بواسطة عناصر الصف الأول (أي صفر). نحتاج إلى برنامج لبناء الصغرى ، نحصل عليه بحذف الصف الأول والعمود k.

--  k-     
deln :: Matrix -> Int -> Matrix
deln matrix k = map (\ r -> (take (k) r)++(drop (k+1) r)) matrix

وهنا الصغير:

--  k-   
minor :: Matrix -> Int -> Double
minor matrix k = det $ deln (drop 1 matrix) k

يرجى ملاحظة: القاصر هو المحدد. نحن نطلق على وظيفة det ، التي لم ننفذها بعد. لتنفيذ det ، سيتعين علينا تكوين المجموع البديل لنواتج العنصر التالي من الصف الأول بواسطة محدد القاصر التالي. لتجنب التعبيرات المرهقة ، دعنا ننشئ دالة منفصلة لتشكيل علامة الجمع:

sgn :: Int -> Double
sgn n = if n `rem` 2 == 0 then 1.0 else (-1.0)

الآن يمكننا حساب المحدد:

--   
det :: Matrix -> Double
det [[a,b],[c,d]] = a*d-b*c
det matrix = sum $ map (\c -> ((matrix !! 0)!!c)*(sgn c)*(minor matrix c))  [0..n]
             where n = length matrix - 1

الكود بسيط للغاية ولا يتطلب أي تعليقات خاصة. لاختبار أداء وظائفنا ، دعنا نكتب الوظيفة الرئيسية:

main = print $ det [[1,2,3],[4,5,6],[7,8,(-9)]]

قيمة هذا المحدد 54 ، والتي يمكن التحقق منها.

طريقة جاوس

سنحتاج إلى بعض وظائف المرافق (التي يمكن استخدامها في مكان آخر). الأول هو تبادل صفين في المصفوفة:

--    
swap :: Matrix -> Int -> Int -> Matrix
swap matrix n1 n2 = map row [0..n]
                    where n=length matrix - 1
                          row k | k==n1 = matrix !! n2
                                | k==n2 = matrix !! n1
                                | otherwise = matrix !! k

كما ترى من الكود أعلاه ، فإن الوظيفة تمر سطراً بسطر. في هذه الحالة ، إذا تم العثور على خط برقم n1 ، يتم استبدال السطر n2 بالقوة (والعكس صحيح). بقيت بقية الخطوط في مكانها.



الدالة التالية تحسب السلسلة r1 المضافة بالسلسلة r2 مضروبًا في العنصر بواسطة الرقم f:

--   r1+f*r2
comb :: Row -> Row -> Double -> Row
comb r1 r2 f = zipWith (\ x y -> x+f*y) r1 r2

كل شيء هنا شفاف للغاية: يتم تنفيذ الإجراءات على صفوف المصفوفة (أي على قوائم [مزدوجة]). لكن الوظيفة التالية تؤدي هذا التحول على المصفوفة (وبطبيعة الحال ، تحصل على مصفوفة جديدة):

--    r1  r2,   f
trans :: Matrix -> Int -> Int -> Double -> Matrix
trans matrix n1 n2 f = map row [0..n]
                       where n=length matrix - 1
                             row k | k==n1 = comb (matrix !! n1) (matrix !! n2) f
                                   | otherwise = matrix !! k

تبحث الدالة getNz عن رقم العنصر الأول غير الصفري في القائمة. تكون مطلوبة عندما يساوي العنصر القطري التالي صفرًا.

--     
getNz :: Row -> Int
getNz xs = if length tmp == 0 then (-1) else snd $ head tmp
           where tmp=dropWhile (\ (x,k) -> (abs x) <= 1.0e-10) $ zip xs [0..]

إذا كانت جميع عناصر القائمة صفرًا ، فستُرجع الدالة -1.



تتحقق وظيفة البحث مما إذا كانت المصفوفة مناسبة للتحول التالي (يجب أن تحتوي على عنصر قطري تالٍ غير صفري). إذا لم يكن كذلك ، يتم تحويل المصفوفة عن طريق تبديل الصفوف.

--        
search :: Matrix -> Int -> Matrix
search matrix k | (abs ((matrix !! k) !! k)) > 1.0e-10 = matrix
                | nz < 0 = matrix  --      
                | otherwise = swap matrix k p 
                           where n   = length matrix
                                 lst = map (\ r -> r !! k) $ drop k matrix
                                 nz  = getNz lst
                                 p   = k + nz

إذا كان من المستحيل العثور على العنصر الأول (غير الصفري) (المصفوفة متدهورة) ، فإن الوظيفة ستعيدها دون تغيير. تشكل الدالة mkzero أصفارًا في العمود التالي من المصفوفة:

--     
mkzero :: Matrix -> Int -> Int -> Matrix
mkzero matrix k p | p>n-1 = matrix
                  | otherwise = mkzero (trans matrix p k (-f)) k (p+1)
                    where n = length matrix
                          f = ((matrix !! p) !! k)/((matrix !! k) !! k)

تشكل وظيفة المثلث الشكل المثلث العلوي للمصفوفة:

--     
triangle :: Matrix -> Int -> Matrix
triangle matrix k | k>=n = matrix
                  | (abs v) <= 1.0e-10 = [[0.0]] --  
                  | otherwise = triangle (mkzero tmp k k1) k1 
                    where n   = length matrix
                          tmp = search matrix k
                          v   = (tmp !! k) !! k --  
                          k1  = k+1

إذا لم يكن من الممكن في المرحلة التالية العثور على العنصر الرئيسي ، فإن الدالة ترجع مصفوفة صفرية من الترتيب الأول. الآن يمكنك تكوين وظيفة العرض لتقليل المصفوفة إلى الشكل المثلثي العلوي:

--  
gauss :: Matrix -> Matrix
gauss matrix = triangle matrix 0 

لحساب المحدد ، علينا ضرب العناصر القطرية. للقيام بذلك ، دعنا ننشئ وظيفة منفصلة:

--   
proddiag :: Matrix -> Double
proddiag matrix = product $ map (\ (r,k) -> r !!k) $ zip matrix [0,1..]

حسنًا ، و "القوس" هو الحساب الفعلي للمحدد:

--  
det :: Matrix -> Double
det matrix = proddiag $ triangle matrix 0

دعنا نتحقق من كيفية عمل هذه الوظيفة:

main = print $ det  [[1,2,3],[4,5,6],[7,8,-9]]

[1 of 1] Compiling Main             ( main.hs, main.o ) 
Linking a.out ...                                                                                 
54.0     

شكرا لمن قرأ حتى النهاية!



يمكن تنزيل الكود هنا