Menguji apakah dua matriks 12x12 memiliki determinan yang sama

11

12×12QJ

det(Q)=det(12IQJ)(1)
J

Saat ini saya melakukan ini dengan perpustakaan armadillo tetapi ternyata terlalu lambat. Masalahnya adalah bahwa saya perlu melakukan ini untuk satu triliun matriks dan ternyata menghitung dua penentu adalah hambatan dari program saya. Karena itu saya punya dua pertanyaan

  1. Apakah ada trik yang bisa saya gunakan untuk menghitung penentu lebih cepat mengingat saya tahu ukurannya? Apakah mungkin ekspansi berantakan untuk matriks yang dapat bekerja dalam kasus ini?12×12

  2. Apakah ada cara lain yang efisien untuk menguji kesetaraan(1)

Edit. Untuk menjawab komentar. Saya perlu menghitung semua graf saling melengkapi yang tidak terhubung sendiri dengan urutan sehingga dan memiliki jumlah pohon rentang yang sama. Motivasi untuk ini dapat ditemukan di pos mathoverflow ini . Sedangkan untuk mesin saya menjalankan ini pada mesin 8 core 3.4GHh secara paralel.13 G ¯ GG13GG¯

Edit. Saya dapat mengurangi waktu berjalan yang diharapkan sebesar 50% dengan membuat program C untuk secara spesifik menghitung determinan matriks . Saran masih diterima.12×12

Jernej
sumber
6
Seberapa lambat terlalu lambat? Berapa lama waktu yang dibutuhkan untuk perangkat keras apa? Apakah triliunan ini independen sehingga Anda dapat menghitung banyak dari penentu ini secara paralel? Jika demikian, seberapa besar mesin yang dapat Anda jalankan? Apa yang menyebabkan masalah ini? Apakah Anda yakin harus menghitung faktor penentu? Q
Bill Barth
3
Seberapa sering (untuk fraksi kasus apa) determinannya sama / berbeda? Jika mereka berbeda sebagian besar waktu, mungkin ada tes yang lebih murah untuk menentukan bahwa mereka mungkin berbeda dan Anda akan memverifikasi bahwa mereka sama hanya jika tes pertama gagal. Sebaliknya jika mereka hampir sama.
Wolfgang Bangerth
1
Seperti yang sudah ditanyakan: dapatkah Anda memberikan beberapa detail dari mana berasal? Mungkin ada pendekatan yang lebih baik daripada faktor penentu komputasi buta. Q
JM
4
Gagasan bahwa kondisi ini harus diuji "untuk satu triliun matriks" menunjukkan 1) bahwa diketahui apriori memiliki beberapa struktur khusus (jika tidak, harapan bahwa kondisi berlaku secara acak sedikit), dan 2) bahwa pendekatan yang lebih baik mungkin untuk mengkarakterisasi semua matriks dengan properti ini (dengan formulasi yang dapat diperiksa secara efisien). QQQ
hardmath
1
@hardmath Ya, adalah matriks integer yang memiliki entri diagonal mulai dari hingga dan sebagai elemen diagonal1 12 - 1Q1121
Jernej

Jawaban:

8

Karena Anda sudah menggunakan C ++ dan matriks Anda pasti positif simetris, saya akan melakukan faktorisasi tidak dipvotivasikan dari dan juga . Di sini saya berasumsi bahwa juga pasti positif, jika tidak, akan memerlukan pivoting untuk stabilitas numerik (juga mungkin meskipun itu tidak pasti positif, pivoting tidak diperlukan, tetapi Anda harus mencobanya). Q 12 I - Q - J 12 I - Q - J L D L TLDLTQ12IQJ12IQJLDLT

Ini lebih cepat daripada faktorisasi LU, dan juga lebih cepat daripada Cholesky karena akar kuadrat dihindari. Penentu hanyalah produk dari elemen-elemen matriks diagonal . Kode untuk melakukan faktorisasi LDL sangat sederhana sehingga Anda dapat menulisnya dalam waktu kurang dari 50 baris C. Halaman Wikipedia di dalamnya menjelaskan algoritme, dan saya memiliki beberapa kode templated sederhana untuk melakukan Cholesky di sini . Anda dapat menyederhanakannya dan memodifikasinya untuk menghindari akar kuadrat untuk mengimplementasikan faktorisasi .L D L TDLDLT

