Memecahkan sistem yang jarang dan sangat tidak terkondisikan

9

Saya bermaksud untuk menyelesaikan Ax = b di mana A adalah matriks persegi atau persegi panjang empat persegi panjang yang rumit, jarang, tidak simetris dan sangat tidak terkondisi (kondisi nomor ~ 1E + 20). Saya telah berhasil menyelesaikan sistem dengan ZGELSS di LAPACK secara akurat. Tetapi seiring dengan meningkatnya derajat kebebasan dalam sistem saya, perlu waktu lama untuk menyelesaikan sistem pada PC dengan ZGELSS karena sparsity tidak dieksploitasi. Baru-baru ini saya mencoba SuperLU (menggunakan penyimpanan Harwell-Boeing) untuk sistem yang sama tetapi hasilnya tidak akurat untuk nomor kondisi> 1E + 12 (Saya tidak yakin apakah ini merupakan masalah numerik dengan pivoting).

Saya lebih cenderung menggunakan pemecah yang sudah dikembangkan. Apakah ada pemecah yang kuat yang dapat memecahkan sistem yang saya sebutkan dengan cepat (yaitu mengeksploitasi sparsity) dan andal (mengingat angka-angka kondisi)?

pengguna1234
sumber
1
Bisakah Anda mengondisasinya? Jika demikian, metode ruang bagian Krylov bisa efektif. Bahkan jika Anda bersikeras pada metode langsung, prakondisi akan membantu mengendalikan kesalahan numerik.
Geoff Oxberry
1
Saya juga membuat pengalaman yang baik dengan memprakondisi cukup banyak seperti yang dijelaskan di sini: en.wikipedia.org/wiki/… Anda dapat melakukan prasyarat dalam aritmatika yang tepat. Matriks saya semuanya padat, jadi tidak bisa mengarahkan Anda ke metode / rutin yang lebih spesifik di sini.
AlexE
11
(machine precision)(condition number)1020__float128
2
Dari mana Anda mendapatkan perkiraan angka kondisi ini? Jika Anda meminta Matlab untuk memperkirakan jumlah kondisi dari sebuah matriks dengan ruang nol, itu mungkin memberi Anda tak terbatas atau kadang-kadang mungkin hanya memberi Anda jumlah yang sangat besar (seperti apa yang Anda miliki). Jika sistem yang Anda lihat memiliki ruang nol dan Anda tahu apa itu, Anda dapat memproyeksikannya dan yang tersisa mungkin memiliki nomor kondisi yang lebih baik. Kemudian Anda dapat menggunakan PETSc atau Trilinos atau apa pun.
Daniel Shapero
3
minAxbperp(N(A))

Jawaban:

13

Ketika Anda menggunakan ZGELSS untuk mengatasi masalah ini, Anda menggunakan dekomposisi nilai singular terpotong untuk mengatur masalah yang sangat tidak terkondisi ini. penting untuk memahami bahwa rutin pustaka ini tidak berusaha menemukan solusi kuadrat terkecil untuk , melainkan berusaha menyeimbangkan menemukan solusi yang meminimalkanterhadap meminimalkan. Ax=bxAxb

Perhatikan bahwa parameter RCOND yang diteruskan ke ZGELSS dapat digunakan untuk menentukan nilai singular mana yang harus dimasukkan dan dikecualikan dari perhitungan solusi. Nilai singular apa pun yang kurang dari RCOND * S (1) (S (1) adalah nilai singular terbesar) akan diabaikan. Anda belum memberi tahu kami bagaimana Anda mengatur parameter RCOND di ZGELSS, dan kami tidak tahu apa-apa tentang tingkat kebisingan dari koefisien dalam matriks Anda atau di sisi kanan , jadi sulit untuk mengatakan apakah Anda telah menggunakan jumlah regularisasi yang tepat. Ab

Anda tampaknya senang dengan solusi yang diregulasi yang Anda peroleh dengan ZGELSS, sehingga tampaknya regularisasi dipengaruhi oleh metode SVD terpotong (yang menemukan solusi minimum antara solusi kuadrat terkecil yang meminimalkan atas ruang solusi yang direntang oleh vektor singular yang terkait dengan nilai singular yang lebih besar dari RCOND * S (1)) memuaskan Anda. xAxb

