Apa yang dimaksud dengan "symplectic" dalam referensi ke integrator numerik, dan apakah odeint SciPy menggunakannya?

Dalam komentar ini saya menulis:

... integrator SciPy default, yang saya asumsikan hanya menggunakan metode symplectic.

di mana saya merujuk ke SciPy's odeint, yang menggunakan metode "non-stiff (Adams)" atau "stiff (BDF) method". Menurut sumber :

def odeint(func, y0, t, args=(), Dfun=None, col_deriv=0, full_output=0,
           ml=None, mu=None, rtol=None, atol=None, tcrit=None, h0=0.0,
           hmax=0.0, hmin=0.0, ixpr=0, mxstep=0, mxhnil=0, mxordn=12,
           mxords=5, printmessg=0):
    """
    Integrate a system of ordinary differential equations.

    Solve a system of ordinary differential equations using lsoda from the
    FORTRAN library odepack.

    Solves the initial value problem for stiff or non-stiff systems
    of first order ode-s::
        dy/dt = func(y, t0, ...)
    where y can be a vector.
    """

Berikut adalah contoh di mana saya menyebarkan orbit satelit di sekitar bumi selama tiga bulan hanya untuk menunjukkan bahwa satelit itu melakukan seperti yang diharapkan.

Saya percaya bahwa integrator non-simplektik memiliki sifat yang tidak diinginkan sehingga mereka cenderung tidak menghemat energi (atau jumlah lain) dan karenanya tidak diinginkan dalam mekanika orbital misalnya. Tapi saya tidak begitu yakin apa yang membuat integrator symplectic symplectic.

Apakah mungkin untuk menjelaskan apa properti itu (yang membuat integrator symplectic symplectic) secara langsung dan (cukup) mudah dimengerti tetapi tidak dengan cara yang tidak akurat? Saya bertanya dari sudut pandang bagaimana integrator berfungsi secara internal , daripada bagaimana kinerjanya dalam pengujian.

Dan apakah kecurigaan saya benar bahwa odeinthanya menggunakan integrator symplectic?

Sebagai aturan praktis yang kuat, Anda seharusnya hanya berharap integrator kotak hitam adalah symplectic jika mengharuskan Anda untuk memisahkan persamaan posisi dan momentum.

origimbo

@origimbo Terima kasih. Ya, dan sepertinya odeintpembungkus Python untuk kode sumber yang cukup lama, mapan, dan direferensikan dengan baik, (pertanyaan diedit, referensi ODEPACK dan LSODA) meskipun saya tentu mengakui menggunakannya dalam mode kotak hitam. Contoh saya yang ditautkan menunjukkan vektor keadaan 6D terdiri dari tiga posisi dan tiga kecepatan.

uhoh

Integrator ODE dalam ODEPACK dan LSODA bukan integrator simpleks .

Brian Borchers

Berikut adalah contoh yang berhasil membandingkan dua pemecah yang sangat sederhana: Euler dan Symplectic Euler: idontgetoutmuch.wordpress.com/2013/08/06/… .

idontgetoutmuch

Buku karya Hairer, Nørsett, dan Wanner memberikan penjelasan yang baik tentang metode symplectic. Lihatlah Gambar 16.1 secara khusus, dan gambar di sini .

Jawaban:

Biarkan saya memulai dengan koreksi. Tidak, odeinttidak memiliki integrator symplectic. Tidak, integrasi symplectic tidak berarti konservasi energi.

Apa yang dimaksud dengan symplectic dan kapan Anda harus menggunakannya?

Pertama-tama, apa arti symplectic? Symplectic berarti bahwa solusinya ada pada manifold symplectic. Manifold symplectic adalah set solusi yang didefinisikan oleh bentuk-2. Detail manifold symplectic mungkin terdengar seperti matematis omong kosong, jadi intinya adalah ada hubungan langsung antara dua set variabel pada manifold tersebut. Alasan mengapa ini penting untuk fisika adalah karena persamaan Hamiltonian secara alami memiliki bahwa solusi berada pada manifold symplectic dalam ruang fase, dengan pemisahan alami menjadi komponen posisi dan momentum. Untuk solusi Hamiltonian sejati, jalur ruang fase itu adalah energi konstan.