Karena Anda juga dapat mengontrol format penyimpanan, Anda dapat lebih mengoptimalkan rutin untuk menyimpan hanya setengah matriks dan mengemasnya dalam array linier untuk memaksimalkan memori lokalitas. Saya juga akan menulis produk kustom titik sederhana dan rutinitas pembaruan peringkat-1 karena ukuran masalahnya sangat kecil, Anda harus membiarkan kompiler menyatukan rutinitas untuk mengurangi overhead panggilan. Karena itu adalah loop ukuran tetap, kompiler harus dapat secara otomatis inline dan membuka gulungan barang bila perlu.

Saya akan menghindari mencoba memainkan trik untuk mengambil keuntungan dari kenyataan bahwa berisi di dalam ekspresi. Sangat mungkin bahwa untuk ukuran masalah sekecil itu, trik-trik ini pada akhirnya menjadi lebih lambat daripada hanya melakukan dua perhitungan penentu terpisah. Tentu saja, satu-satunya cara untuk memverifikasi klaim ini adalah dengan mencobanya.Q12IQJQ

Victor Liu
sumber
1
Saya mendukung rekomendasi penerapan , alias Cholesky bebas-root, karena Armadillo tampaknya tidak memiliki cara untuk mengambil keuntungan dari kepastian positif / dominasi diagonal. LDLT
hardmath
5

Tanpa beberapa informasi tentang konstruksi matriks simetrik nyata pasti positif ini, saran-saran yang akan dibuat adalah kebutuhan yang cukup terbatas.12×12

Saya mengunduh paket Armadillo dari Sourceforge dan melihat dokumentasi. Cobalah untuk meningkatkan kinerja komputasi secara terpisah dan , di mana adalah matriks peringkat satu dari semua yang ada, dengan menetapkan misalnya . Dokumentasi mencatat bahwa ini adalah default untuk matriks hingga ukuran , jadi karena kelalaian saya menganggap opsi adalah default untuk kasus .det(Q)det(12IQJ)Jdet(Q,slow=false)4×4slow=true12×12

