Bagaimana cara menghitung residu secara numerik?

15

Saya perlu menghitung integral berikut: f ( E ) = T r

12πiCf(E)dE
Di manahadalah matriks (satu partikel kinetik dan energi potensial dinyatakan dalam basis),Gadalah matriks yang bergantung padaE(fungsi satu-partikel banyak-tubuh Green) dan integral kontur adalah setengah lingkaran kiri. Integrandf(E)memiliki kutub pada sumbu nyata negatif dan mahal untuk mengevaluasi. Apa cara paling efektif untuk menghitung integral seperti itu?
f(E)=Tr((h+E)G(E))
hGEf(E)

Inilah penelitian saya sejauh ini:

1) Saya menggunakan integrasi Gaussian, jalur integrasi saya adalah persegi panjang. Saya memperbaiki sisi kiri dan kanan (yaitu lebar) dan bermain dengan tinggi (di atas dan di bawah sumbu nyata) sehingga untuk urutan integrasi yang diberikan saya mendapatkan akurasi tertinggi. Misalnya untuk pesanan 20, jika ketinggiannya terlalu besar, akurasinya turun (jelas), tetapi jika terlalu kecil, ia juga turun (teori saya adalah bahwa ia membutuhkan lebih banyak titik di sekitar kutub saat ketinggian naik ke 0). Saya puas dengan ketinggian optimal 0,5 untuk fungsi saya.

2) Lalu saya mengatur sisi kanan persegi panjang pada E0, biasanya E0 = 0, tetapi bisa juga E0 = -0.2 atau yang serupa.

3) Saya mulai memindahkan sisi kiri dari persegi panjang ke kiri dan untuk setiap langkah saya melakukan integrasi urutan konvergensi untuk memastikan integral saya sepenuhnya terkonvergensi untuk setiap persegi panjang. Dengan menambah lebar, saya akhirnya mendapatkan nilai konvergen dalam batas setengah lingkaran kiri tak hingga.

Perhitungannya sangat lambat dan juga tidak terlalu akurat untuk lebar yang besar. Salah satu perbaikannya adalah dengan mempartisi lebar yang panjang menjadi "elemen" dan menggunakan integrasi Gaussian pada setiap elemen (seperti dalam FE).

Pilihan lain adalah mengintegrasikan lingkaran kecil di sekitar masing-masing kutub dan merangkumnya. Masalah:

a) Bagaimana cara menemukan kutub fungsi numerik ? Itu harus kuat. Satu-satunya yang saya tahu adalah bahwa mereka berada di poros nyata negatif. Untuk beberapa dari mereka (tapi tidak semua) saya juga tahu tebakan awal yang cukup bagus. Apakah ada metode yang berfungsi untuk fungsi analitik f ( E ) ? Atau apakah itu tergantung pada bentuk aktual dari f ( E ) ?f(E)f(E)f(E)

b) Setelah kita tahu kutub, skema numerik apa yang terbaik untuk mengintegrasikan lingkaran kecil di sekitarnya? Haruskah saya menggunakan integrasi Gaussian pada lingkaran? Atau haruskah saya menggunakan distribusi poin yang seragam?

Pilihan lain mungkin sekali saya tahu kutub berkat a), mungkin ada beberapa cara semi-analitik untuk mendapatkan Residu tanpa perlu integrasi yang kompleks. Tetapi untuk sekarang saya akan senang hanya mengoptimalkan integrasi kontur.

OndČej Čertík
sumber
1
Sudahkah Anda memeriksa buku "Metode Numerik untuk Pembalikan Transformasi Laplace" oleh Cohen (2007)? IIRC, Robert Piessens (dari ketenaran QUADPACK) juga menangani topik ini.
GertVdE

Jawaban:

7

Saya dapat menawarkan saran untuk pertanyaan pertama Anda: Jika Anda tahu kutub Anda berada di suatu tempat di sepanjang sumbu nyata, Anda bisa melokalkannya cukup efisien menggunakan interpolasi / pendekatan Rasional . Ini sama dengan menemukan polinomial dan q ( x ) sedemikian rupap(x)q(x)

f(x)p(x)q(x)

untuk dalam beberapa interval. Tiang-tiang f ( x ) kemudian harus cocok dengan akar q ( x )xf(x) q(x) .

Interpolasi / aproksimasi rasional bisa menjadi hal yang rumit, tetapi saya baru-baru ini ikut menulis makalah tentang algoritma yang stabil untuk menghitungnya menggunakan SVD. Makalah ini berisi kode Matlab yang mengimplementasikan algoritme, dan versi yang lebih luas tersedia sebagai fungsi ratinterpdalam proyek Chebfun , di mana saya adalah salah satu pengembangnya.

Untuk pertanyaan kedua Anda, makalah ini mungkin bermanfaat.

Pedro
sumber
Terima kasih untuk semua tipsnya! Berikut adalah kode netlib.org/toms/579 dari makalah Bengt Fornberg. Sayangnya, ada beberapa bug numerik, karena ini adalah output yang saya dapatkan: gist.github.com/2942970#file_output . Jadi saya harus mengimplementasikan ulang atau men-debug-nya. Tautan Chebfun memberi saya 404 (saya mencobanya beberapa bulan yang lalu dengan hasil yang sama, jadi mungkin itu tidak bekerja dari Amerika Serikat).
Ondřej Čertík
@ OndřejČertík: Saya sendiri belum pernah menggunakan kode TOMS 579, jadi saya tidak tahu harus memberi tahu Anda tentang kesalahan apa. Sedangkan untuk beranda Chebfun, dapatkah Anda mencoba "googling" dan melihat apakah itu berfungsi?
Pedro
Google menemukan beranda Chebfun dan menampilkan versi dalam tembolok. Tetapi ketika saya mengklik halaman tersebut, inilah yang saya dapatkan: pastehtml.com/view/c1ts4h3ct.html
Ondřej Čertík
Coba peramban lain? Atau dari ISP yang berbeda. Situs web berfungsi dengan baik dari sini (di AS.)
Costis
Saya mencoba Firefox dan Chrome. Jadi itu harus oleh ISP saya. Aneh.
Ondřej Čertík