Integrator symplectic adalah integrator yang solusinya berada pada manifold symplectic. Karena kesalahan diskritisasi, ketika memecahkan sistem Hamiltonian itu tidak mendapatkan lintasan yang benar pada manifold. Sebaliknya, yang lintasan itu sendiri terganggu untuk pesanan dari lintasan yang benar. Lalu ada pergeseran linear karena kesalahan numerik lintasan ini dari waktu ke waktu. Integrator normal cenderung mengalami pergeseran kuadratik (atau lebih), dan tidak memiliki jaminan global yang baik tentang jalur ruang fase ini (hanya lokal). $\mathcal{O}(\Delta t^n)$ $n$

Ini cenderung berarti bahwa integrator symplectic cenderung menangkap pola lama lebih baik daripada integrator normal karena kurangnya penyimpangan dan ini hampir menjamin periodisitas. Notebook ini menampilkan properti tersebut dengan baik pada masalah Kepler . Gambar pertama menunjukkan apa yang saya bicarakan dengan sifat periodik dari solusi.

Ini diselesaikan dengan menggunakan integrator symplectic urutan ke - 6 dari Kahan dan Li dari DifferentialEquations.jl . Anda dapat melihat bahwa enerji tidak dilestarikan dengan tepat, tetapi variasinya bergantung pada seberapa jauh solusi berlipat ganda dari manifold yang sebenarnya. Tetapi karena solusi numerik itu sendiri berada pada manifold symplectic, cenderung hampir secara berkala (dengan beberapa penyimpangan numerik linier yang dapat Anda lihat), membuatnya menjadi sangat baik untuk integrasi jangka panjang. Jika Anda melakukan hal yang sama dengan RK4, Anda bisa mendapatkan bencana:

Anda dapat melihat bahwa masalahnya adalah tidak ada periodisitas sejati dalam solusi numerik dan karenanya lembur cenderung melayang.

Ini menyoroti alasan sebenarnya untuk memilih integrator symplectic: integrator symplectic baik pada integrasi lama pada masalah yang memiliki properti symplectic (sistem Hamiltonian) . Jadi mari kita bahas beberapa hal. Perhatikan bahwa Anda tidak selalu memerlukan integrator symplectic bahkan pada masalah symplectic. Untuk kasus ini, metode Runge-Kutta ke-5 adaptif dapat dilakukan dengan baik. Inilah Tsit5:

Perhatikan dua hal. Pertama, ia mendapatkan akurasi yang cukup baik sehingga Anda tidak dapat melihat pergeseran sebenarnya dalam plot ruang fase. Namun, di sisi kanan Anda dapat melihat bahwa ada penyimpangan energi ini, dan jika Anda melakukan integrasi yang cukup lama, metode ini tidak akan berfungsi sebaik metode solusi dengan properti periodik. Tapi itu menimbulkan pertanyaan, bagaimana cara efisiensi-vs vs hanya mengintegrasikan sangat akurat? Yah, ini agak kurang pasti. Dalam DiffEqBenchmarks.jl Anda dapat menemukan beberapa tolok ukur menyelidiki pertanyaan ini. Misalnya, notebook inimelihat kesalahan energi vs runtime pada sistem persamaan Hamiltonian dari model Boson quadruple dan menunjukkan bahwa jika Anda menginginkan akurasi yang sangat tinggi, maka bahkan untuk waktu integrasi yang cukup lama, lebih efisien jika hanya menggunakan RK or Runge-Kutta Nystrom yang tinggi ( Metode RKN). Ini masuk akal karena untuk memenuhi sifat simptektik integrator memberikan beberapa efisiensi dan cukup banyak harus diperbaiki langkah waktu (ada beberapa penelitian membuat kemajuan ke yang terakhir tetapi tidak terlalu jauh).

Selain itu, perhatikan dari kedua notebook ini bahwa Anda juga dapat mengambil metode standar dan memproyeksikannya kembali ke solusi yang berlipat ganda setiap langkah (atau setiap beberapa langkah). Inilah yang dilakukan contoh-contoh menggunakan panggilan balik ManifoldProjection DifferentialEquations.jl . Anda melihat bahwa jaminan hukum konservasi ditegakkan tetapi dengan biaya tambahan untuk menyelesaikan sistem implisit setiap langkah. Anda juga dapat menggunakan pemecah ODE yang sepenuhnya implisit atau matriks massa tunggal untuk menambah persamaan konservasi, tetapi hasil akhirnya adalah bahwa metode ini lebih mahal secara komputasi sebagai pengorbanan.