Apa yang slow=true mungkin dilakukan adalah berputar sebagian atau penuh dalam mendapatkan bentuk eselon baris, dari mana determinan mudah ditemukan. Namun Anda tahu sebelumnya bahwa matriks adalah positif pasti, jadi pivot tidak diperlukan untuk stabilitas (setidaknya secara dugaan untuk sebagian besar perhitungan Anda. Tidak jelas apakah paket Armadillo melempar pengecualian jika pivot menjadi terlalu kecil, tetapi ini harus menjadi fitur paket aljabar linear numerik yang masuk akal. EDIT: Saya menemukan kode Armadillo yang mengimplementasikan dalam file header , menggunakan template C ++ untuk fungsionalitas substansial. Pengaturan ini tampaknya tidak mempengaruhi bagaimanaQdetinclude\armadillo_bits\auxlib_meat.hppslow=false12×12penentu akan dilakukan karena perhitungan akan "dilempar ke atas tembok" ke LAPACK (atau ATLAS) pada saat itu tanpa indikasi bahwa pivot tidak diperlukan; lihat det_lapackdan doanya dalam file itu.

Poin lainnya adalah mengikuti rekomendasi mereka untuk membangun paket Armadillo yang menghubungkan ke penggantian kecepatan tinggi untuk BLAS dan LAPACK, jika Anda memang menggunakannya; lihat Sec. 5 file Armadillo README.TXT untuk detailnya. [Penggunaan BLAS atau LAPACK versi 64-bit khusus juga disarankan untuk kecepatan pada mesin 64-bit saat ini.]

Reduksi baris ke bentuk eselon pada dasarnya adalah eliminasi Gaussian, dan memiliki kompleksitas aritmatika . Untuk kedua matriks ini, berarti dua kali lipat yang berfungsi, atau . Operasi-operasi ini mungkin menjadi "penghambat" dalam pemrosesan Anda, tetapi ada sedikit harapan bahwa tanpa struktur khusus di (atau beberapa hubungan yang diketahui di antara triliun kasus uji yang memungkinkan amortisasi) pekerjaan dapat dikurangi menjadi .23n3+O(n2)43n3+O(n2)QO(n2)

Sebagai perbandingan, ekspansi oleh kofaktor umum matriks melibatkanoperasi perkalian (dan kira-kira sebanyak penambahan / pengurangan), jadi untuk perbandingan ( vs ) jelas mendukung penghapusan lebih dari kofaktor.n×nn!n=1212!=47900160023n3=1152

Pendekatan lain yang membutuhkan bekerja akan mengurangi ke bentuk tridiagonal dengan transformasi Householder, yang juga menempatkan ke dalam bentuk tridiagonal. Komputasi dan selanjutnya dapat dilakukan dalam operasi . [Efek dari pembaruan peringkat satu dalam penentu kedua dapat dinyatakan sebagai faktor skalar yang diberikan dengan menyelesaikan satu sistem tridiagonal.]43n3+O(n2)Q12IQdet(Q)det(12IQJ)O(n)J

Menerapkan perhitungan independen semacam itu mungkin bermanfaat sebagai pemeriksaan atas hasil panggilan yang berhasil (atau gagal) ke detfungsi Armadillo .

Kasus Khusus: Seperti yang disarankan oleh Komentar Jernej, anggaplah bahwa mana seperti sebelumnya adalah matriks (peringkat 1) dari semua yang ada dan adalah matriks diagonal nonsingular (positif). Memang untuk aplikasi yang diusulkan dalam teori grafik ini akan menjadi matriks bilangan bulat. Maka rumus eksplisit untuk adalah:Q=DJJD=diag(d1,,dn)det(Q)

det(Q)=(i=1ndi)(1i=1ndi1)

Sebuah sketsa pembuktiannya memberikan peluang untuk menggambarkan penerapan yang lebih luas, yaitu setiap kali memiliki determinan yang diketahui dan sistem dengan cepat diselesaikan. Mulailah dengan memperhitungkan:DDv=(11)T

det(DJ)=det(D)det(ID1J)

Sekarang kembali menjadi peringkat 1, yaitu . Perhatikan bahwa penentu kedua sederhana:D1J(d11dn1)T(11)

f(1)=det(ID1J)

di mana adalah polinomial karakteristik dari . Sebagai matriks peringkat 1, harus memiliki (setidaknya) faktor untuk memperhitungkan ruang nolnya. Nilai eigen "hilang" adalah , seperti yang dapat dilihat dari perhitungan:f(x)D1Jf(x)n1xdi1

D1J(d11dn1)T=(di1)(d11dn1)T

Karena itu karakteristik polinomial , dan adalah seperti yang ditunjukkan di atas untuk , .f(x)=xn1(xdi1)f(1)det(ID1J)1di1

Perhatikan juga bahwa jika , maka , sebuah matriks diagonal yang determinannya hanyalah produk dari entri diagonalnya.Q=DJ12IQJ=12ID+JJ=12ID

hardmath
sumber
Hm .. sebenarnya adalah mana adalah matriks adjacency dari jadi saya berpikir bahwa hasil ini mungkin tidak benar. Khususnya itu akan menyiratkan bahwa jumlah pohon spanning dari grafik ditentukan oleh urutan derajatnya yang tidak berlaku. QDAAGG
Jernej
Entri off-diagonal biasanya akan berisi 0 dan -1. The dekomposisi disarankan oleh Victor mengambil keuntungan dari simetri dan mengurangi jangka terkemuka di count operasi dari untuk . Ada pendekatan integer yang tepat, tetapi itu mungkin tidak diperlukan untuk matriks dan entri ukuran sederhana Anda. Jika saya mengerti konstruksinya, adalah pasti positif untuk alasan yang sama adalah. QLDLT23n313n312IQJQ
hardmath
@Jernej: Jika Anda yakin sesuatu yang saya nyatakan tidak benar, saya telah membuat ruang obrolan berdasarkan pada Pertanyaan ini tempat diskusi dapat di utas tanpa komentar yang tidak perlu di sini.
hardmath
1

Jika Anda memiliki cara terstruktur untuk menghitung grafik yang ingin Anda hitung determinannya, mungkin Anda bisa menemukan pembaruan peringkat rendah yang memindahkan Anda dari satu grafik ke grafik lainnya.

Jika demikian, maka Anda dapat menggunakan lemma penentu matriks untuk menghitung dengan murah determinan grafik berikutnya yang akan dihitung menggunakan pengetahuan Anda tentang determinan grafik saat ini.

Yaitu, untuk matriks dan vektor : Ini dapat digeneralisasi jika U dan V adalah matriks dan adalah : Au,v

det(A+uvT)=(1+vTA1u)det(A)
n×mAn×n
det(A+UVT)=det(Im+VTA1U)det(A)

Untuk menghitung invers secara efisien, Anda dapat menggunakan rumus Sherman-Morrison untuk mendapatkan invers dari matriks berikutnya dari matriks saat ini:

(A+uvT)1=A1A1uvTA11+vTA1u

Richard
sumber