Pertanyaan Anda dapat diformulasikan ulang sebagai "Bagaimana saya dapat secara efisien memperoleh solusi kuadrat terkecil yang terregulasi untuk masalah kuadrat linear terkecil yang besar, jarang, dan sangat tidak terkondisi ini?"

Rekomendasi saya adalah menggunakan metode berulang (seperti CGLS atau LSQR) untuk meminimalkan masalah kuadrat terkecil yang diatur secara eksplisit.

minAxb2+α2x2

di mana parameter regularisasi disesuaikan sehingga masalah kuadrat terkecil teredam dengan baik dan sehingga Anda senang dengan solusi yang diregulasi yang dihasilkan. α

Brian Borchers
sumber
Permintaan maaf saya karena tidak menyebutkan ini sejak awal. Masalah yang dipecahkan adalah persamaan Helmholtz dari akustik menggunakan FEM. Sistem ini dikondisikan dengan buruk karena dasar gelombang bidang yang digunakan untuk mendekati solusi.
user1234
Dari mana koefisien dan berasal? Apakah mereka mengukur data? "tepat" nilai-nilai dari desain beberapa objek (yang dalam praktiknya tidak bisa mesin untuk toleransi yang 15 digit ...)? Ab
Brian Borchers
1
Matriks A dan b dibentuk menggunakan formulasi lemah Helmholtz PDE, lihat: asadl.org/jasa/resource/1/jasman/v119/i3/…
user1234
9

Jed Brown telah menunjukkan hal ini dalam komentar untuk pertanyaan, tetapi sebenarnya tidak banyak yang dapat Anda lakukan dalam presisi ganda biasa jika jumlah kondisi Anda besar: dalam kebanyakan kasus, Anda mungkin tidak akan mendapatkan satu digit akurasi dalam solusi Anda dan, lebih buruk lagi, Anda bahkan tidak tahu karena Anda tidak dapat secara akurat mengevaluasi residu yang sesuai dengan vektor solusi Anda. Dengan kata lain: ini bukan pertanyaan tentang solver linier mana yang harus Anda pilih - tidak ada solver linier yang dapat melakukan sesuatu yang berguna untuk matriks seperti itu.

Situasi semacam ini biasanya terjadi karena Anda memilih basis yang tidak sesuai. Misalnya, Anda mendapatkan matriks yang dikondisikan buruk jika Anda memilih fungsi sebagai dasar metode Galerkin. (Ini mengarah ke matriks Hilbert, yang terkenal buruk kondisinya). Solusi dalam kasus-kasus seperti ini bukanlah dengan bertanya pemecah mana yang dapat memecahkan sistem linier, tetapi tanyakan apakah ada basis yang lebih baik yang dapat digunakan. Saya akan mendorong Anda untuk melakukan hal yang sama: pikirkan tentang merumuskan kembali masalah Anda sehingga Anda tidak berakhir dengan matriks semacam ini.1,x,x2,x3,...

Wolfgang Bangerth
sumber
Ketika mendiskritisasi masalah yang ditimbulkan buruk untuk PDE, misalnya persamaan panas terbalik, pasti kita akan berakhir dengan persamaan matriks yang tidak dikondisikan. Ini bukan kasus yang bisa kita pecahkan dengan merumuskan ulang persamaan atau memilih pemecah matriks yang efisien atau meningkatkan presisi dalam angka floating point. Dalam hal ini [yaitu masalah invers akustik], metode regularisasi diperlukan.
tqviet
7

Cara termudah / tercepat untuk menyelesaikan masalah yang tidak terkondisikan adalah untuk meningkatkan presisi perhitungan (dengan kekerasan). Cara lain (namun tidak selalu mungkin) adalah merumuskan kembali masalah Anda.

Anda mungkin perlu menggunakan presisi empat kali lipat (34 angka desimal). Meskipun 20 digit akan hilang dalam suatu kursus (karena nomor kondisi) Anda masih akan mendapatkan 14 digit yang benar.

Jika itu menarik, sekarang solver jarang quad-presisi tersedia di MATLAB juga.

(Saya penulis kotak alat yang disebutkan).

Pavel Holoborodko
sumber