Jadi untuk meringkas, kelas masalah di mana Anda ingin meraih integrator symplectic adalah mereka yang memiliki solusi pada manifold symplectic (sistem Hamilton) di mana Anda tidak ingin menginvestasikan sumber daya komputasi untuk memiliki yang sangat tepat (toleransi <1e-12) solusi dan tidak membutuhkan energi yang tepat / dll. konservasi. Ini menyoroti bahwa ini semua tentang properti integrasi jangka panjang, jadi Anda tidak harus berbondong-bondong ke mereka semua mau tak mau seperti yang ditunjukkan beberapa literatur. Tetapi mereka masih merupakan alat yang sangat penting di banyak bidang seperti Astrofisika di mana Anda memiliki integrasi lama yang Anda butuhkan untuk menyelesaikan cukup cepat tanpa memiliki akurasi yang absurd.

Di mana saya menemukan integrator symplectic? Integrator symplectic seperti apa yang ada?

Pada umumnya ada dua kelas integrator symplectic. Ada integrator Runge-Kutta symplectic (yang ditunjukkan pada contoh di atas) dan ada metode Runge-Kutta implisit yang memiliki properti symplectic. Seperti @origimbo menyebutkan, integrator Runge-Kutta symplectic mengharuskan Anda memberi mereka struktur yang dipartisi sehingga mereka dapat menangani posisi dan momentum secara terpisah. Namun, bertentangan dengan komentar tersebut, metode Runge-Kutta implisit bersifat simplektik tanpa memerlukan ini, tetapi sebaliknya membutuhkan penyelesaian sistem nonlinier. Ini tidak terlalu buruk karena jika sistemnya tidak kaku, sistem nonlinear ini dapat diselesaikan dengan iterasi fungsional atau percepatan Anderson, tetapi metode simptektik RK mungkin masih lebih disukai untuk efisiensi (itu '

Yang mengatakan, odeint tidak memiliki metode dari salah satu keluarga ini , jadi itu bukan pilihan yang baik jika Anda mencari integrator symplectic. Di Fortran, situs Hairer memiliki set kecil yang dapat Anda gunakan . Mathematica memiliki beberapa bawaan . Para pemecah ODL GSL memiliki integrator titik Gaussian RK implisit yang IIRC bersifat symplectic, tapi itu satu-satunya alasan untuk menggunakan metode GSL.

Tetapi set integrator symplectic paling komprehensif dapat ditemukan di DifferentialEquations.jl di Julia (ingat ini digunakan untuk notebook di atas). Daftar metode Runge-Kutta symplectic yang tersedia ditemukan di halaman ini dan Anda akan melihat bahwa metode titik tengah implisit juga symplectic (metode Trapesium Runge-Kutta implisit dianggap "hampir symplectic" karena dapat dibalik). Tidak hanya memiliki serangkaian metode terbesar, tetapi juga sumber terbuka (Anda dapat melihat kode dan pengujiannya dalam bahasa tingkat tinggi) dan memiliki banyak tolok ukur . Notebook pengantar yang bagus untuk menggunakannya untuk menyelesaikan masalah fisik adalah notebook tutorial ini. Tapi tentu saja disarankan Anda memulai dengan paket melalui tutorial ODE pertama .

Secara umum Anda dapat menemukan analisis terperinci dari suite persamaan diferensial numerik di posting blog ini . Ini cukup terperinci, tetapi karena harus mencakup banyak topik, masing-masing kurang detail dari ini, jadi jangan ragu untuk memintanya diperluas dengan cara apa pun.

Chris Rackauckas
sumber

Dengan jawaban ini saya sepertinya telah mencapai Stack Exchange Jackpot! Ini adalah jawaban yang sempurna bagi saya, karena saya segera memahaminya dan bagian-bagian yang tidak membuat saya ingin membaca lebih lanjut. Saya sangat menghargai waktu yang Anda gunakan untuk mencari jawaban ini secara menyeluruh serta menyertakan tautan bermanfaat dan instruktif lainnya.

uhoh

Sebelum masuk ke perincian matematika, secara kasar kita bisa mengatakan bahwa symplectic berarti pelestarian volume , bukan?

Miguel

FTR, alasan urutan ke-5 Runge-Kutta yang adaptif berkinerja jauh lebih baik daripada RK4 di sini bukan karena memiliki urutan yang lebih tinggi tetapi memilih ukuran langkah yang lebih cocok. Alasan RK4 berkinerja sangat buruk adalah terutama karena ukuran langkahnya terlalu tinggi di perigee; pemecah yang sama dengan setengah ukuran langkah akan memberikan solusi yang jauh lebih baik. (Hanya, itu akan menghabiskan banyak waktu menyelesaikan orbit dengan halus di sekitar apogee, di mana ini tidak diperlukan sama sekali.)

leftaroundabout

eksposisi yang sangat baik. Sebagai pertanyaan sampingan: OP dimulai dengan referensi ke Python - apakah ada tutorial / paket Python yang direkomendasikan di sepanjang garis contoh Julia yang ditautkan?

Quetzalcoatl

Satu-satunya paket Python yang saya tahu untuk jenis integrator ini adalah diffeqpy , di mana itu tidak didokumentasikan pada README tetapi Anda dapat mengakses semua metode yang sama dan menulis ulang ini dalam Python menggunakan paket itu.

Chris Rackauckas

$\mathbf{p}$ $\mathbf{q}$ $\mathcal{H(\mathbf{p},\mathbf{q})}$

\frac{d q}{d t} = + \frac{\partial H}{\partial hal}

$\frac{d\mathbf{q}}{dt}=+\frac{\partial \mathcal{H}}{\partial\mathbf{p}}$

\frac{d hal}{d t} = - \frac{\partial H}{\partial q} .

$\frac{d\mathbf{p}}{dt}=-\frac{\partial \mathcal{H}}{\partial\mathbf{q}}.$

H

$\mathcal{H}$

hal (t), q (t) = ϕ_{t} (hal (t_{0}), q (t_{0}))

$\mathbf{p}(t),\mathbf{q}(t)= \phi_t(\mathbf{p}(t_0),\mathbf{q}(t_0))$

d p \land d q

$d\mathbf{p}\wedge d\mathbf{q}$

p

$p$

q

$q$ adalah satu dimensi Anda dapat menganggap ini sebagai mengatakan bahwa area di dalam kurva tertutup pada ruang fase dilestarikan. Ini memastikan semua jenis properti stabilitas yang bagus, karena "bola" lintasan harus tetap "dekat" satu sama lain.

Dalam hal angka, integrator symplectic bertindak dengan cara yang sama, juga melestarikan area / dua bentuk ini. Pada gilirannya ini berarti bahwa ada "Hamiltonian numerik" yang dikonservasi (yang mungkin tidak [dibaca 'tidak'] sama dengan yang persis). Perhatikan bahwa stabilitas tidak sama dengan keakuratan, sehingga sebagian besar keuntungan dari metode simpleks datang ketika mengintegrasikan untuk waktu yang sangat lama (misalnya metode Anda dapat dengan cepat menempatkan satelit di sisi Bumi yang salah, sementara tidak pernah membiarkannya membusuk menjadi saya t).

origimbo
sumber

Terima kasih untuk ini! Sekarang saya akan menggunakan kata-kata di atas nilai gaji saya. n-bola lintasan lebih berisiko ketika mereka dekat bifurkasi seperti yang ada dalam simulasi 3-tubuh. lih. Doedel et al. 2007, Int. J. Bifurcation dan Chaos, v 17, no. 8 (2007) 2625–2677 Bagaimana saya melakukannya? Juga ieec.cat/hosted/web-libpoint/papers/…

uhoh

Kecuali jika pembaca mengetahui rincian matematisnya, penyebutan stabilitas itu menyesatkan, karena konservasi volume tidak berarti bahwa lintasan individual tetap dekat.

Miguel

@Miguel Saya pikir ini adalah salah satu situasi di mana pembaca yang tidak mengikuti detail matematika dibodohi dengan cara baik, tetapi dalam hal troika numerik biasa tentang akurasi, stabilitas dan efisiensi komputasi, saya berpendapat bahwa menekankan stabilitas manfaatnya bermanfaat. Saya senang menerima saran untuk menulis ulang jika Anda dapat memikirkan sesuatu yang lebih baik yang tidak sengaja tidak akurat.

origimbo

^{2}

$^2$

^{2}

$^2$

@Miguel: Tapi gumpalan partikel dibiarkan membelah menjadi dua bagian atau lebih. Its Total Volume hanya harus tetap konstan.

Wolfgang Bangerth