Difference between revisions of "Oscillating one-dimensional systems"

From ccitonlinewiki
Jump to: navigation, search
(= 4.3.7 Metode Runge-Kutta Orde 4)
 
(151 intermediate revisions by 32 users not shown)
Line 1: Line 1:
== Studi kasus dan Terjemahan ==
+
<- back to [[Studi kasus komputasi teknik]]
 +
 
 +
== Knowledge base ==
 +
 
 +
 
 +
 
 +
== Studi kasus ==
  
 
[[File:1d oscillating dynamic system 1.png]]
 
[[File:1d oscillating dynamic system 1.png]]
Line 120: Line 126:
 
====4.3.3 Memprogram Metode Numerik; Kasus Khusus====
 
====4.3.3 Memprogram Metode Numerik; Kasus Khusus====
  
 +
Program sederhana untuk (4.47) - ( 4.48) mengikuti ide yang sama seperti di bagian 4.2.3:
 +
 +
[[File:4.3.3.fadhli.JPG|500px]]
 +
 +
(Lihat file osc_FE.py.)
  
 +
Karena kita sudah tahu solusi yang tepat sebagai u(t) = Xo cos ωt , kami beralasan sebagai berikut untuk menemukan interval simulasi yang sesuai [0,T] dan juga berapa poin kita harus memilih. Solusinya memiliki periode P = 2π/ω. (Periode P adalah waktunya perbedaan antara dua puncak u(t) ~ cos ωt curve). Simulasi untuk tiga periode fungsi cosinus, T = 3P, dan memilih Δt sehingga ada 20 Interval per periode menghasilkan Δt = P/20 dan total Nt  = T/ Δt = t interval. Sisanya dari program ini adalah pengodean langsung dari skema Forward Euler.
 +
 +
Gambar 4.16 menunjukkan perbandingan antara solusi numerik dan tepat solusi persamaan diferensial. Yang mengejutkan kami, solusi numeriknya terlihat salah. Apakah perbedaan ini disebabkan oleh kesalahan pemrograman atau masalah dengan metode Forward Euler?
  
 
Pertama-tama, bahkan sebelum mencoba menjalankan program, Anda harus menghitung dua langkah dalam putaran waktu dengan kalkulator sehingga Anda memiliki beberapa hasil antara untuk dibandingkan. Menggunakan X0 = 2. Dt = 0: 157079632679, dan ω = 2, kita mendapatkan u1 = 2, v = -1,25663706, u2 = 1,80260791, dan v2 = 2,51327412. Perhitungan semacam itu menunjukkan bahwa program itu tampaknya benar. (Kemudian, kita dapat menggunakan nilai-nilai tersebut untuk membangun tes unit dan fungsi tes yang sesuai.)
 
Pertama-tama, bahkan sebelum mencoba menjalankan program, Anda harus menghitung dua langkah dalam putaran waktu dengan kalkulator sehingga Anda memiliki beberapa hasil antara untuk dibandingkan. Menggunakan X0 = 2. Dt = 0: 157079632679, dan ω = 2, kita mendapatkan u1 = 2, v = -1,25663706, u2 = 1,80260791, dan v2 = 2,51327412. Perhitungan semacam itu menunjukkan bahwa program itu tampaknya benar. (Kemudian, kita dapat menggunakan nilai-nilai tersebut untuk membangun tes unit dan fungsi tes yang sesuai.)
Line 131: Line 145:
  
 
[[File:Simulation with different steps.PNG|500px]]
 
[[File:Simulation with different steps.PNG|500px]]
 
  
 
==== '''4.3.4 Sebuah Penyelesaian dari Metode Numerik ''' ====
 
==== '''4.3.4 Sebuah Penyelesaian dari Metode Numerik ''' ====
Line 137: Line 150:
 
Dalam skema Forward Euler,
 
Dalam skema Forward Euler,
  
kita dapat mengganti un pada persamaan terakhir dengan nilai unC1 yang baru dihitung dari
+
[[File:4.3.4.(1).JPG|500px]]
 +
 
 +
kita dapat mengganti u^n pada persamaan terakhir dengan nilai u^n+1 yang baru dihitung dari
 
persamaan pertama:
 
persamaan pertama:
 +
 +
[[File:4.3.4.(2).JPG|500px]]
  
 
Sebelum membenarkan perbaikan ini secara matematis, mari kita coba pada contoh sebelumnya. Hasilnya muncul pada Gambar 4.18. Kita melihat bahwa amplitudo tidak tumbuh, tetapi
 
Sebelum membenarkan perbaikan ini secara matematis, mari kita coba pada contoh sebelumnya. Hasilnya muncul pada Gambar 4.18. Kita melihat bahwa amplitudo tidak tumbuh, tetapi
Line 148: Line 165:
 
Mari kita tafsirkan skema yang disesuaikan secara matematis. Pertama kami memesan (4,49) - (4,50)
 
Mari kita tafsirkan skema yang disesuaikan secara matematis. Pertama kami memesan (4,49) - (4,50)
 
sedemikian rupa sehingga perbedaan pendekatan terhadap derivatif menjadi transparan:
 
sedemikian rupa sehingga perbedaan pendekatan terhadap derivatif menjadi transparan:
(4,51)
 
  
 +
[[File:4.3.4.(10).JPG|500px]]
 +
 +
[[File:4.3.4.(3).JPG|500px]]
  
 (4,52)
 
 
Kami menafsirkan (4,51) sebagai persamaan diferensial sampel pada titik mesh tn, karena
 
Kami menafsirkan (4,51) sebagai persamaan diferensial sampel pada titik mesh tn, karena
 
kami memiliki vn di sisi kanan. Sisi kiri kemudian perbedaan maju atau
 
kami memiliki vn di sisi kanan. Sisi kiri kemudian perbedaan maju atau
 
Meneruskan perkiraan Euler ke turunan u0
 
Meneruskan perkiraan Euler ke turunan u0
 
, lihat Gambar 4.2. Di samping itu,
 
, lihat Gambar 4.2. Di samping itu,
kami menginterpretasikan (4,52) sebagai persamaan diferensial sampel pada titik mesh tnC1, karena kami miliki di sisi kanan. Dalam hal ini, perbedaan aproksimasi pada
+
kami menginterpretasikan (4,52) sebagai persamaan diferensial sampel pada titik mesh tnC1, karena kami miliki di sisi kanan.  
 +
 
 +
[[File:4.3.4.(4).jpeg]]
 +
 
 +
Dalam hal ini, perbedaan aproksimasi pada
 
sisi kiri adalah perbedaan ke belakang,
 
sisi kiri adalah perbedaan ke belakang,
  
 +
[[File:4.3.4.(5).jpeg]]
  
  
Line 164: Line 187:
 
Gambar 4.19 mengilustrasikan perbedaan mundur. Kesalahan dalam perbedaan mundur sebanding dengan t, sama seperti untuk perbedaan maju (tetapi konstanta proporsionalitas dalam istilah kesalahan memiliki tanda yang berbeda). Diskretisasi yang dihasilkan
 
Gambar 4.19 mengilustrasikan perbedaan mundur. Kesalahan dalam perbedaan mundur sebanding dengan t, sama seperti untuk perbedaan maju (tetapi konstanta proporsionalitas dalam istilah kesalahan memiliki tanda yang berbeda). Diskretisasi yang dihasilkan
 
metode untuk (4,52) sering disebut sebagai skema Backward Euler.
 
metode untuk (4,52) sering disebut sebagai skema Backward Euler.
 +
 
Untuk meringkas, gunakan perbedaan maju untuk persamaan pertama dan mundur
 
Untuk meringkas, gunakan perbedaan maju untuk persamaan pertama dan mundur
 
Perbedaan untuk hasil persamaan kedua dalam metode yang jauh lebih baik daripada hanya menggunakan
 
Perbedaan untuk hasil persamaan kedua dalam metode yang jauh lebih baik daripada hanya menggunakan
 
maju perbedaan dalam kedua persamaan.
 
maju perbedaan dalam kedua persamaan.
 +
 
Cara standar untuk mengekspresikan skema ini dalam fisika adalah dengan mengubah urutan
 
Cara standar untuk mengekspresikan skema ini dalam fisika adalah dengan mengubah urutan
 
persamaan,
 
persamaan,
  
 +
[[File:4.3.4.(6).jpeg]]
  
 
dan terapkan perbedaan maju ke (4,53) dan perbedaan mundur ke (4,54):
 
dan terapkan perbedaan maju ke (4,53) dan perbedaan mundur ke (4,54):
  
=== 4.3.7 Metode Runge-Kutta Orde 4 ===
+
[[File:4.3.4.(7).jpg]]
 +
 
 +
Artinya, pertama kecepatan v diperbarui dan kemudian posisi u, menggunakan kecepatan yang paling baru dihitung. Tidak ada perbedaan antara (4,55) - (4,56) dan (4,49) -
 +
(4,50) sehubungan dengan akurasi, jadi urutan persamaan diferensial asli
 +
tidak apa-apa. Skema (4.55) - (4.56) berada di bawah nama Semi-implisit
 +
Euler4 atau Euler-Cromer. Implementasi (4.55) - (4.56) ditemukan dalam file
 +
osc_EC.py. Inti dari kode itu seperti
 +
 
 +
[[File:4.3.4.(8).jpg]]
 +
 
 +
 
 +
[[File:4.3.4.(9).jpg]]
 +
 
 +
 
 +
 
 +
 
 +
 
 +
==== 4.3.5 Metode Runge-Kutta orde 2 (atau metode Heun) ====
 +
Sebuah metode yang cukup populer digunakan untuk menyelesaikan persamaan diferensial biasa (ODE) vector dan skalar adalah metode Runge-Kutta Order (RK2), atau biasa dikenal dengan metode Heun. Ide dasar pada metode ini, yang pertama untuk ODE skalar, adalah dengan membentuk aproksimasi perbedaan terpusat (centered difference) terhadap turunan antara dua titik waktu yang didefinisikan sebagai berikut:
 +
 
 +
[[File:Kolab1.JPG]]
 +
 
 +
Formula dari centered difference tersebut dapat digambarkan melalui Gambar 4.20. Error pada aproksimasi centered difference ini proporsional terhadap nilai ∆t2, 1 order lebih tinggi dibandingkan dengan pendekatan forward and backward difference, yang berarti nilai jika kita memiliki sebuah nilai ∆t, maka error nya akan berkurang secara effektif dengan menggunakan centered difference karena nilai error tersebut berkurang dengan faktor 4, daripada faktor 2.
 +
 
 +
[[File:Kolab2.JPG]]
 +
 
 +
Permasalahan yang ada pada skema centered difference semacam ini untuk persamaan ODE secara umum, u’=f(u,t) adalah kita mendapatkan
 +
 
 +
[File:Kolab3.JPG]]
 +
 
 +
Yang mana ini akan menyulitkan karena kita tidak mengetahui berapa nilai un+1/2. Namun demikian, ktia dapat mengaproksimasi nilai f diantara dua level waktu dengan menggunakan rata-rata aritmatik dari nilai f tesebut pada saat tn dan tn+1  :
 +
 
 +
[[FIle:Kolab4.JPG]]
 +
 
 +
Kemudian hasilnya adalah :
 +
[[File:results435.jpg]]
 +
Dimana berupa persamaan aljabar nonlinear untuk
 +
[[File:f435.jpeg]]
 +
dan bukan fungsi linear dari u.
 +
sehingga untuk menyelesaikan fungsi
 +
[[File:f4351.jpg]]
 +
tanpa menyelesaikannya dengan persamaan nonlinear, dapat diprediksi [[File:f4352.jpg]]
 +
menggunakan persamaan Forward Euler:
 +
[[File:f4353.jpg]]
 +
Sehingga dapat digunakan metode
 +
[[File:f4354.jpg]]
 +
metode tersebut dapat diaplikasikan untuk ODEs scalar dan vector.
 +
 
 +
Untuk system osilasi dengan
 +
[[File:f4355.jpg]]
 +
 
 +
Pada file osc_Heun.py terdapat implementasinya. File tersebut menjalankan simulasi untuk 10 period dengan 20 kali langkah per periode.
 +
 
 +
Solusi Numerical dan eksak yang berkaitan dengan ini terdapat di fig. 4.21. dapat diliat bahwa amplitude meningkat namun tidak sebanyak pada metode forward euler. Bgaimanapun juga, metode forward euler adalah yang terbaik.
 +
Perlu diingat juga bahwa metode forward euler memberikan prediksi yang lebih baik, seperti contohnya untuk persoalan pertumbuhan/peluruhan, atau SIR mode. Akan tetapi metode orde 2 runge-kutta atau metod heun juga bisa dipertimbangkan. Meskipun untuk menyelesaikan persoalan osilasi, metode euler sudah terbaik.
 +
 
 +
 
 +
 
 +
==== 4.3.6 Perangkat Lunak untuk Menyelesaikan ODEs ====
 +
 
 +
Terdapat banyak metode yang dapat digunakan untuk menyelesaikan ODEs, dan alangkah baiknya kita memilih akses yang mudah untuk mengimplementasikannya ke berbagai metode, terutama metode  adaptif yang canggih dan kompleks yang dapat menyesuaikan nilai Δt secara otomatis untuk mendapatkan nilai akurasi yang ditentukan. Phyton Odespy3 merupakan salah satu perangkat yang dapat memberikan akses yang mudah ke berbagai metode numerik untuk menyelesaikan ODEs.
 +
 
 +
Salah satu contoh termudah dalam penggunaan Odespy adalah untuk menyelesaikan masalah u’ = u, u(0) = 2, untuk 100 time steps sampai t = 4:
 +
 
 +
import odespy
 +
 
 +
def f(u, t):
 +
    return u
 +
 
 +
method = odespy.Heun  #or, e.g., odespy.ForwardEuler
 +
solver  = method(f)
 +
solver.set_initial _condition(2)
 +
time_points = np.linspace(0, 4, 101)
 +
u. t = solver.solve (time_points)
 +
 
 +
Dengan kata lain, kalian mendefinisikan sebuah fungsi f(u, t), menginisialisasi sebuah objek penyelesaian Odespy, mengatur kondisi awal, menghitung titik waktu pengumpulan dimana anda menginginkan solusinya, dan bertanya mengenai solusinya. Variabel arrays u dan t dapat dibuat menjadi sebuah grafik secara langsung, yaitu: plot(t,u).
 +
 
 +
Fitur menarik yang dimiliki oleh Odespy ialah parameter permasalahan dapat menjadi sebuah argumen pada fungsi f(u, t) penggunanya. Sebagai contoh, apabila permasalahan ODE kita adalah u’ = -au + b, dengan 2 parameter yaitu a dan b, kita dapat menuliskan fungsi f kita menjadi
 +
 
 +
def f(u, t, a, b):
 +
    return -a*u + b
 +
 
 +
Sebagai tambahan, permasalahan yang bergantung pada argumen a dan b dapat ditransfer ke fungsi ini bila kita mengumpulkan nilainya dalam sebuah daftar atau tuple ketika membuat sebuah pemecahan Odespy dan menggunakan argumen f_args:
 +
 
 +
a = 2
 +
b = 1
 +
solver = method(f, f_args=[a, b])
 +
 
 +
Hal ini merupakan sebuah fitur yang baik karena parameter permasalahan haruslah selain sebagai sebuah variabel global – sekarang dapat menjadi sebuah argument dalam fungsi kita secara alami.
 +
 
 +
Menggunakan Odespy untuk menyelesaikan osilasi ODEs seperti u” + ω2u = 0, diformulasikan sebagai sebuah sistem u’ = v dan v’ = -ω2u, dilakukan sebagai berikut. Kita tentukan sebuah nilai time steps per periode dan hitung time steps yang diasosiasikan serta waktu akhir simulasi (T), cantumkan sebuah nilai periode untuk disimulasikan:
 +
 
 +
Import odespy
 +
 
 +
# Define the ODE system
 +
# u’ = v
 +
# v’ = -omega**2*u
 +
 
 +
def f(sol, t, omega=2):
 +
      u, v = sol
 +
      return [v, -omega**2*u]
 +
 
 +
#Set and compute problem dependent parameters
 +
omega = 2
 +
X_0 = 1
 +
number_of_periods = 40
 +
time_intervals_per_period = 20
 +
from numpy import pi, linspace, cos
 +
P = 2*pi/omega #length of one period  # length of one period
 +
dt = P/time_intervals_per_period      # time step
 +
T = number_of_periods*P                # final simulation time
 +
 
 +
# Create Odespy solver object
 +
odespy_method = odespy.RK2
 +
solver = odespy_method(f, f_args=[omega])
 +
 
 +
# The initial condition for the system is collected in a list
 +
Solver.set_initial_condition([X_0, 0])
 +
 
 +
# Compute the desired time points where we want the solution
 +
N_t = int(round(T/dt)) # no of time intervals
 +
Time_points = linspace(0, T, N_t+1)
 +
 
 +
# Solve the ODE problem
 +
sol, t = solver.solve(time_points)
 +
 
 +
# Note: sol contains both displacement and velocity
 +
# extract original variables
 +
u = sol[:,0]
 +
v = sol[:,1]
 +
 
 +
Dua pernyataan terakhir menjadi penting karena dua fungsi u dan v di dalam sistem ODE tersebut tergabung bersama dalam sebuah array di dalam pemecahan Odespy. Solusi pada sistem ODE ditunjukan sebagai array 2 dimesi dimana kolom pertama (sol[:,0]) disimpan sebagai u dan kolom kedua (sol[:,1]) disimpan sebagai v. Mengeplot u dan v merupakan sebuah masalah dalam menjalankan plot(t, u, t, v).
 +
 
 +
Catatan
 +
 
 +
Di dalam fungsi tersebut kita menuliskan f(sol, t, omega) dibandingkan menulis f (u, t, omega) untuk mengindikasikan bahwa solusi pada f adalah solusi pada waktu t dimana nilai u dan t tergabung bersama: sol = [u,v]. Kita dapat juga menggunakan u sebagai argumen:
 +
 
 +
def f(u, t, omega=2):
 +
      u, v = u
 +
      return [v, -omega**2*u]
 +
 
 +
Ini hanya berarti kita mendefinisikan ulang nama u pada fungsi tersebut untuk merata-ratakan solusi pada waktu t untuk komponen pertama pada sistem ODE tersebut.
 +
 
 +
Untuk beralih ke metode numerik lain, tinggal substitusikan RK2 dengan nama yang sesuai dari metode yang diinginkan. Mengetik pydoc odespy pada terminal window memunculkan daftar dari metode yang dijalankan. Cara yang sangat sederhana dalam memilih metode ini menyarankan penambahan yang jelas dari kode diatas: kita dapat menentukan daftar metode, menjalankan semua metode, dan membandingkan setiap kurva u pada sebuah plot. Sebagaimana odespy juga mengandung skema Euler-Cromer, kita menulis kembali sistem ini dengan v’ = -w2u sebagai ODE pertama dan u’ = v sebagai ODE kedua, karena ini adalah pilihan standar ketika menggunakan metode Euler-Cromer (juga pada odespy):
 +
 
 +
def f(u, t, omega=2):
 +
v, u = u
 +
return [-omega**2*u, v]
 +
 
 +
Perubahan persamaan ini juga mempengaruhi kondisi awal: komponen pertama adalah nol dan yang kedua adalah X_0 maka kita perlu melewati daftar [0, X_0] untuk solver.set_ initial_condition.
 +
 
 +
Kode osc_odespy.py mengandung detail:
 +
 
 +
def compare(odespy_methods,
 +
omega,
 +
X_0,
 +
number_of_periods,
 +
time_intervals_per_period=20):
 +
from numpy import pi, linspace, cos
 +
P = 2*pi/omega                  # length of one period
 +
dt = P/time_intervals_per_period
 +
T = number_of_periods*P
 +
# If odespy_methods is not a list, but just the name of
 +
# a single Odespy solver, we wrap that name in a list
 +
# so we always have odespy_methods as a list
 +
if type(odespy_methods) != type([]):
 +
odespy_methods = [odespy_methods]
 +
# Make a list of solver objects
 +
solvers = [method(f, f_args=[omega]) for method in
 +
odespy_methods]
 +
for solver in solvers:
 +
solver.set_initial_condition([0, X_0])
 +
# Compute the time points where we want the solution
 +
dt = float(dt) # avoid integer division
 +
N_t = int(round(T/dt))
 +
time_points = linspace(0, N_t*dt, N_t+1)
 +
legends = []
 +
for solver in solvers:
 +
sol, t = solver.solve(time_points)
 +
v = sol[:,0]
 +
u = sol[:,1]
 +
# Plot only the last p periods
 +
p = 6
 +
m = p*time_intervals_per_period # no time steps to plot
 +
plot(t[-m:], u[-m:])
 +
hold(’on’)
 +
legends.append(solver.name())
 +
xlabel(’t’)
 +
# Plot exact solution too
 +
plot(t[-m:], X_0*cos(omega*t)[-m:], ’k--’)
 +
legends.append(’exact’)
 +
legend(legends, loc=’lower left’)
 +
axis([t[-m], t[-1], -2*X_0, 2*X_0])
 +
title(’Simulation of %d periods with %d intervals per period’
 +
% (number_of_periods, time_intervals_per_period))
 +
savefig(’tmp.pdf’); savefig(’tmp.png’)
 +
show()
 +
 
 +
Fitur baru pada kode ini adalah kemampuan untuk mem-plot hanya periode p terakhir, yang memperbolehkan kita untuk menjalankan long time simulations dan melihat hasil akhir tanpa plot yang berantakan dengan terlalu banyak periode. Syntax t[-m:] mem-plot elemen m terakhir dalam t (indeks negatif dalam hitungan susunan/daftar Pyhton dari akhir).
 +
 
 +
Kita bisa membandingkan metode Heun (atau setara metode RK2) dengan skema Euler-Crome:
 +
 
 +
compare(odespy_methods=[odespy.Heun, odespy.EulerCromer],
 +
omega=2, X_0=2, number_of_periods=20,
 +
time_intervals_per_period=20)
 +
 
 +
Gambar 4.22 menunjukkan bagaimana metode Heun (garis biru dengan piringan kecil) memiliki error yang cukup besar pada amplitude dan fase sesudah setelah periode 14-20 (kiri atas), namun menggunakan sebanyak tiga kali langkah waktu membuat kurvanya hampir sama (kanan atas). Akan tetapi setelah periode 194-200 error tersebut telah berkembang (kiri bawah), tetapi dapat cukup dikurangi dengan mengurangi separuh langkah waktu (kanan bawah).
 +
 
 +
Dengan semua metode di Odespy, sekarang menjadi mudah untuk mulai menjelajahi metode-metode lain, seperti perbedaan mundur (backward differences) bukannya perbedaan maju (forward differences) yang digunakan dalam skema Forward Euler. Latihan 4.17 mengatasi permasalahan tersebut.
 +
 
 +
Odespy berisi metode adaptif yang cukup canggih di mana pengguna "dijamin" untuk mendapatkan solusi dengan akurasi yang ditentukan. Tidak ada jaminan matematis, tetapi error untuk sebagian besar kasus tidak akan menyimpang secara signifikan dari toleransi pengguna yang mencerminkan keakuratan. Metode yang sangat populer dari jenis ini adalah metode Runge-Kutta-Fehlberg, yang menjalankan metode Runge-Kutta orde 4 dan menggunakan metode Runge-Kutta orde 5 untuk memperkirakan error sehingga dapat disesuaikan untuk menjaga error di bawah toleransi. Metode ini juga dikenal luas sebagai ode45, karena itulah nama fungsi yang mengimplementasikan metode ini di Matlab. Kita dapat dengan mudah menguji metode Runge-Kutta-Fehlberg segera setelah kita tahu nama Odespy yang sesuai, yaitu RKFehlberg:
 +
 
 +
compare(odespy_methods=[odespy.EulerCromer, odespy.RKFehlberg],
 +
omega=2, X_0=2, number_of_periods=200,
 +
time_intervals_per_period=40)
 +
 
 +
[[File:oscillating17-2.png]]
 +
 
 +
Perhatikan bahwa argumen time_intervals_per_period mengacu pada titik waktu di mana kami ingin solusinya. Poin-poin ini juga yang digunakan untuk perhitungan numerik dalam pemecah odespy.EulerCromer, sedangkan pemecah odespy.RKFehlberg akan menggunakan satu set titik waktu yang tidak diketahui karena interval waktu disesuaikan ketika metode berjalan. Orang dapat dengan mudah melihat titik-titik yang sebenarnya digunakan oleh metode karena ini tersedia sebagai himpunan solver.t_all (tetapi merencanakan atau memeriksa titik-titik membutuhkan modifikasi di dalam metode perbandingan).
 +
 
 +
Gambar 4.23 menunjukkan contoh komputasi di mana metode Runge-Kutta-Fehlberg jelas lebih unggul daripada skema Euler-Cromer dalam simulasi yang lama, tetapi perbandingannya tidak terlalu adil karena metode Runge-Kutta_Fehlberg berlaku sekitar dua kali lebih banyak langkah waktu dalam hal perhitungan ini dan melakukan lebih banyak pekerjaan per langkah waktu. Ini adalah tugas yang cukup rumit untuk membandingkan dua metode yang sangat berbeda dalam cara yang wajar sehingga pekerjaan komputasi versus akurasi dilaporkan secara ilmiah dengan baik.
 +
 
 +
[[File:oscillating18-2.png]]
 +
 
 +
==== 4.3.7 Metode Runge-Kutta Orde 4 ====
 
Metode Runge-Kutta Orde 4 adalah metode yang sering digunakan secara luas untuk menyelesaikan ODEs, karena menghasilkan data dengan tingkat akurasi yang tinggi bahkan dalam time step yang tidak terlalu kecil.
 
Metode Runge-Kutta Orde 4 adalah metode yang sering digunakan secara luas untuk menyelesaikan ODEs, karena menghasilkan data dengan tingkat akurasi yang tinggi bahkan dalam time step yang tidak terlalu kecil.
  
Gambar 4.24
+
[[File:1-.PNG]]
  
 
Algoritma; Pertama-tama kita nyatakan algoritma 4-stage
 
Algoritma; Pertama-tama kita nyatakan algoritma 4-stage
  
Rumus 4.59
+
[[File:2-.PNG]]
  
 
Dimana
 
Dimana
Rumus 4.60
 
Rumus 4.61
 
Rumus 4.62
 
  
Aplikasi; Kita bisa menjalankan simulasi seperti pada Figs. 4.16, 4.18, dan 4.21, untuk 40 periode. 10 periode terakhir ditunjukan melalui  Fig. 2.24. Hasil yang ditunjukan terlihat impresif sebagaimana penggunaan metode Euler-Cromer.
+
[[File:3-.PNG]]
 +
 
 +
[[File:4-.PNG]]
 +
 
 +
[[File:5-.PNG]]
  
Implementasi; Tingkatan dalam metode Runge-Kutta orde-4 bisa dengan mudah diimplementasikan sebagai modifikasi dari osc_Heun.py code. Sebagai alternatif, salah satu dapat menggunakan osc_odespy.py code dengan menyediakan argumen odespy_methods-[odespy.RK4] untuk membandingkan fungsi.
 
  
Derivasi; Derivasi dari metode Runge-Kutta orde-4 dapat disajikan dengan cara pedagogis yang menyatukan banyak elemen fundamental dari teknik diskritisasi numerik dan bisa menggambarkan banyak aspek “numerical thinking ”ketika membangun perkiraan metode solusi.
+
'''Aplikasi'''; Kita bisa menjalankan simulasi seperti pada Figs.  4.16, 4.18, dan 4.21, untuk 40 periode. 10 periode terakhir ditunjukan melalui  Fig. 2.24. Hasil yang ditunjukan terlihat impresif sebagaimana penggunaan metode Euler-Cromer.
  
Kita mulai dengan mengintegrasikan general ODE (gambar rumus1) dari waktu ke waktu, mulai dari (gambar rumus2) sampai (gambar rumus 3),
 
  
(Persamaan1)
+
'''Implementasi'''; Tingkatan dalam metode Runge-Kutta orde-4 bisa dengan mudah diimplementasikan sebagai modifikasi dari osc_Heun.py code. Sebagai alternatif, salah satu dapat menggunakan osc_odespy.py code dengan menyediakan argumen odespy_methods-[odespy.RK4] untuk membandingkan fungsi.
Tujuan dari komputasi (persamaan2), ketika (persamaan3) pada saat ini lebih dikenal dengan nilai (u tanda miring). Tantangan mengintegralkan muncul ketika integrand mengandung u yang tidak diketahuai antara tn sampai tn+1.
+
 
 +
 
 +
'''Derivasi'''; Derivasi dari metode Runge-Kutta orde-4 dapat disajikan dengan cara pedagogis yang menyatukan banyak elemen fundamental dari teknik diskritisasi numerik dan bisa menggambarkan banyak aspek “numerical thinking ”ketika membangun perkiraan metode solusi.
 +
 
 +
Kita mulai dengan mengintegrasikan general ODE [[File:6-.PNG]] dari waktu ke waktu, mulai dari tn sampai t(n_1),
 +
 
 +
[[File:9-.PNG]]
 +
 
 +
Tujuan dari komputasi [[File:10-.PNG]], ketika [[File:11-.PNG]] pada saat ini lebih dikenal dengan nilai ''u''. Tantangan mengintegralkan muncul ketika integrand mengandung ''u'' yang tidak diketahuai antara tn sampai t(n+1).
  
 
Integral tersebut dapat diperkirakan dengan menggunakan Simpson’s rule yang telah terkenal
 
Integral tersebut dapat diperkirakan dengan menggunakan Simpson’s rule yang telah terkenal
(Persamaan 4)
+
 
Permasalahan dengan persamaan ini adalah kita tidak mengetahui nilai dari (persamaan5) dan (pers6) karena hanya u^n yang tersedia dan hanya f^n yang dapat dihitung.
+
[[File:12-.PNG]]
Untuk melanjutkana, idenya dalah menggunakan berbagai perkiraan untuk (persamaan7) dan (pers8) berdasarkan penggunaan skema yang telah diketahui untuk ODE dalam interval (pers9) dan (pers10). Mari kita bagi persamaan integral menjadi empat suku.
+
 
persamaan
+
Permasalahan dengan persamaan ini adalah kita tidak mengetahui nilai dari [[File:13-.PNG]] dan [[File:14-.PNG]] karena hanya u^n yang tersedia dan hanya f^n yang dapat dihitung.
 +
 
 +
Untuk melanjutkan, idenya dalah menggunakan berbagai perkiraan untuk [[File:15-.PNG]] dan [[File:16-.PNG]] berdasarkan penggunaan skema yang telah diketahui untuk ODE dalam interval [[File:17-.PNG]] dan [[File:18-.PNG]]. Mari kita bagi persamaan integral menjadi empat suku.
 +
 
 +
 
 +
[[File:19-.PNG]]
 +
 
 +
Dimana [[File:C01.JPG|40px]], [[File:C02.JPG|40px]], dan [[File:C03.JPG|40px]] adalah pendekatan untuk [[File:C04.JPG|40px]] dan [[File:C05.JPG|40px]] yang dapat digunakan pada perhitungan. Untuk [[File:C01.JPG|40px]] dapat menggunakan pendekatan untuk [[File:C06.JPG|40px]] berdasarkan tahap Forward Euler pada size [[File:C14(2).JPG|27px]]
 +
 
 +
 
 +
[[File:4-63.JPG|400px]]
 +
 
 +
 
 +
Persamaan ini mempermudah prediksi [[File:C04.JPG|40px]], sehingga untuk [[File:C02.JPG|40px]] kita dapat mencoba metode Backward Euler untuk memperkirakan [[File:C06.JPG|40px]]
 +
 
 +
 
 +
[[File:4-64.JPG|400px]]
 +
 
 +
 
 +
Dengan [[File:C02.JPG|40px]] sebagai pendekatan untuk [[File:C04.JPG|40px]], pada akhirnya bentuk akhir dari [[File:C03.JPG|40px]] dapat menggunakan metode midpoint (atau central difference, juga disebut metode Crank-Nicholson) untuk memperkirakan [[File:C15(2).JPG|30px]].
 +
 
 +
 
 +
[[File:4-65.JPG|400px]]
 +
 
 +
 
 +
Kita telah menggunakan metode Forward dan Backward Euler, juga centered difference approximation pada konteks Simpsons rule. Diharapkan kombinasi dari metode ini dapat menghasilkan overall time stepping dari [[File:C07.JPG|20px]] ke [[File:C08.JPG|40px]] yang lebih akurat dibandingkan individual steps (yang memiliki error proportional dengan [[File:C09.JPG|20px]] dan [[File:C10(2).JPG|25px]]). Hal ini benar bahwa: error numerik yang terjadi seperti [[File:C11(2).JPG|40px]] Untuk konstanta ''C'', artinya error lebih cepat mendekati nol ketika time step size dikurangi, dibandingkan  dengan metode Forward Euler [[File:C12(2).JPG|80px]], metode Euler-Cromer [[File:C12(2).JPG|80px]],atau Runge Kutta orde 2, atau metode Heuns [[File:C13(2).JPG|80px]].
 +
 
 +
Perhatikan bahwa Metode Runge-Kutta Orde 4 sepenuhnya eksplisit jadi tidak diperlukan untuk menyelesaikannya dengan persamaan aljabar baik secara linier maupun non linier, terlepas dari apa yang terlihat pada ''f''. Namun nilai kesetabilannya kondisional dan bergantung pada nilai ''f'' tersebut. Ada sebuah bagian besar  dari metode implisit Runge-Kutta yang nilai kesetabilannya tidak kondisional. namun diperlukan solusi dari persamaan aljabar yang melibatkan nilai ''f'' pada setiap "''time step''". Odespy dapat dimanfaatkan untuk mendukung penyelesaian dari banyak metode Runge-Katta yang eksplisit. Tetapi belum bisa digunakan untuk metode Runge-Katta yang implisit.
 +
 
 +
==== 4.3.8 Efek Lain : Damping, Nonlinearity, dan external force ====
 +
 
 +
Model permasalahan u’’ + ω2u = 0 adalah model matematika yang paling simple untuk oscilating system. Namun, Model ini lebih banyak membutuhkan metode numerik, seperti yang sudah kita lihat, dan sangat berguna untuk menjadi tolak ukur untuk mengevaluasi kinerja dari metode numerik.
 +
 
 +
Dalam Pengaplikasian dikehidupan nyata lebih banyak melibatkan efek fisika, yang mengarahkan ke persamaan diferensial dengan ketentuan yang lebih banyak dan juga lebih kompleks. biasanya, memiliki kekuatan redaman f (u ') dan pegas s (u). Kedua gaya ini tergantung pada nonlinear dari uraiannya, u’ atau u. sebagai tambahan, gaya lingkungan F(t) jufga bekerja pada sistem. Contohnya, pendulum klasik memiliki “pegas” nonlinear atau mengembalikan gaya s(u) ~ sin (u), dan gaya tahan dari udara pada pendulum menyebabkan terjadinya gaya redam f(u’) ~ |u’|u. Contoh dari gaya lingkungan adalah getaran dari tanah (seperti gempa) dan juga seperti ombak atau angina.
 +
 
 +
Dengan tiga jenis gaya yang bekerja pada sistem : F(t), f(u’), dan s(u). maka dapat ditulis persamaan F(t) – f(u’) – s(u). Tanda mines didepan f dan s menunjukan bahwa fungsi ini didefinisikan sebagai gaya yang melawan gerakan. Sebagai Contoh, Pegas yang terpasang pada roda mobil dikombinasikan dengan beberapa perdeam yang efektif. Masing-masing memiliki gaya redam f(u’) = bu’ yang bekerja melawan kecepatan pegas u’. gaya fisika yang sesuai dapat dtulis –f: -bu’, yang menunjuk ke bawah saat pegas diregangkan (dan poin u’ ke atas), sedangkan -f bertindak ke atas saat pegas dikompresi (dan poin u’ ke bawah).
 +
 
 +
Gambar 4.25 menunjukan contoh dari massa m terpasang dengan pegas nonlinear dan dashpot, dan bersubyek pada gaya lingkungan F(t). Namun, model umum yang kita miliki dapat juga digunakan pada pendulum pada gambar 4.26 dengan s (u) = m g sin θ dan f (u ̇) = 1/2  C_D Aϱ(θ|) ̇θ| (Dimana CD = 0.4, A adalah area perpotongan dari body dan ϱ adalah densitas udara)
 +
 
 +
[[File:Gambar425.png]]
 +
 
 +
Gambar 4.25 Sistem Oscillating General
 +
 
 +
Hukum Newton kedua untuk sistem yang dapat ditulis dengan akselerasi waktu massa pada sisi kiri dan gaya pada sisi kanan:
 +
 
 +
[[File:438rumus1.png]]
 +
 
 +
Bagaimanapun persamaan ini lebih umum disusun ulang menjadi
 +
 
 +
[[File:438rumus2.png]]
 +
 
 +
Karena persamaan diferensial adalah orde 2, disebabkan oleh istilah u^'', kita membutuhkan dua kondisi awal:
 +
 
 +
[[File:438rumus3.png]]
 +
 
 +
[[File:gambar426.png]]
 +
 
 +
Gambar 4.26 Sebuah pendulum dengan gaya
 +
 
 +
Catat bahwa dengan pilihan [[File:438rumus4.png]] kita memperoleh kembali persamaan diferensial biasa [[File:438rumus5.png]]
 +
 
 +
Bagaimana kita bisa menyelesaikan (4.66)? sebagaimana persamaan diferensial biasa yang simpel [[File:438rumus6.png]] kita mulai dengan menulis ulang persamaan diferensial biasa orde 2 sebagai sebuah sistem dari dua persamaan diferensial biasa orde 1:
 +
 
 +
[[File:438rumus8.png]]
 +
 +
Kondisi awal menjadi
 +
 
 +
[[File:438rumus9.png]]
 +
 
 +
Setiap metode dari sebuah sistem persamaan diferensial biasa orde 1 dapat digunakan untuk menyelesaikan [[File:438rumus10.png]]
 +
 
 +
'''The Euler-Cromer scheme'''
 +
 
 +
Sebuah pilihan atraktif dari sebuah implementasi, akurasi dan efisiensi sudut pandang adalah skema Euler-Cromer dimana kita mengambil sebuah perbedaan kedepan pada (4.68) dan perbedaan kebelakang pada (4.69):
 +
 
 +
[[File:438rumus11.png]]
 +
 
 +
Kita dapat dengan mudah menyelesaikan [[File:438rumus12.png]] yang tidak diketahui:
 +
 
 +
[[File:438rumus13.png]]
 +
 
 +
 
 +
'''kata kata dalam perintah ODEs'''
 +
 +
Perintah ODE dalam sistem ODE penting untuk model yang diperluas (4.68) - (4.69). Bayangkan kita menulis persamaan untuk u’ terlebih dahulu dan kemudian untuk v’. Metode Euler-Cromer akan menggunakan  forward difference untuk u^n+1 dan kemudian backward difference akan menggunakannya untuk v^n+1. Yang Terkhir akan menyebabkan persamaan nonlinear algebraic untuk v^n+1
 +
 
 +
[[File:(4.3.8) 1.png]]
 +
 
 +
jika f(v) adalah fungsi nonlinear dari v Ini akan membutuhkan metode numerik untuk persamaan aljabar nonlinier untuk mencari v^n+1 saat memperbarui v^n+1 melalui forward difference memberikan persamaan untuk v^n+1 itu linear dan sepele untuk dipecahkan dengan tangan.
 +
 
 +
File osc_EC_general.pymemiliki fungsi Euler Cromer yang mengimplementasikan metode ini:
 +
 
 +
[[File:(4.3.8) 2.png]]
 +
 
 +
[[File:(4.3.8) 3.png]]
 +
 
 +
Metode Runge-Kutta orde ke 4
 +
 
 +
Metode RK4 hanya mengevaluasi sisi kanan  sistem ODE,
 +
 
 +
[[File:(4.3.8) 4.png]]
 +
 
 +
untuk nilai-nilai u, v, dan t yang diketahui, maka metode ini sangat sederhana untuk digunakan terlepas dari bagaimana fungsi s(u) dan f(v)dipilih.
 +
 +
 
 +
 
  
  
 
==== 4.3.9 ilustrasi redaman linier ====
 
==== 4.3.9 ilustrasi redaman linier ====
  
Kami menganggap sistem rekayasa dengan pegas linier, s (u) = kx, dan peredam kental, di mana gaya peredaman adalah porpotional terhadap u ', f (u') = bu ', untuk beberapa konstanta b> 0. Pilihan ini dapat memodelkan sistem pegas vertikal di dalam mobil (tetapi insinyur sering suka menggambarkan sistem tersebut dengan massa bergerak horizontal seperti yang digambarkan pada Gambar 4.25). kita dapat memilih nilai-nilai sederhana untuk konstanta untuk mengilustrasikan efek dasar redaman (dan kegembiraan selanjutnya). Memilih osilasi sebagai fungsi u (t) = cos t sederhana dalam kasus undamped, kita dapat menetapkan m = 1, k = 1, b = 0,3, Uo = 1, Vo = 0. Fungsi berikut mengimplementasikan kasus ini:
+
Kami menganggap sistem rekayasa dengan pegas linier, s(u) = kx, dan peredam kental, di mana gaya peredaman adalah porpotional terhadap u', f(u') = bu', untuk beberapa konstanta b > 0. Pilihan ini dapat memodelkan sistem pegas vertikal di dalam mobil (tetapi insinyur sering suka menggambarkan sistem tersebut dengan massa bergerak horizontal seperti yang digambarkan pada Gambar 4.25). kita dapat memilih nilai-nilai sederhana untuk konstanta untuk mengilustrasikan efek dasar redaman (dan kegembiraan selanjutnya). Memilih osilasi sebagai fungsi u(t) = cos t sederhana dalam kasus tak teredam, kita dapat menetapkan m = 1, k = 1, b = 0,3, Uo = 1, Vo = 0. Fungsi berikut mengimplementasikan kasus ini:
  
 +
[[File:Wafi_439-1.png|500px]]
  
 +
Fungsi plot_u adalah kumpulan plot untuk merencanakan u(t), atau bagian darinya. Gambar 4.27 menunjukkan efek dari bu': kita memiliki osilasi dengan (perkiraan) periode 2π, seperti yang diharapkan, tetapi amplitudo teredam secara efisien.
  
Fungsi plot_u adalah kumpulan pernyataan plot untuk merencanakan u (t), atau bagian darinya. Gambar 4.27 menunjukkan efek dari istilah bu: kita memiliki osilasi dengan (perkiraan) periode 2π, seperti yang diharapkan, tetapi amplitudo teredam secara efisien.
+
[[File:Kania 439-2.png|500px]]
  
  
Ceritakan tentang kebingungan dengan masalah berskala
+
'''Komentar mengenai pekerjaan dengan masalah berskala'''
Daripada menetapkan b = 0,3 dan m - k = Uo = 1 sebagai nilai fisik yang cukup "tidak mungkin", akan lebih baik untuk skala persamaan mu '' + bu '+ ku = 0. Ini berarti bahwa kita memperkenalkan independen tanpa tergantung dan tidak tergantung variabel:
 
  
 +
Alih-alih menetapkan b = 0,3 dan m = k = Uo = 1 sebagai nilai fisik yang “tidak mungkin”, akan lebih baik untuk skala persamaan mu" + bu' + ku = 0. Ini mengartikan bahwa kita memasukan variabel independen dan dependen yang tak berdimensi :
  
 +
[[File:Kania_439-3.png|200px]]
  
Di mana tc dan uc adalah ukuran karakteristik waktu dan perpindahan, masing-masing, sehingga ¯t dan ¯u memiliki ukuran tipikal mereka di sekitar kesatuan. Dalam masalah ini, kita dapat memilih uc = Uo dan tc = √ (m / k). Ini memberikan masalah skala (atau tanpa dimensi) berikut untuk kuantitas tak berdimensi ¯u (¯t):
+
Di mana tc dan uc adalah ukuran karakteristik waktu dan perpindahan, sehingga [[File:Kania_439-5.png|15px]] dan [[File:Kania_439-6.png|20px]] memiliki ukuran tipikal mereka didekat kesatuan. Dalam masalah ini, kita dapat memilih [[File:Kania_439-7.png|70px]] dan [[File:Kania_439-8.png|80px]]. Ini memberikan masalah yang berskala (atau tanpa dimensi) berikut untuk kuantitas tak berdimensi [[File:Kania_439-9.png|40px]]:
  
 +
[[File:Kania_439-4.png|600px]]
  
 +
Faktnya adalah hanya ada satu parameter fisik di kasus ini: angka β. Menyelesaikan masalah ini begitu juga terkait dengan masalah utama dengan parameter yaitu m = k = Uo = 1 dan b = β. Tetapi untuk menyelsaikan masalah dengan satuan lebih umum: jika kita memdapatkan solusi ¯u(¯t;β), kita dapat menemukan solusi fisik pada kasus ini, dikarenakan :
  
Fakta yang mencolok adalah bahwa hanya ada satu parameter fisik dalam masalah ini: angka β yang tidak disebutkan. Memecahkan masalah ini terkait dengan penyelesaian masalah asli (dengan dimensi) dengan parameter m = k = Uo = 1 dan b = β. Namun, menyelesaikan masalah tak berdimensi lebih umum: jika kita memiliki solusiu (¯t; β), kita dapat menemukan solusi fisik dari berbagai masalah karena
+
[[File:439rumus.png|200px]]
  
 +
Selama β didapat, kita dapat menemukan u untuk Uo , k, dan m dengan rumus diatas, dengan begitu pengerjaan simulasi dapat dipersingkat waktu. Ini menunjukkan pengerjaan dengan skala atau masalah satuan.
  
 
Selama β diperbaiki, kita dapat menemukan u untuk Uo, k dan m dari rumus di atas! Dengan cara ini, simulasi yang memakan waktu dapat dilakukan hanya sekali, tetapi masih memberikan banyak solusi. Ini menunjukkan kekuatan bekerja dengan masalah skala atau tanpa dimensi.
 
  
 
==== 4.3.10. Ilustrasi Redaman Linier Dengan Eksitasi Sinusoidal ====
 
==== 4.3.10. Ilustrasi Redaman Linier Dengan Eksitasi Sinusoidal ====
 
Sekarang kita akan memperluas contoh sebelumnya untuk menambah beberapa gaya osilasi eksternal pada sistem: F (t) = Asin (wt). Mengendarai mobil di jalan dengan lonjakan sinusoidal mungkin memberikan eksitasi eksternal pada sistem pegas di mobil (w terkait dengan kecepatan mobil).
 
Sekarang kita akan memperluas contoh sebelumnya untuk menambah beberapa gaya osilasi eksternal pada sistem: F (t) = Asin (wt). Mengendarai mobil di jalan dengan lonjakan sinusoidal mungkin memberikan eksitasi eksternal pada sistem pegas di mobil (w terkait dengan kecepatan mobil).
 +
 +
    from math import pi,sin
 +
    w = 3
 +
    A = 0.5
 +
    F = lambda t: A*sin(w*t)
 +
 +
kita dapatkan grafik pada gambar 4.28 .Perbedaan yang mencolok dari Gambar 4.27 adalah bahwa osilasi dimulai sebagai sinyal ''cos t'' teredam tanpa banyak pengaruh gaya eksternal, tetapi kemudian osilasi bebas dari sistem yang tidak teredam ''(cos t) u’’ + u =  0'' mati dan gaya eksternal ''0: 5 sin.(3t)''  menimbulkan osilasi dengan periode yang lebih pendek ''2phi/3''. Dianjurkan untuk menggunakan beberapa nilai A yang lebih besar dan beralih dari sinus ke acosinus dalam F dan mengamati efeknya. Jika mencarinya di dalam buku fisika, Anda dapat menemukan solusi analitik yang tepat untuk masalah persamaan diferensial dalam kasus ini.
 +
 +
====4.3.11. Sistem pegas-massa dengan gesekan luncur====
 +
 +
Sebuah benda dengan massa ''m'' bekerja pada sebuah pegas dengan kekakuan ''k'' saat meluncur pada sebuah bidang permukaan. Benda tersebut mengalami gaya gesek ''f(u')'' disebabkan terjadi kontak antara benda dengan bidang permukaan seperti terlihat pada Gambar 4.30. Gaya gesek ''f(u')'' dapat dimodelkan dengan gesekan Coulomb sebagai berikut:
 +
 +
[[File:Eq4.3.11.1.png|180px|center]]
 +
 +
Dimana ''μ'' adalah koefisien gesek, dan mg merupakan gaya normal pada bidang permukaan benda yang bergerak. Formula ini dapat juga ditulis sebagai ''f(u') = μmg sign (u')'', dengan syarat fungsi signum sign (x) didefinisikan nol untuk ''x'' = 0 (numpy.sign mempunyai sifat ini). Untuk memastikan bahwa signum dari definisi ''f'' benar, ingat bahwa gaya fisis aktual adalah ''-f'' dan positif (misal ''f''<0) ketika gaya tersebut bekerja berlawanan dengan benda yang bergerak dengan kecepatan ''u'''<0.
 +
 +
[[File:606px-1d_oscillating_dynamic_system_29.1.png|600px|thumb|center|Gambar 4.30 Sketsa dari sebuah subjek sistem osilasi dinamis untuk gesekan luncur dan gaya pegas satu dimensi.]]
 +
 +
Gaya pegas nonlinear diambil sebagai:
 +
 +
[[File:Eq4.3.11.2.png|160px|center]]
 +
 +
Yang mana nilai ''-ku'' diperkirakan untuk nilai ''u'' yang kecil, namun stabil pada ±''k/α'' untuk nilai ±''αu'' yang besar. Berikut adalah plot dengan ''k''=1000 dan ''u'' ∈ [-0.1,0.1] untuk tiga nilai ''α'':
 +
 +
[[File:591px-1d_oscillating_dynamic_system_30.1.png|center]]
 +
 +
Jika tidak ada gaya eksitasi eksternal yang bekerja pada benda, maka persamaan gerak yang kita dapatkan adalah:
 +
 +
[[File:Eq4.3.11.3.png|300px|center]]
 +
 +
Mari kita simulasikan situasi dimana sebuah benda dengan massa 1 kg meluncur pada bidang permukaan dengan ''μ'' = 0.4, terikat pada pegas dengan kekakuan ''k'' = 1000 kg/s^2. Perpindahan awal benda adalah 10 cm, dan parameter ''α'' dalam ''s(u)'' diatur pada 60 1/m.
 +
 +
Dengan menggunakan fungsi EulerCromer dari kode osc_EC_general, kita dapat menulis fungsi sliding_friction untuk menyelesaikan masalah ini:
 +
 +
Def    sliding_friction():
 +
        from numpy import tanh, sign
 +
       
 +
        f = lambda v: mu*m*g*sign(v)
 +
        alpha = 60.0
 +
        s = lambda u: k/alpha*tanh(alpha*u)
 +
        F = lambda t: 0
 +
       
 +
        g = 9.81
 +
        mu = 0.4
 +
        m = 1
 +
        k = 1000
 +
       
 +
        U_0 = 0.1
 +
        V_0 = 0
 +
       
 +
        T = 2
 +
        dt = T/5000
 +
       
 +
        u, v, t = EulerCromer(f=f, s=s, F=F, m=m, T=T,
 +
                              U_0=U_0, V_0=V_0, dt=dt)
 +
        plot_u(u, t)
 +
 +
Setelah menjalankan fungsi sliding_friction memberi kita hasil seperti pada Gambar. 4.31 dengan ''s(u)= -k/α tanh⁡(αu)'', (kiri) dan versi linierisasi ''s(u)=ku'' (kanan).
 +
 +
[[File:Photoagh.png|600px|thumb|center|Gambar 4.31 Efek pegas nonlinear (kiri) dan linier (kanan) pada gesekan luncur]]
 +
 +
====4.3.12. Metode finite diference; Undamped, Linear Case====
 +
 +
Selanjutnya kita akan membahas metode numerik untuk ODE orde kedua
 +
 +
  u^''+ω^2 u=0,  u(0)=U_0,u^' (0)=0,t∈(0,T]
 +
 +
tanpa menulis ulang ODE sebagai sistem ODE orde pertama. Motivasi utama untuk "metode solusi lain" adalah bahwa prinsip-prinsip diskritisasi menghasilkan skema yang sangat baik, dan yang lebih penting, pemikiran seputar diskritisasi bisa digunakan kembali ketika memecahkan persamaan diferensial parsial.
 +
Gagasan utama dari metode numerik ini adalah untuk memperkirakan urutan kedua turunan u'' dengan selisih terbatas. Sementara ada beberapa pilihan perbedaan perkiraan untuk derivatif orde pertama, ada satu rumus yang mendominasi untuk turunan orde kedua:
 +
 +
[[File:Persamaan4.74.jpg]]
 +
 +
Error dalam perkiraan tersebut proporsional terhadap ∆t^2. Biarkan ODE valid di beberapa titik waktu yang berubah ubah t_n,
 +
 +
  u^'' (t_n )+ ω^2 u (t_n )=0
 +
 +
Selanjutnya memasukkan rumus perkiraan  (4.74) diatas, sehingga di dapatkan
 +
 +
[[File:Persamaan4.75.jpg]]
 +
 +
Sekarang diasumsikan bahwa u^(n-1) dan u^n sudah dihitung, dan u^(n+1) adalah yang baru
 +
tidak diketahui. Memecahkan sehubungan dengan u^(n+1)
 +
 +
[[File:Persamaan4.76.jpg]]
 +
 +
Masalah besar muncul ketika kita ingin memulai skema. Kita tahu bahwa u^0 = U_0, tetapi menerapkan (4,76) untuk n=0 untuk menghitung u^1
 +
 +
[[File:Persamaan4.77.jpg]]
 +
 +
Dimana kita tidak mengetahui U-1. Kondisi awal U’ (0) = 0 dapat membantu kiti untuk menghilangkan U-1 dan kondisi ini bagaimanapun juga harus dimasukkan dalam beberapa cara. Untuk tujuan ini, kami mendiskritasikan u’(0) = 0 dengan perbedaan terpusat,
 +
 +
 +
[[File:Persamaan4.78.jpg]]
 +
 +
Oleh karena itu, u-1 = u1, dan kita dapat menggunakan relasi ini untuk menghilangkan u1 di persamaan 4.77
 +
 +
[[File:Persamaan4.79.png]]
 +
 +
Dengan U0 = U0 dan u1 dihitung dari persamaan (4.78), kita dapat menghitung u2, u3, dan seterusnya dari persamaan (4.76). Latihan 4.19 meminta Anda untuk mengeksplorasi bagaimana langkah-langkah di atas diubah seandainya  kita memiliki kondisi awal bukan nol u’ (0) = V0
 +
 +
Kita dapat memperkirakan kondisi awal U’(0) dengan menggunakan Forward difference
 +
 +
[[File:Persamaan4.80.jpg]]
 +
 +
Mengarah pada u1 = u0 . lalu kita dapat menggunakan persamaan (4.76) untuk langkah selanjutnya . Walaupun forward difference memiliki kesalahan proporsional ke ∆t . dimana centered difference yang kita gunakan memiliki error proporsional ke ∆t2. Yang dimana  kompatibel dengan akurasi (erro yang enunjukan ∆t2) yang digunakan dalam diskritisasi persamaan diferensial.
 +
Metode untuk ODE orde kedua yang dijelaskan di atas berjalan di bawah nama metode Störmer atau integrasi Verlet 7. Ternyata metode ini secara matematis setara dengan skema Euler-Cromer  Atau lebih tepatnya, rumus umum (4,76) setara dengan rumus Euler-Cromer.
 +
  
 
====4.3.13 Metode finite diference; damping linier====
 
====4.3.13 Metode finite diference; damping linier====
Line 250: Line 723:
 
[[File:4.80a.png]]
 
[[File:4.80a.png]]
  
dan memasukkan perkiraan perbedaan terhingga pada u" dan u' hasil dalam
+
Dan memasukkan perkiraan perbedaan terhingga pada u" dan u' hasil dalam
  
  
Line 256: Line 729:
  
  
dimana F" adalah notasi pendek untuk F(t). Persamaan (4.81) adalah linear dalam
+
Dimana F" adalah notasi pendek untuk F(t). Persamaan (4.81) adalah linear dalam
 
u^(n+1) tak diketahui kita dapat dengan mudah memecahkan untuk kuantitas ini:
 
u^(n+1) tak diketahui kita dapat dengan mudah memecahkan untuk kuantitas ini:
  
Line 262: Line 735:
 
[[File:4.82.png]]
 
[[File:4.82.png]]
  
== Artikel 1 Hasil diskusi : judul .. ==
 
  
== Artikel 2 Hasil diskusi : judul ..===
+
Pada kasus tanpa redaman, kita membutuhkan formula khusus untuk u1. kondisi awal U`(0) = 0 menyatakan bahwa u-1 = u1, dan dengan persamaan (4.82) untuk n = 0, kita mendapatkan.
 +
 
 +
[[File:4.8.3casees.JPG]]
 +
 
 +
Pada kasus yang lebih unun dengan sebuah bentuk redaman nonlinier f(u`),
 +
 
 +
 
 +
[[File:4.8.3.2.1.2.JPG]]
 +
 
 +
 
 +
Kita mendapatkan
 +
 
 +
[[File:4.8.3.2.1.2.1.JPG]]
 +
 
 +
Dimana sebauh persamaan ajabar non linier untuk un+1 bahwa harus diseleseikan dengan metode numerik. Skema lebih bagus diperoleh dari penggunaan "backward difference" untuk u`,
 +
 
 +
[[File:4.8.3.2.1.2.1.2.JPG]]
 +
 
 +
Karena pada bagian redaman akan lebih diketahui, yang hanya melibatkan un dan un-1, dan kita dapat dengan mudah menyelesaikan untuk un+1.
 +
Kelemahan dari backward difference dibandingkan dengan centered difference (4.80) adalah ini mengurangi urutan akurasi dalam skema keseluruhan dari ∆t2 ke ∆t. Pada kenyataanya, skema Euler-Cromer mengevaluasi istilah redaman nonlinear sebagai f(vn), saat menghitung vn+1, dan ini setara dengan menggunakan backward difference di atas. Akibatnya, kenyamanan skema Euler-Cromer untuk redaman nonlinier datang dengan konsekuensi menurunkan akurasi keseluruhan skema dari urutan kedua ke urutan pertama pada ∆t. Menggunakan trik yang sama dalam skema beda hingga {finite difference} untuk persamaan diferensial orde kedua, yaitu, menggunakan backward difference dalam f(u’), membuat skema ini sama bagus dan akuratnya seperti skema Euler-Cromer pada kasus nonlinier umum mu”+f(u’)+s(u) = F.
 +
 
 +
 
 +
=='''Hasil tugas kaloborasi Oscillating 1-D Dynamic Artikel 1'''==
 +
 
 +
=== Pembagian Tema ===
 +
 
 +
4.3.1 [http://air.eng.ui.ac.id/index.php?title=Oldy_Fahlovi Oldy Fahlovvi], [http://air.eng.ui.ac.id/index.php?title=Muchalis_Zikramansyah_Masuku Muchalis Zikramansyah Masuku], [http://air.eng.ui.ac.id/index.php?title=AHMAD_ZIKRI Ahmad Zikri], [http://air.eng.ui.ac.id/index.php?title=Muhammad_Irfan_Dzaky Muhammad Irfan Dzaky]
 +
 
 +
4.3.2 Andhika, faturahman
 +
 
 +
4.3.4 iqbal & Alghi & Adam, aji suryadi
 +
 
 +
4.3.5 Shabrina & Edo, jerry, raihan
 +
 
 +
4.3.6 ronald & Desy & yophie, ardi
 +
 
 +
4.3.7 Kania & Chandra, evi & Dieter
 +
 
 +
4.3.9 Wafirul & fajri & keni, maha
 +
 
 +
4.3.10 [http://air.eng.ui.ac.id/index.php?title=Daniel_Meino_Soedira Daniel Meino Soedira] & [http://air.eng.ui.ac.id/index.php?title=Paskal_Rachman Paskal Rachman],
 +
[http://air.eng.ui.ac.id/index.php?title=Joko.triwardono Joko Triwardono] & [http://air.eng.ui.ac.id/index.php?title=Aghnia_Ilmiah_Nurhudan Aghnia Ilmiah Nurhudan]
 +
 
 +
4.3.11 Bambang ali.
 +
 
 +
4.3.12 Bagus, maheka, adzanna, Adinda
 +
 
 +
4.3.13 Harry, wisnu, Ichwan, fadli
 +
 
 +
=== Artikel Kolaborasi : ''1-D OSCILLATING SYSTEM'' arranged by [http://air.eng.ui.ac.id/index.php?title=Oldy_Fahlovi Oldy Fahlovvi], [http://air.eng.ui.ac.id/index.php?title=Muchalis_Zikramansyah_Masuku Muchalis Zikramansyah Masuku], [http://air.eng.ui.ac.id/index.php?title=AHMAD_ZIKRI Ahmad Zikri], [http://air.eng.ui.ac.id/index.php?title=Muhammad_Irfan_Dzaky Muhammad Irfan Dzaky]===
 +
 
 +
Berikut ini kami lampirkan tugas kolaborasi tentang ''1-D OSCILLATING SYSTEM'' dalam bentuk slideshow.
 +
<gallery mode="slideshow">
 +
File:Artikel Kolaborasi - Komputasi Teknik-01.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-02.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-03.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-04.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-05.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-06.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-07.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-08.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-09.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-10.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-11.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-12.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-13.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-14.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-15.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-16.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-17.jpg
 +
File:Artikel Kolaborasi - Komputasi Teknik-18.jpg
 +
</gallery>
 +
 
 +
=== Artikel 1 Hasil diskusi : OSCILLATING 1-D DYNAMIC SYSTEM with 1 MASS, 3 SPRING and 1 DAMPING ===
 +
 
 +
<gallery mode="slideshow">
 +
File:ANN_matlab_evi_001.png
 +
File:ANN_matlab_evi_002.png
 +
File:ANN_matlab_evi_003.png
 +
File:ANN_matlab_evi_004.png
 +
File:ANN_matlab_evi_005.png
 +
File:ANN_matlab_evi_006.png
 +
File:ANN_matlab_evi_007.png
 +
File:ANN_matlab_evi_008.png
 +
File:ANN_matlab_evi_014.png
 +
File:ANN_matlab_evi_015.png
 +
File:ANN_matlab_evi_011.png
 +
File:ANN_matlab_evi_012.png
 +
File:ANN_matlab_evi3.png
 +
</gallery>
 +
 
 +
=== Tema 4.3.13  Linear Damping Oscillation ===
 +
<gallery mode="slideshow">
 +
File:Wisnu 12346798.png
 +
File:Wisnu 123467989.png
 +
File:Wisnu 12346798910.png
 +
File:Wisnu 12346798435435.png
 +
File:Wisnu 1234679843fdsaf4.png
 +
File:Wisnu 12346798fdsaf4.png
 +
File:qwerwqerqwerq_caseII_1.png
 +
File:qwerwqerqwerq_caseII_2.png
 +
File:24-04-2020-1-tugas komtek.png
 +
File:2020-04-24 23 12 57-Spyder (Python 3.7).png
 +
File:2020-04-24 23 13 22-Spyder (Python 3.7).png
 +
File:2020-04-24 23 13 53-Spyder (Python 3.7).png
 +
File:2020-04-24 23 13 53-Spyder (Python 3.7).png
 +
File:hasil-24-04-2020.png
 +
File:2020-04-24 23 47 29-Book1 - Excel.png
 +
</gallery>
 +
 
 +
=== Artikel Kolaborasi - Sistem Osilasi Satu Dimensi ===
 +
 
 +
Berikut adalah tugas Artikel Kolaborasi kelompok kami mengenai ''Penggunaan Perangkat Lunak Python untuk Menyelesaikan ODEs pada Sistem Mekanik Berosilasi'' dengan anggota sebagai berikut:
 +
 
 +
1. Ardy Lefran Lololau
 +
 
 +
2. I Gusti Agung Ayu Desy Wulandari
 +
 
 +
3. Ronald Akbar
 +
 
 +
4. Yophie Dikaimana
 +
 
 +
<gallery mode="slideshow">
 +
File:Artikel_4.3.6_1.jpg
 +
File:Artikel_4.3.6_2.jpg
 +
File:Artikel_4.3.6_3.jpg
 +
File:Artikel_4.3.6_4.jpg
 +
File:Artikel_4.3.6_5.jpg
 +
File:Artikel_4.3.6_6.jpg
 +
File:Artikel_4.3.6_7.jpg
 +
File:Artikel_4.3.6_8.jpg
 +
File:Artikel_4.3.6_9.jpg
 +
</gallery>
 +
 
 +
Selanjutnya berikut ini adalah hasil diskusi kelompok kami mengenai sistem osilasi satu dimensi dengan menggunakan metode Artificial Neural Network (ANN)
 +
 
 +
<gallery mode="slideshow">
 +
File:ANN_final_1.jpg
 +
File:ANN_final_2.jpg
 +
File:ANN_final_3.jpg
 +
File:ANN_final_4.jpg
 +
File:ANN_final_5.jpg
 +
File:ANN_final_6.jpg
 +
File:ANN_final_7.jpg
 +
File:ANN_final_8.jpg
 +
File:ANN_final_9.jpg
 +
File:ANN_final_10.jpg
 +
</gallery>
 +
 
 +
=== Tugas kolaborasi 4.3.4. ARTIKEL OSCILLATING ORDE 1 DUMPING PADA SISTEM SPRING MOBIL ===
 +
 
 +
<gallery mode="slideshow">
 +
 
 +
File:AK_1.jpg
 +
File:AK_2.jpg
 +
File:AK_3.jpg
 +
File:AK_4.jpg
 +
File:AK_5.jpg
 +
File:AK_6.jpg
 +
File:AK_7.jpg
 +
File:AK_8.jpg
 +
File:AK_9.jpg
 +
 
 +
</gallery>
 +
 
 +
=== Tugas kolaborasi 4.3.5 sistem osilasi satu dimensi runge - kutta : studi kasus shock breaker motor ===
 +
[[File:Collab13 (5).jpg]]
 +
[[File:Collab13 (6).jpg]]
 +
[[File:Collab13 (7).jpg]]
 +
[[File:Collab13 (8).jpg]]
 +
[[File:Collab13 (9).jpg]]
 +
[[File:Collab13 (10).jpg]]
 +
[[File:Collab13 (11).jpg]]
 +
[[File:Collab13 (12).jpg]]
 +
[[File:Collab13 (13).jpg]]
 +
[[File:Collab13 (14).jpg]]
 +
[[File:Collab13 (15).jpg]]
 +
[[File:Collab13 (16).jpg]]
 +
[[File:Collab13 (1).jpg]]
 +
[[File:Collab13 (2).jpg]]
 +
[[File:Collab13 (3).jpg]]
 +
[[File:Collab13 (4).jpg]]
 +
 
 +
 
 +
=== Artikel Tugas 4.3.9 Illustrasi redaman linear Kania Dyah Nastiti, Mohamad wafirul Hadi, Maha Hidayatullah Akbar, Fajri Octadiansyah ===
 +
[[File:Halaman 1 artikell.png]]
 +
[[File:Halaman 2 artikell.png]]
 +
 
 +
[[File:Halaman 3 artikell.png]]
 +
 
 +
[[File:Halaman 4 artikell.png]]
 +
 
 +
[[File:Halaman 5 artikell.png]]
 +
 
 +
[[File:Halaman 6 artikell.png]]
 +
 
 +
[[File:Halaman 7 artikell.png]]
 +
 
 +
[[File:Halaman 8 artikell.png]]
 +
 
 +
[[File:Halaman 9 artikell.png]]
 +
 
 +
[[File:Halaman 10 artikell.png]]
 +
 
 +
[[File:Halaman 11 artikell.png]]
 +
 
 +
[[File:Halaman 12 artikell.png]]
 +
 
 +
[[File:Halaman 13 artikell.png]]
 +
 
 +
 
 +
=== Artikel 4.3.10. Osilasi Pegas dengan Peredam ===
 +
 
 +
Tugas Kolaborasi Artikel Komputasi Teknik
 +
 
 +
Anggota Kelompok :
 +
 
 +
[http://air.eng.ui.ac.id/index.php?title=Aghnia_Ilmiah_Nurhudan 1. Aghnia Ilmiah Nurhudan]
 +
 
 +
[http://air.eng.ui.ac.id/index.php?title=Daniel_Meino_Soedira 2. Daniel Meino Soedira]
 +
 
 +
[http://air.eng.ui.ac.id/index.php?title=Joko.triwardono 3. Joko Triwardono]
 +
 
 +
[http://air.eng.ui.ac.id/index.php?title=Paskal_Rachman 4. Paskal Rachman]
 +
 
 +
Berikut merupakan hasil dari Tugas Kolaborasi Artikel Komputasi Teknik:
 +
 
 +
<gallery mode="slideshow">
 +
File:OSILASI_PEGAS_DENGAN_PEREDAM_Page_1.jpg
 +
File:OSILASI PEGAS DENGAN PEREDAM-1.jpg
 +
File:OSILASI PEGAS DENGAN PEREDAM-2.jpg
 +
File:Osilasi_daniel3.jpg
 +
File:Osilasi_daniel4.jpg
 +
File:OSILASI_PEGAS_DENGAN_PEREDAM_Page_6.jpg
 +
File:OSILASI_PEGAS_DENGAN_PEREDAM_Page_7.jpg
 +
File:Osi7_joko.jpg
 +
File:Osi8_joko.jpg
 +
</gallery>
 +
 
 +
<!-- [[File:OSILASI_PEGAS_DENGAN_PEREDAM_Page_1.jpg|700px]]
 +
 
 +
[[File:OSILASI PEGAS DENGAN PEREDAM-1.jpg|700px]]
 +
 
 +
[[File:OSILASI PEGAS DENGAN PEREDAM-2.jpg|700px]]
 +
 
 +
[[File:Osilasi_daniel3.jpg|700px]]
 +
 
 +
[[File:Osilasi_daniel4.jpg|700px]]
 +
 
 +
[[File:OSILASI_PEGAS_DENGAN_PEREDAM_Page_6.jpg|700px]]
 +
 
 +
[[File:OSILASI_PEGAS_DENGAN_PEREDAM_Page_7.jpg|500px]]
 +
 
 +
[[File:Osi7_joko.jpg]]
 +
 
 +
[[File:Osi8_joko.jpg]] -->
 +
 
 +
=== Artikel 4.3.12. Metode finite diference; Undamped, Linear Case ===
 +
 
 +
[[File:Smoosilasi1.jpg|500 px]]
 +
 
 +
[[File:Smoosilasi2.jpg|500 px]]
 +
 
 +
[[File:Smoosilasi3.jpg|500 px]]
 +
 
 +
[[File:Smoosilasi4.jpg|500 px]]
 +
 
 +
[[File:Smoosilasi5.jpg|500 px]]
 +
 
 +
 
 +
=== Artikel: Penggunaan ANN untuk menyelesaikan sistem osilasi 1D ===
 +
Oleh: Adhika Satyadharma 1906323956
 +
 
 +
 
 +
'''1. Pendahuluan'''
 +
 
 +
Dalam kasus ini, persamaan osilasi 1D akan dicoba diselesaikan dengan menggunakan Artificial Neural Network (ANN). Untuk lebih jelasnya, persamaan yang akan diselesaikan adalah d2x/dt2 + w2x = 0. Adapun untuk menyelesaikan persamaan ini, setidaknya diperlukan 3 buah input yaitu 2 boundary condition (x(t=0) dan dx/dt (t=0)), dan nilai w. Untuk outputnya sendiri, solusi yang ditargetkan dalam kasus ini adalah solusi persamaan defrential tersebut dalam rentang dt.  
 +
 
 +
 
 +
'''2. Detail Teknis'''
 +
 
 +
'''2.1 Initial Thinking'''
 +
 
 +
Karena output yang diinginkan adalah solusi diskrit dalam rentang dt, hal ini memunculkan pertanyaan berapakah nilai dt yang tepat digunakan dan berapa banyak interval yang akan digunakan? Mempertimbangkan agar data dari setiap contoh kasus yang digunakan untuk training dapat dinyatakan dengan format dan ukuran yang sama serta bahwa tidak ada nilai dt yang general, maka nilai dt dibiarkan bebas, jumlah interval diset dalam range 40 hinnga 50 data dan total waktu diasumsikan sebesar 5s.
 +
 
 +
 
 +
Untuk lebih jelasnya berikut variasi dari input yang digunakan:
 +
 
 +
0<= x(t=0) <=1
 +
 
 +
-x(t=0)*w <= dx/dt(t=0) <= x(t=0)*w
 +
 
 +
0 <= w <= 10
 +
 
 +
T = 5
 +
 
 +
40 <= n <= 50
 +
 
 +
dt = T/n
 +
 
 +
 
 +
'''2.2 Training Data Generation'''
 +
 
 +
Data training yang akan digunakan merupakan solusi analitis dari persamaan defrential tersebut yaitu x = x(t=0)*cos(wt+c), c=arcsin(x(t=0)/w/(dx/dt(t=0))). Data solusi analitis ini dihasilkan dengan excel dimana input-input yang diperlukan diset secara random dalam range masing-masing. Secara total 600 contoh kasus dihasilkan dari rumus ini, 500 contoh kasus akan digunakan untuk training ANN dan 100 untuk testing.
 +
 
 +
 
 +
'''2.3. Kode ANN'''
 +
 
 +
Kode ANN yang digunakan dalam kasus ini mengikuti panduan dari https://iamtrask.github.io/2015/07/12/basic-python-network/. Codingan ini dilakukan dalam bahasa phython, yang dieksekusi oleh aplikasi Spyder. Mengenai detail ANN yang digunakan, terdapat 10 hidden neuron digunakan dalam ANN ini. Hanya saja karena dalam kodingan ini tidak terdapat nilai bias dalam setiap hidden neuron, untuk mengkompensasinya digunakanlah input ke-5 yang selalu bernilai 1.
 +
 
 +
Coding ANN yang dilakukan terdapat 2 jenis yaitu training dan testing. Codingan untuk yang training mempunyai format seperti berikut:
 +
 
 +
# Define Sigmoid Function
 +
# Read Data (Input, and Output)
 +
# Normalize
 +
# Define Weights
 +
# Set Iteration
 +
# Forward Propagation
 +
# Backwards Propagation
 +
# Error Calculation
 +
# Backup Data
 +
# Set Break Condition
 +
# Export Weights
 +
 
 +
 
 +
Adapun untuk coding mengenai masalah testing, format yang digunakan adalah sebagai berikut:
 +
 
 +
# Read Data (Input, Analytical Output, Weights)
 +
# Normalize
 +
# Calculate ANN
 +
# Unormalize
 +
# Export Results
 +
# Calculate Error + Print Error
 +
 
 +
(Kode tidak dilampirkan karena cukup panjang dan terdiri dari berberapa file)
 +
 
 +
 
 +
'''3. Hasil'''
 +
 
 +
Sebelum masuk ke hasil, perlu diketahui bahwa training ANN ini belum mencapi hasil yang memuaskan. Error yang dihasilkan dari ANN ini masih cukup besar. Dari hasil proses iterasinya, menurut kami apabila dilakukan lebih banyak iterasi hasilnya dapat lebih baik, namun setelah kami perkirakan (secara kasar) hal ini dapat memakan waktu yang sangat-sangat lama. Karena kami sendiri kurang pengalaman dalam ANN, kami kurang mengetahui cara yang paling efektif untuk mempercepat proses iterasinya sehingga hasilnya akan konvergen.
 +
 
 +
Dari hasil testing sendiri, secara rata-rata terdapat perbedaan nilai 0.123 antara hasil analitis dan ANN. Adapun setelah berberapa contoh kasus data pada testing diplot, dapat terlihat bahwa ada kasus-kasus dimana hasil antara ANN dan analitis cukup mirip, dan ada juga yang berbeda jauh.
 +
 
 +
 
 +
[[File:ANN-Osilasi1D_1312351_1_v1.png]]
 +
 
 +
=== Artikel 4.3.11. Sistem pegas-massa dengan gesekan luncur ===
 +
 
 +
[[File:file1.jpg]]
 +
[[File:file2.jpg]]
 +
[[File:file3.jpg]]
 +
[[File:file4.jpg]]
 +
 
  
 
== Artikel .... Hasil diskusi  : judul ...==
 
== Artikel .... Hasil diskusi  : judul ...==

Latest revision as of 21:18, 2 May 2020

<- back to Studi kasus komputasi teknik

Contents

Knowledge base

Studi kasus

1d oscillating dynamic system 1.png

1d oscillating dynamic system 2.png

1d oscillating dynamic system 3.png

1d oscillating dynamic system 4.png

1d oscillating dynamic system 5.png

1d oscillating dynamic system 6.png

1d oscillating dynamic system 7.png

1d oscillating dynamic system 8.png

1d oscillating dynamic system 9.png

1d oscillating dynamic system 10.png

1d oscillating dynamic system 11.png

1d oscillating dynamic system 12.png

1d oscillating dynamic system 13.png

1d oscillating dynamic system 14.png

1d oscillating dynamic system 15.png

1d oscillating dynamic system 16.png

1d oscillating dynamic system 17.png

1d oscillating dynamic system 18.png

1d oscillating dynamic system 19.png

1d oscillating dynamic system 20.png

1d oscillating dynamic system 21.png

1d oscillating dynamic system 22.png

1d oscillating dynamic system 23.png

1d oscillating dynamic system 24.png

1d oscillating dynamic system 25.png

1d oscillating dynamic system 26.png

1d oscillating dynamic system 27.png

1d oscillating dynamic system 28.png

1d oscillating dynamic system 29.png

1d oscillating dynamic system 30.png

1d oscillating dynamic system 31.png

1d oscillating dynamic system 32.png

1d oscillating dynamic system 33.png

1d oscillating dynamic system 34.png

1d oscillating dynamic system 35.png

Ref. Linge S, Langtangen HP, Programming for Computations - A Gentle Introduction to Numerical Simulations with Python

Terjemahan

4.3.1 Penurunan Model yang Sederhana

alt text
Banyak sistem keteknikan (engineering) berkaitan dengan osilasi, dan persamaan diferensial merupakan kunci utama untuk memahami, memprediksi, dan mengontrol osilasi. Kita mulai dengan model paling sederhana yang berkaitan dengan dinamika penting dari sistem osilasi. suatu benda dengan massa m melekat/dikaitkan pada pegas dan bergerak sepanjang garis tanpa gesekan, lihat Gambar 4.15 di samping untuk sketsa (rolling wheels menunjukkan “tidak ada gesekan”). Ketika pegas diregangkan (atau dikompresi), gaya pegas menarik (atau mendorong) bodi (penampang m) kembali dan bekerja "melawan" gerakan. Lebih tepatnya, misalkan x (t) adalah posisi bodi pada sumbu x, dimana bodi bergerak. Pegas tidak direntangkan ketika x= 0, sehingga gaya adalah nol, dan x= 0 karenanya posisi keseimbangan bodi. Gaya pegas adalah -kx, dimana k adalah konstanta yang diukur. Kami berasumsi bahwa tidak ada gaya lain (mis., Tidak ada gesekan). Hukum Newton ke-2 F=ma kemudian memiliki F=-kx dan a=x ̈ ,

Az 4.41.png

yang dapat ditulis ulang sebagai:

Az 4.42.png

dengan memperkenalkan ω=√(k/m)   (yang sangat umum).

Persamaan (4.42) adalah persamaan diferensial orde kedua, dan oleh karena itu kita memerlukan dua kondisi awal, satu pada posisi x(0) dan satu pada kecepatan x’(0). Di sini kita memilih bodi untuk berhenti, tetapi menjauh dari posisi setimbang:

Az 4.42a.png

Solusi tepat untuk Pers. (4.42) dengan kondisi awal ini adalah x(t)=X0 cosωT. Ini dapat dengan mudah diverifikasi dengan mensubsitusikan ke Pers. (4.42) dan memeriksa kondisi awal. Solusinya mengatakan bahwa sistem massa pegas berosilasi bolak-balik seperti yang dijelaskan oleh kurva kosinus.

Persamaan diferensial (4.42) muncul dalam banyak konteks lainnya. Contoh klasik adalah pendulum sederhana yang berosilasi bolak-balik. Buku-buku fisika berasal, dari hukum gerak kedua Newton, itu diperoleh:

Az 4.42b.png

dimana m adalah massa bodi di ujung pendulum dengan panjang L, g adalah percepatan gravitasi, dan ϴ merupakan sudut yang dibuat pendulum dengan vertikal. Mempertimbangkan sudut kecil ϴ, sin ϴ ≈ ϴ, dan kita dapatkan Pers. (4.42) dengan x = ϴ, ω=√(g/L) , x(0)=Θ, dan x’(0)=0, jika Θ merupakan sudut awal dan pendulum diam di t=0.


4.3.2 Solusi Numerik

Kita telah melihat metode numerik untuk mengendalikan turunan orde kedua, dan beberapa pilihan lainnya merupakan tambahan, akan tetapi kita mengetahu cara menyelesaikan persamaan turunan orde pertama dan bahkan sistem-sistem pada persamaan orde pertama. Dengan hanya sedikit, tetapi cukup umum, cara yang dapat kita tuliskan pada persamaan 4.42 sebagai sebuah sistem orde pertama dari 2 persamaan turunan. Kita memperkenalkan u=x dan v=x^'=u' sebagai 2 fungsi baru yang tidak diketahui. Dua persamaan yang sesuai muncul dari definisi v=u' dan persamaan asal (4.42):

Eviii4.43.JPG

(memperlihatkan bahwa kita dapat menggunakan u"=v') untuk menghilangkan turunan orde kedua dari hokum kedua newton). Selanjutnya kita dapat menerapkan metode forward euler untuk persamaan 4.43 dan 4.44, seperti yang sudah dilakukan pada section 4.2.2:

Eviii4.45.JPG

Sehingga menghasilkan skema komputasi sebagai berikut,

Eviii4.47.JPG


4.3.3 Memprogram Metode Numerik; Kasus Khusus

Program sederhana untuk (4.47) - ( 4.48) mengikuti ide yang sama seperti di bagian 4.2.3:

4.3.3.fadhli.JPG

(Lihat file osc_FE.py.)

Karena kita sudah tahu solusi yang tepat sebagai u(t) = Xo cos ωt , kami beralasan sebagai berikut untuk menemukan interval simulasi yang sesuai [0,T] dan juga berapa poin kita harus memilih. Solusinya memiliki periode P = 2π/ω. (Periode P adalah waktunya perbedaan antara dua puncak u(t) ~ cos ωt curve). Simulasi untuk tiga periode fungsi cosinus, T = 3P, dan memilih Δt sehingga ada 20 Interval per periode menghasilkan Δt = P/20 dan total Nt = T/ Δt = t interval. Sisanya dari program ini adalah pengodean langsung dari skema Forward Euler.

Gambar 4.16 menunjukkan perbandingan antara solusi numerik dan tepat solusi persamaan diferensial. Yang mengejutkan kami, solusi numeriknya terlihat salah. Apakah perbedaan ini disebabkan oleh kesalahan pemrograman atau masalah dengan metode Forward Euler?

Pertama-tama, bahkan sebelum mencoba menjalankan program, Anda harus menghitung dua langkah dalam putaran waktu dengan kalkulator sehingga Anda memiliki beberapa hasil antara untuk dibandingkan. Menggunakan X0 = 2. Dt = 0: 157079632679, dan ω = 2, kita mendapatkan u1 = 2, v = -1,25663706, u2 = 1,80260791, dan v2 = 2,51327412. Perhitungan semacam itu menunjukkan bahwa program itu tampaknya benar. (Kemudian, kita dapat menggunakan nilai-nilai tersebut untuk membangun tes unit dan fungsi tes yang sesuai.)

Simulation of an Oscillating System.PNG

Langkah selanjutnya adalah mengurangi delta t parameter diskritisasi dan melihat apakah hasilnya menjadi lebih akurat. Gambar 4.17 menunjukkan solusi numerik dan tepat untuk kasus delta t = P / 40; P / 160; P / 2000. Hasilnya jelas menjadi lebih baik, dan resolusi terakhir memberikan grafik yang tidak dapat dibedakan secara visual. Namun demikian, resolusi terakhir melibatkan 6000 interval komputasi secara total, yang dianggap cukup banyak. Namun, ini bukan masalah pada laptop modern, karena perhitungan hanya membutuhkan sepersekian detik.

Meskipun 2000 interval per periode osilasi tampaknya cukup untuk solusi numerik yang akurat, grafik kanan bawah pada Gambar 4.17 menunjukkan bahwa jika kita meningkatkan waktu simulasi, di sini hingga 20 periode, ada sedikit pertumbuhan amplitudo, yang menjadi signifikan dari waktu ke waktu. . Kesimpulannya adalah bahwa metode Forward Euler memiliki masalah mendasar dengan amplitudo yang tumbuh, dan bahwa diperlukan delta yang sangat kecil untuk mencapai hasil yang memuaskan. Semakin lama simulasi, semakin kecil Delta t. Sudah pasti saatnya untuk mencari metode numerik yang lebih efektif!

Simulation with different steps.PNG

4.3.4 Sebuah Penyelesaian dari Metode Numerik

Dalam skema Forward Euler,

4.3.4.(1).JPG

kita dapat mengganti u^n pada persamaan terakhir dengan nilai u^n+1 yang baru dihitung dari persamaan pertama:

4.3.4.(2).JPG

Sebelum membenarkan perbaikan ini secara matematis, mari kita coba pada contoh sebelumnya. Hasilnya muncul pada Gambar 4.18. Kita melihat bahwa amplitudo tidak tumbuh, tetapi fase tidak sepenuhnya benar. Setelah 40 periode (Gbr. 4.18 kanan) kita melihat signifikan perbedaan antara solusi numerik dan tepat. Penurunan t menurun kesalahan. Misalnya, dengan 2000 interval per periode, kami hanya melihat fase kecil kesalahan bahkan setelah 50.000 periode (!). Kita dapat menyimpulkan bahwa perbaikan tersebut menghasilkan metode numerik yang sangat baik! Mari kita tafsirkan skema yang disesuaikan secara matematis. Pertama kami memesan (4,49) - (4,50) sedemikian rupa sehingga perbedaan pendekatan terhadap derivatif menjadi transparan:

4.3.4.(10).JPG

4.3.4.(3).JPG

Kami menafsirkan (4,51) sebagai persamaan diferensial sampel pada titik mesh tn, karena kami memiliki vn di sisi kanan. Sisi kiri kemudian perbedaan maju atau Meneruskan perkiraan Euler ke turunan u0 , lihat Gambar 4.2. Di samping itu, kami menginterpretasikan (4,52) sebagai persamaan diferensial sampel pada titik mesh tnC1, karena kami miliki di sisi kanan.

4.3.4.(4).jpeg

Dalam hal ini, perbedaan aproksimasi pada sisi kiri adalah perbedaan ke belakang,

4.3.4.(5).jpeg


Gambar 4.19 mengilustrasikan perbedaan mundur. Kesalahan dalam perbedaan mundur sebanding dengan t, sama seperti untuk perbedaan maju (tetapi konstanta proporsionalitas dalam istilah kesalahan memiliki tanda yang berbeda). Diskretisasi yang dihasilkan metode untuk (4,52) sering disebut sebagai skema Backward Euler.

Untuk meringkas, gunakan perbedaan maju untuk persamaan pertama dan mundur Perbedaan untuk hasil persamaan kedua dalam metode yang jauh lebih baik daripada hanya menggunakan maju perbedaan dalam kedua persamaan.

Cara standar untuk mengekspresikan skema ini dalam fisika adalah dengan mengubah urutan persamaan,

4.3.4.(6).jpeg

dan terapkan perbedaan maju ke (4,53) dan perbedaan mundur ke (4,54):

4.3.4.(7).jpg

Artinya, pertama kecepatan v diperbarui dan kemudian posisi u, menggunakan kecepatan yang paling baru dihitung. Tidak ada perbedaan antara (4,55) - (4,56) dan (4,49) - (4,50) sehubungan dengan akurasi, jadi urutan persamaan diferensial asli tidak apa-apa. Skema (4.55) - (4.56) berada di bawah nama Semi-implisit Euler4 atau Euler-Cromer. Implementasi (4.55) - (4.56) ditemukan dalam file osc_EC.py. Inti dari kode itu seperti

4.3.4.(8).jpg


4.3.4.(9).jpg



4.3.5 Metode Runge-Kutta orde 2 (atau metode Heun)

Sebuah metode yang cukup populer digunakan untuk menyelesaikan persamaan diferensial biasa (ODE) vector dan skalar adalah metode Runge-Kutta Order (RK2), atau biasa dikenal dengan metode Heun. Ide dasar pada metode ini, yang pertama untuk ODE skalar, adalah dengan membentuk aproksimasi perbedaan terpusat (centered difference) terhadap turunan antara dua titik waktu yang didefinisikan sebagai berikut:

Kolab1.JPG

Formula dari centered difference tersebut dapat digambarkan melalui Gambar 4.20. Error pada aproksimasi centered difference ini proporsional terhadap nilai ∆t2, 1 order lebih tinggi dibandingkan dengan pendekatan forward and backward difference, yang berarti nilai jika kita memiliki sebuah nilai ∆t, maka error nya akan berkurang secara effektif dengan menggunakan centered difference karena nilai error tersebut berkurang dengan faktor 4, daripada faktor 2.

Kolab2.JPG

Permasalahan yang ada pada skema centered difference semacam ini untuk persamaan ODE secara umum, u’=f(u,t) adalah kita mendapatkan

[File:Kolab3.JPG]]

Yang mana ini akan menyulitkan karena kita tidak mengetahui berapa nilai un+1/2. Namun demikian, ktia dapat mengaproksimasi nilai f diantara dua level waktu dengan menggunakan rata-rata aritmatik dari nilai f tesebut pada saat tn dan tn+1  :

Kolab4.JPG

Kemudian hasilnya adalah : Results435.jpg Dimana berupa persamaan aljabar nonlinear untuk File:F435.jpeg dan bukan fungsi linear dari u. sehingga untuk menyelesaikan fungsi F4351.jpg tanpa menyelesaikannya dengan persamaan nonlinear, dapat diprediksi F4352.jpg menggunakan persamaan Forward Euler: F4353.jpg Sehingga dapat digunakan metode F4354.jpg metode tersebut dapat diaplikasikan untuk ODEs scalar dan vector.

Untuk system osilasi dengan F4355.jpg

Pada file osc_Heun.py terdapat implementasinya. File tersebut menjalankan simulasi untuk 10 period dengan 20 kali langkah per periode.

Solusi Numerical dan eksak yang berkaitan dengan ini terdapat di fig. 4.21. dapat diliat bahwa amplitude meningkat namun tidak sebanyak pada metode forward euler. Bgaimanapun juga, metode forward euler adalah yang terbaik. Perlu diingat juga bahwa metode forward euler memberikan prediksi yang lebih baik, seperti contohnya untuk persoalan pertumbuhan/peluruhan, atau SIR mode. Akan tetapi metode orde 2 runge-kutta atau metod heun juga bisa dipertimbangkan. Meskipun untuk menyelesaikan persoalan osilasi, metode euler sudah terbaik.


4.3.6 Perangkat Lunak untuk Menyelesaikan ODEs

Terdapat banyak metode yang dapat digunakan untuk menyelesaikan ODEs, dan alangkah baiknya kita memilih akses yang mudah untuk mengimplementasikannya ke berbagai metode, terutama metode adaptif yang canggih dan kompleks yang dapat menyesuaikan nilai Δt secara otomatis untuk mendapatkan nilai akurasi yang ditentukan. Phyton Odespy3 merupakan salah satu perangkat yang dapat memberikan akses yang mudah ke berbagai metode numerik untuk menyelesaikan ODEs.

Salah satu contoh termudah dalam penggunaan Odespy adalah untuk menyelesaikan masalah u’ = u, u(0) = 2, untuk 100 time steps sampai t = 4:

import odespy

def f(u, t):

   return u

method = odespy.Heun #or, e.g., odespy.ForwardEuler solver = method(f) solver.set_initial _condition(2) time_points = np.linspace(0, 4, 101) u. t = solver.solve (time_points)

Dengan kata lain, kalian mendefinisikan sebuah fungsi f(u, t), menginisialisasi sebuah objek penyelesaian Odespy, mengatur kondisi awal, menghitung titik waktu pengumpulan dimana anda menginginkan solusinya, dan bertanya mengenai solusinya. Variabel arrays u dan t dapat dibuat menjadi sebuah grafik secara langsung, yaitu: plot(t,u).

Fitur menarik yang dimiliki oleh Odespy ialah parameter permasalahan dapat menjadi sebuah argumen pada fungsi f(u, t) penggunanya. Sebagai contoh, apabila permasalahan ODE kita adalah u’ = -au + b, dengan 2 parameter yaitu a dan b, kita dapat menuliskan fungsi f kita menjadi

def f(u, t, a, b):

   return -a*u + b

Sebagai tambahan, permasalahan yang bergantung pada argumen a dan b dapat ditransfer ke fungsi ini bila kita mengumpulkan nilainya dalam sebuah daftar atau tuple ketika membuat sebuah pemecahan Odespy dan menggunakan argumen f_args:

a = 2 b = 1 solver = method(f, f_args=[a, b])

Hal ini merupakan sebuah fitur yang baik karena parameter permasalahan haruslah selain sebagai sebuah variabel global – sekarang dapat menjadi sebuah argument dalam fungsi kita secara alami.

Menggunakan Odespy untuk menyelesaikan osilasi ODEs seperti u” + ω2u = 0, diformulasikan sebagai sebuah sistem u’ = v dan v’ = -ω2u, dilakukan sebagai berikut. Kita tentukan sebuah nilai time steps per periode dan hitung time steps yang diasosiasikan serta waktu akhir simulasi (T), cantumkan sebuah nilai periode untuk disimulasikan:

Import odespy

  1. Define the ODE system
  2. u’ = v
  3. v’ = -omega**2*u

def f(sol, t, omega=2):

     u, v = sol
     return [v, -omega**2*u]
  1. Set and compute problem dependent parameters

omega = 2 X_0 = 1 number_of_periods = 40 time_intervals_per_period = 20 from numpy import pi, linspace, cos P = 2*pi/omega #length of one period # length of one period dt = P/time_intervals_per_period # time step T = number_of_periods*P # final simulation time

  1. Create Odespy solver object

odespy_method = odespy.RK2 solver = odespy_method(f, f_args=[omega])

  1. The initial condition for the system is collected in a list

Solver.set_initial_condition([X_0, 0])

  1. Compute the desired time points where we want the solution

N_t = int(round(T/dt)) # no of time intervals Time_points = linspace(0, T, N_t+1)

  1. Solve the ODE problem

sol, t = solver.solve(time_points)

  1. Note: sol contains both displacement and velocity
  2. extract original variables

u = sol[:,0] v = sol[:,1]

Dua pernyataan terakhir menjadi penting karena dua fungsi u dan v di dalam sistem ODE tersebut tergabung bersama dalam sebuah array di dalam pemecahan Odespy. Solusi pada sistem ODE ditunjukan sebagai array 2 dimesi dimana kolom pertama (sol[:,0]) disimpan sebagai u dan kolom kedua (sol[:,1]) disimpan sebagai v. Mengeplot u dan v merupakan sebuah masalah dalam menjalankan plot(t, u, t, v).

Catatan

Di dalam fungsi tersebut kita menuliskan f(sol, t, omega) dibandingkan menulis f (u, t, omega) untuk mengindikasikan bahwa solusi pada f adalah solusi pada waktu t dimana nilai u dan t tergabung bersama: sol = [u,v]. Kita dapat juga menggunakan u sebagai argumen:

def f(u, t, omega=2):

     u, v = u
     return [v, -omega**2*u]

Ini hanya berarti kita mendefinisikan ulang nama u pada fungsi tersebut untuk merata-ratakan solusi pada waktu t untuk komponen pertama pada sistem ODE tersebut.

Untuk beralih ke metode numerik lain, tinggal substitusikan RK2 dengan nama yang sesuai dari metode yang diinginkan. Mengetik pydoc odespy pada terminal window memunculkan daftar dari metode yang dijalankan. Cara yang sangat sederhana dalam memilih metode ini menyarankan penambahan yang jelas dari kode diatas: kita dapat menentukan daftar metode, menjalankan semua metode, dan membandingkan setiap kurva u pada sebuah plot. Sebagaimana odespy juga mengandung skema Euler-Cromer, kita menulis kembali sistem ini dengan v’ = -w2u sebagai ODE pertama dan u’ = v sebagai ODE kedua, karena ini adalah pilihan standar ketika menggunakan metode Euler-Cromer (juga pada odespy):

def f(u, t, omega=2): v, u = u return [-omega**2*u, v]

Perubahan persamaan ini juga mempengaruhi kondisi awal: komponen pertama adalah nol dan yang kedua adalah X_0 maka kita perlu melewati daftar [0, X_0] untuk solver.set_ initial_condition.

Kode osc_odespy.py mengandung detail:

def compare(odespy_methods, 
omega, 
X_0, 
number_of_periods, 
time_intervals_per_period=20): 
from numpy import pi, linspace, cos 
P = 2*pi/omega                  # length of one period 
dt = P/time_intervals_per_period 
T = number_of_periods*P
# If odespy_methods is not a list, but just the name of 
# a single Odespy solver, we wrap that name in a list 
# so we always have odespy_methods as a list 
if type(odespy_methods) != type([]): 
odespy_methods = [odespy_methods] 
# Make a list of solver objects 
solvers = [method(f, f_args=[omega]) for method in 
odespy_methods] 
for solver in solvers: 
solver.set_initial_condition([0, X_0]) 
# Compute the time points where we want the solution 
dt = float(dt) # avoid integer division 
N_t = int(round(T/dt)) 
time_points = linspace(0, N_t*dt, N_t+1) 
legends = [] 
for solver in solvers: 
sol, t = solver.solve(time_points) 
v = sol[:,0] 
u = sol[:,1] 
# Plot only the last p periods 
p = 6 
m = p*time_intervals_per_period # no time steps to plot 
plot(t[-m:], u[-m:]) 
hold(’on’) 
legends.append(solver.name()) 
xlabel(’t’) 
# Plot exact solution too 
plot(t[-m:], X_0*cos(omega*t)[-m:], ’k--’) 
legends.append(’exact’) 
legend(legends, loc=’lower left’) 
axis([t[-m], t[-1], -2*X_0, 2*X_0]) 
title(’Simulation of %d periods with %d intervals per period’ 
% (number_of_periods, time_intervals_per_period)) 
savefig(’tmp.pdf’); savefig(’tmp.png’) 
show()

Fitur baru pada kode ini adalah kemampuan untuk mem-plot hanya periode p terakhir, yang memperbolehkan kita untuk menjalankan long time simulations dan melihat hasil akhir tanpa plot yang berantakan dengan terlalu banyak periode. Syntax t[-m:] mem-plot elemen m terakhir dalam t (indeks negatif dalam hitungan susunan/daftar Pyhton dari akhir).

Kita bisa membandingkan metode Heun (atau setara metode RK2) dengan skema Euler-Crome:

compare(odespy_methods=[odespy.Heun, odespy.EulerCromer], 
omega=2, X_0=2, number_of_periods=20, 
time_intervals_per_period=20)

Gambar 4.22 menunjukkan bagaimana metode Heun (garis biru dengan piringan kecil) memiliki error yang cukup besar pada amplitude dan fase sesudah setelah periode 14-20 (kiri atas), namun menggunakan sebanyak tiga kali langkah waktu membuat kurvanya hampir sama (kanan atas). Akan tetapi setelah periode 194-200 error tersebut telah berkembang (kiri bawah), tetapi dapat cukup dikurangi dengan mengurangi separuh langkah waktu (kanan bawah).

Dengan semua metode di Odespy, sekarang menjadi mudah untuk mulai menjelajahi metode-metode lain, seperti perbedaan mundur (backward differences) bukannya perbedaan maju (forward differences) yang digunakan dalam skema Forward Euler. Latihan 4.17 mengatasi permasalahan tersebut.

Odespy berisi metode adaptif yang cukup canggih di mana pengguna "dijamin" untuk mendapatkan solusi dengan akurasi yang ditentukan. Tidak ada jaminan matematis, tetapi error untuk sebagian besar kasus tidak akan menyimpang secara signifikan dari toleransi pengguna yang mencerminkan keakuratan. Metode yang sangat populer dari jenis ini adalah metode Runge-Kutta-Fehlberg, yang menjalankan metode Runge-Kutta orde 4 dan menggunakan metode Runge-Kutta orde 5 untuk memperkirakan error sehingga dapat disesuaikan untuk menjaga error di bawah toleransi. Metode ini juga dikenal luas sebagai ode45, karena itulah nama fungsi yang mengimplementasikan metode ini di Matlab. Kita dapat dengan mudah menguji metode Runge-Kutta-Fehlberg segera setelah kita tahu nama Odespy yang sesuai, yaitu RKFehlberg:

compare(odespy_methods=[odespy.EulerCromer, odespy.RKFehlberg], 
omega=2, X_0=2, number_of_periods=200, 
time_intervals_per_period=40)

Oscillating17-2.png

Perhatikan bahwa argumen time_intervals_per_period mengacu pada titik waktu di mana kami ingin solusinya. Poin-poin ini juga yang digunakan untuk perhitungan numerik dalam pemecah odespy.EulerCromer, sedangkan pemecah odespy.RKFehlberg akan menggunakan satu set titik waktu yang tidak diketahui karena interval waktu disesuaikan ketika metode berjalan. Orang dapat dengan mudah melihat titik-titik yang sebenarnya digunakan oleh metode karena ini tersedia sebagai himpunan solver.t_all (tetapi merencanakan atau memeriksa titik-titik membutuhkan modifikasi di dalam metode perbandingan).

Gambar 4.23 menunjukkan contoh komputasi di mana metode Runge-Kutta-Fehlberg jelas lebih unggul daripada skema Euler-Cromer dalam simulasi yang lama, tetapi perbandingannya tidak terlalu adil karena metode Runge-Kutta_Fehlberg berlaku sekitar dua kali lebih banyak langkah waktu dalam hal perhitungan ini dan melakukan lebih banyak pekerjaan per langkah waktu. Ini adalah tugas yang cukup rumit untuk membandingkan dua metode yang sangat berbeda dalam cara yang wajar sehingga pekerjaan komputasi versus akurasi dilaporkan secara ilmiah dengan baik.

Oscillating18-2.png

4.3.7 Metode Runge-Kutta Orde 4

Metode Runge-Kutta Orde 4 adalah metode yang sering digunakan secara luas untuk menyelesaikan ODEs, karena menghasilkan data dengan tingkat akurasi yang tinggi bahkan dalam time step yang tidak terlalu kecil.

1-.PNG

Algoritma; Pertama-tama kita nyatakan algoritma 4-stage

2-.PNG

Dimana

3-.PNG

4-.PNG

5-.PNG


Aplikasi; Kita bisa menjalankan simulasi seperti pada Figs. 4.16, 4.18, dan 4.21, untuk 40 periode. 10 periode terakhir ditunjukan melalui Fig. 2.24. Hasil yang ditunjukan terlihat impresif sebagaimana penggunaan metode Euler-Cromer.


Implementasi; Tingkatan dalam metode Runge-Kutta orde-4 bisa dengan mudah diimplementasikan sebagai modifikasi dari osc_Heun.py code. Sebagai alternatif, salah satu dapat menggunakan osc_odespy.py code dengan menyediakan argumen odespy_methods-[odespy.RK4] untuk membandingkan fungsi.


Derivasi; Derivasi dari metode Runge-Kutta orde-4 dapat disajikan dengan cara pedagogis yang menyatukan banyak elemen fundamental dari teknik diskritisasi numerik dan bisa menggambarkan banyak aspek “numerical thinking ”ketika membangun perkiraan metode solusi.

Kita mulai dengan mengintegrasikan general ODE 6-.PNG dari waktu ke waktu, mulai dari tn sampai t(n_1),

9-.PNG

Tujuan dari komputasi 10-.PNG, ketika 11-.PNG pada saat ini lebih dikenal dengan nilai u. Tantangan mengintegralkan muncul ketika integrand mengandung u yang tidak diketahuai antara tn sampai t(n+1).

Integral tersebut dapat diperkirakan dengan menggunakan Simpson’s rule yang telah terkenal

12-.PNG

Permasalahan dengan persamaan ini adalah kita tidak mengetahui nilai dari 13-.PNG dan 14-.PNG karena hanya u^n yang tersedia dan hanya f^n yang dapat dihitung.

Untuk melanjutkan, idenya dalah menggunakan berbagai perkiraan untuk 15-.PNG dan 16-.PNG berdasarkan penggunaan skema yang telah diketahui untuk ODE dalam interval 17-.PNG dan 18-.PNG. Mari kita bagi persamaan integral menjadi empat suku.


19-.PNG

Dimana C01.JPG, C02.JPG, dan C03.JPG adalah pendekatan untuk C04.JPG dan C05.JPG yang dapat digunakan pada perhitungan. Untuk C01.JPG dapat menggunakan pendekatan untuk C06.JPG berdasarkan tahap Forward Euler pada size C14(2).JPG


4-63.JPG


Persamaan ini mempermudah prediksi C04.JPG, sehingga untuk C02.JPG kita dapat mencoba metode Backward Euler untuk memperkirakan C06.JPG


4-64.JPG


Dengan C02.JPG sebagai pendekatan untuk C04.JPG, pada akhirnya bentuk akhir dari C03.JPG dapat menggunakan metode midpoint (atau central difference, juga disebut metode Crank-Nicholson) untuk memperkirakan C15(2).JPG.


4-65.JPG


Kita telah menggunakan metode Forward dan Backward Euler, juga centered difference approximation pada konteks Simpsons rule. Diharapkan kombinasi dari metode ini dapat menghasilkan overall time stepping dari C07.JPG ke C08.JPG yang lebih akurat dibandingkan individual steps (yang memiliki error proportional dengan C09.JPG dan C10(2).JPG). Hal ini benar bahwa: error numerik yang terjadi seperti C11(2).JPG Untuk konstanta C, artinya error lebih cepat mendekati nol ketika time step size dikurangi, dibandingkan dengan metode Forward Euler C12(2).JPG, metode Euler-Cromer C12(2).JPG,atau Runge Kutta orde 2, atau metode Heuns C13(2).JPG.

Perhatikan bahwa Metode Runge-Kutta Orde 4 sepenuhnya eksplisit jadi tidak diperlukan untuk menyelesaikannya dengan persamaan aljabar baik secara linier maupun non linier, terlepas dari apa yang terlihat pada f. Namun nilai kesetabilannya kondisional dan bergantung pada nilai f tersebut. Ada sebuah bagian besar dari metode implisit Runge-Kutta yang nilai kesetabilannya tidak kondisional. namun diperlukan solusi dari persamaan aljabar yang melibatkan nilai f pada setiap "time step". Odespy dapat dimanfaatkan untuk mendukung penyelesaian dari banyak metode Runge-Katta yang eksplisit. Tetapi belum bisa digunakan untuk metode Runge-Katta yang implisit.

4.3.8 Efek Lain : Damping, Nonlinearity, dan external force

Model permasalahan u’’ + ω2u = 0 adalah model matematika yang paling simple untuk oscilating system. Namun, Model ini lebih banyak membutuhkan metode numerik, seperti yang sudah kita lihat, dan sangat berguna untuk menjadi tolak ukur untuk mengevaluasi kinerja dari metode numerik.

Dalam Pengaplikasian dikehidupan nyata lebih banyak melibatkan efek fisika, yang mengarahkan ke persamaan diferensial dengan ketentuan yang lebih banyak dan juga lebih kompleks. biasanya, memiliki kekuatan redaman f (u ') dan pegas s (u). Kedua gaya ini tergantung pada nonlinear dari uraiannya, u’ atau u. sebagai tambahan, gaya lingkungan F(t) jufga bekerja pada sistem. Contohnya, pendulum klasik memiliki “pegas” nonlinear atau mengembalikan gaya s(u) ~ sin (u), dan gaya tahan dari udara pada pendulum menyebabkan terjadinya gaya redam f(u’) ~ |u’|u. Contoh dari gaya lingkungan adalah getaran dari tanah (seperti gempa) dan juga seperti ombak atau angina.

Dengan tiga jenis gaya yang bekerja pada sistem : F(t), f(u’), dan s(u). maka dapat ditulis persamaan F(t) – f(u’) – s(u). Tanda mines didepan f dan s menunjukan bahwa fungsi ini didefinisikan sebagai gaya yang melawan gerakan. Sebagai Contoh, Pegas yang terpasang pada roda mobil dikombinasikan dengan beberapa perdeam yang efektif. Masing-masing memiliki gaya redam f(u’) = bu’ yang bekerja melawan kecepatan pegas u’. gaya fisika yang sesuai dapat dtulis –f: -bu’, yang menunjuk ke bawah saat pegas diregangkan (dan poin u’ ke atas), sedangkan -f bertindak ke atas saat pegas dikompresi (dan poin u’ ke bawah).

Gambar 4.25 menunjukan contoh dari massa m terpasang dengan pegas nonlinear dan dashpot, dan bersubyek pada gaya lingkungan F(t). Namun, model umum yang kita miliki dapat juga digunakan pada pendulum pada gambar 4.26 dengan s (u) = m g sin θ dan f (u ̇) = 1/2 C_D Aϱ(θ|) ̇θ| (Dimana CD = 0.4, A adalah area perpotongan dari body dan ϱ adalah densitas udara)

Gambar425.png

Gambar 4.25 Sistem Oscillating General

Hukum Newton kedua untuk sistem yang dapat ditulis dengan akselerasi waktu massa pada sisi kiri dan gaya pada sisi kanan:

438rumus1.png

Bagaimanapun persamaan ini lebih umum disusun ulang menjadi

438rumus2.png

Karena persamaan diferensial adalah orde 2, disebabkan oleh istilah u^, kita membutuhkan dua kondisi awal:

438rumus3.png

Gambar426.png

Gambar 4.26 Sebuah pendulum dengan gaya

Catat bahwa dengan pilihan 438rumus4.png kita memperoleh kembali persamaan diferensial biasa 438rumus5.png

Bagaimana kita bisa menyelesaikan (4.66)? sebagaimana persamaan diferensial biasa yang simpel 438rumus6.png kita mulai dengan menulis ulang persamaan diferensial biasa orde 2 sebagai sebuah sistem dari dua persamaan diferensial biasa orde 1:

438rumus8.png

Kondisi awal menjadi

438rumus9.png

Setiap metode dari sebuah sistem persamaan diferensial biasa orde 1 dapat digunakan untuk menyelesaikan 438rumus10.png

The Euler-Cromer scheme

Sebuah pilihan atraktif dari sebuah implementasi, akurasi dan efisiensi sudut pandang adalah skema Euler-Cromer dimana kita mengambil sebuah perbedaan kedepan pada (4.68) dan perbedaan kebelakang pada (4.69):

438rumus11.png

Kita dapat dengan mudah menyelesaikan 438rumus12.png yang tidak diketahui:

438rumus13.png


kata kata dalam perintah ODEs

Perintah ODE dalam sistem ODE penting untuk model yang diperluas (4.68) - (4.69). Bayangkan kita menulis persamaan untuk u’ terlebih dahulu dan kemudian untuk v’. Metode Euler-Cromer akan menggunakan forward difference untuk u^n+1 dan kemudian backward difference akan menggunakannya untuk v^n+1. Yang Terkhir akan menyebabkan persamaan nonlinear algebraic untuk v^n+1

(4.3.8) 1.png

jika f(v) adalah fungsi nonlinear dari v Ini akan membutuhkan metode numerik untuk persamaan aljabar nonlinier untuk mencari v^n+1 saat memperbarui v^n+1 melalui forward difference memberikan persamaan untuk v^n+1 itu linear dan sepele untuk dipecahkan dengan tangan.

File osc_EC_general.pymemiliki fungsi Euler Cromer yang mengimplementasikan metode ini:

(4.3.8) 2.png

(4.3.8) 3.png

Metode Runge-Kutta orde ke 4

Metode RK4 hanya mengevaluasi sisi kanan sistem ODE,

(4.3.8) 4.png

untuk nilai-nilai u, v, dan t yang diketahui, maka metode ini sangat sederhana untuk digunakan terlepas dari bagaimana fungsi s(u) dan f(v)dipilih.



4.3.9 ilustrasi redaman linier

Kami menganggap sistem rekayasa dengan pegas linier, s(u) = kx, dan peredam kental, di mana gaya peredaman adalah porpotional terhadap u', f(u') = bu', untuk beberapa konstanta b > 0. Pilihan ini dapat memodelkan sistem pegas vertikal di dalam mobil (tetapi insinyur sering suka menggambarkan sistem tersebut dengan massa bergerak horizontal seperti yang digambarkan pada Gambar 4.25). kita dapat memilih nilai-nilai sederhana untuk konstanta untuk mengilustrasikan efek dasar redaman (dan kegembiraan selanjutnya). Memilih osilasi sebagai fungsi u(t) = cos t sederhana dalam kasus tak teredam, kita dapat menetapkan m = 1, k = 1, b = 0,3, Uo = 1, Vo = 0. Fungsi berikut mengimplementasikan kasus ini:

Wafi 439-1.png

Fungsi plot_u adalah kumpulan plot untuk merencanakan u(t), atau bagian darinya. Gambar 4.27 menunjukkan efek dari bu': kita memiliki osilasi dengan (perkiraan) periode 2π, seperti yang diharapkan, tetapi amplitudo teredam secara efisien.

Kania 439-2.png


Komentar mengenai pekerjaan dengan masalah berskala

Alih-alih menetapkan b = 0,3 dan m = k = Uo = 1 sebagai nilai fisik yang “tidak mungkin”, akan lebih baik untuk skala persamaan mu" + bu' + ku = 0. Ini mengartikan bahwa kita memasukan variabel independen dan dependen yang tak berdimensi :

Kania 439-3.png

Di mana tc dan uc adalah ukuran karakteristik waktu dan perpindahan, sehingga Kania 439-5.png dan Kania 439-6.png memiliki ukuran tipikal mereka didekat kesatuan. Dalam masalah ini, kita dapat memilih Kania 439-7.png dan Kania 439-8.png. Ini memberikan masalah yang berskala (atau tanpa dimensi) berikut untuk kuantitas tak berdimensi Kania 439-9.png:

Kania 439-4.png

Faktnya adalah hanya ada satu parameter fisik di kasus ini: angka β. Menyelesaikan masalah ini begitu juga terkait dengan masalah utama dengan parameter yaitu m = k = Uo = 1 dan b = β. Tetapi untuk menyelsaikan masalah dengan satuan lebih umum: jika kita memdapatkan solusi ¯u(¯t;β), kita dapat menemukan solusi fisik pada kasus ini, dikarenakan :

439rumus.png

Selama β didapat, kita dapat menemukan u untuk Uo , k, dan m dengan rumus diatas, dengan begitu pengerjaan simulasi dapat dipersingkat waktu. Ini menunjukkan pengerjaan dengan skala atau masalah satuan.


4.3.10. Ilustrasi Redaman Linier Dengan Eksitasi Sinusoidal

Sekarang kita akan memperluas contoh sebelumnya untuk menambah beberapa gaya osilasi eksternal pada sistem: F (t) = Asin (wt). Mengendarai mobil di jalan dengan lonjakan sinusoidal mungkin memberikan eksitasi eksternal pada sistem pegas di mobil (w terkait dengan kecepatan mobil).

    from math import pi,sin
    w = 3
    A = 0.5
    F = lambda t: A*sin(w*t)

kita dapatkan grafik pada gambar 4.28 .Perbedaan yang mencolok dari Gambar 4.27 adalah bahwa osilasi dimulai sebagai sinyal cos t teredam tanpa banyak pengaruh gaya eksternal, tetapi kemudian osilasi bebas dari sistem yang tidak teredam (cos t) u’’ + u = 0 mati dan gaya eksternal 0: 5 sin.(3t) menimbulkan osilasi dengan periode yang lebih pendek 2phi/3. Dianjurkan untuk menggunakan beberapa nilai A yang lebih besar dan beralih dari sinus ke acosinus dalam F dan mengamati efeknya. Jika mencarinya di dalam buku fisika, Anda dapat menemukan solusi analitik yang tepat untuk masalah persamaan diferensial dalam kasus ini.

4.3.11. Sistem pegas-massa dengan gesekan luncur

Sebuah benda dengan massa m bekerja pada sebuah pegas dengan kekakuan k saat meluncur pada sebuah bidang permukaan. Benda tersebut mengalami gaya gesek f(u') disebabkan terjadi kontak antara benda dengan bidang permukaan seperti terlihat pada Gambar 4.30. Gaya gesek f(u') dapat dimodelkan dengan gesekan Coulomb sebagai berikut:

Eq4.3.11.1.png

Dimana μ adalah koefisien gesek, dan mg merupakan gaya normal pada bidang permukaan benda yang bergerak. Formula ini dapat juga ditulis sebagai f(u') = μmg sign (u'), dengan syarat fungsi signum sign (x) didefinisikan nol untuk x = 0 (numpy.sign mempunyai sifat ini). Untuk memastikan bahwa signum dari definisi f benar, ingat bahwa gaya fisis aktual adalah -f dan positif (misal f<0) ketika gaya tersebut bekerja berlawanan dengan benda yang bergerak dengan kecepatan u'<0.

Gambar 4.30 Sketsa dari sebuah subjek sistem osilasi dinamis untuk gesekan luncur dan gaya pegas satu dimensi.

Gaya pegas nonlinear diambil sebagai:

Eq4.3.11.2.png

Yang mana nilai -ku diperkirakan untuk nilai u yang kecil, namun stabil pada ±k/α untuk nilai ±αu yang besar. Berikut adalah plot dengan k=1000 dan u ∈ [-0.1,0.1] untuk tiga nilai α:

591px-1d oscillating dynamic system 30.1.png

Jika tidak ada gaya eksitasi eksternal yang bekerja pada benda, maka persamaan gerak yang kita dapatkan adalah:

Eq4.3.11.3.png

Mari kita simulasikan situasi dimana sebuah benda dengan massa 1 kg meluncur pada bidang permukaan dengan μ = 0.4, terikat pada pegas dengan kekakuan k = 1000 kg/s^2. Perpindahan awal benda adalah 10 cm, dan parameter α dalam s(u) diatur pada 60 1/m.

Dengan menggunakan fungsi EulerCromer dari kode osc_EC_general, kita dapat menulis fungsi sliding_friction untuk menyelesaikan masalah ini:

Def    sliding_friction():
       from numpy import tanh, sign
       
       f = lambda v: mu*m*g*sign(v)
       alpha = 60.0
       s = lambda u: k/alpha*tanh(alpha*u)
       F = lambda t: 0
       
       g = 9.81
       mu = 0.4
       m = 1
       k = 1000
       
       U_0 = 0.1
       V_0 = 0
       
       T = 2
       dt = T/5000
       
       u, v, t = EulerCromer(f=f, s=s, F=F, m=m, T=T, 
                             U_0=U_0, V_0=V_0, dt=dt)
       plot_u(u, t)

Setelah menjalankan fungsi sliding_friction memberi kita hasil seperti pada Gambar. 4.31 dengan s(u)= -k/α tanh⁡(αu), (kiri) dan versi linierisasi s(u)=ku (kanan).

Gambar 4.31 Efek pegas nonlinear (kiri) dan linier (kanan) pada gesekan luncur

4.3.12. Metode finite diference; Undamped, Linear Case

Selanjutnya kita akan membahas metode numerik untuk ODE orde kedua

 u^+ω^2 u=0,   u(0)=U_0,u^' (0)=0,t∈(0,T]

tanpa menulis ulang ODE sebagai sistem ODE orde pertama. Motivasi utama untuk "metode solusi lain" adalah bahwa prinsip-prinsip diskritisasi menghasilkan skema yang sangat baik, dan yang lebih penting, pemikiran seputar diskritisasi bisa digunakan kembali ketika memecahkan persamaan diferensial parsial. Gagasan utama dari metode numerik ini adalah untuk memperkirakan urutan kedua turunan u dengan selisih terbatas. Sementara ada beberapa pilihan perbedaan perkiraan untuk derivatif orde pertama, ada satu rumus yang mendominasi untuk turunan orde kedua:

Persamaan4.74.jpg

Error dalam perkiraan tersebut proporsional terhadap ∆t^2. Biarkan ODE valid di beberapa titik waktu yang berubah ubah t_n,

 u^ (t_n )+ ω^2 u (t_n )=0

Selanjutnya memasukkan rumus perkiraan (4.74) diatas, sehingga di dapatkan

Persamaan4.75.jpg

Sekarang diasumsikan bahwa u^(n-1) dan u^n sudah dihitung, dan u^(n+1) adalah yang baru tidak diketahui. Memecahkan sehubungan dengan u^(n+1)

Persamaan4.76.jpg

Masalah besar muncul ketika kita ingin memulai skema. Kita tahu bahwa u^0 = U_0, tetapi menerapkan (4,76) untuk n=0 untuk menghitung u^1

Persamaan4.77.jpg

Dimana kita tidak mengetahui U-1. Kondisi awal U’ (0) = 0 dapat membantu kiti untuk menghilangkan U-1 dan kondisi ini bagaimanapun juga harus dimasukkan dalam beberapa cara. Untuk tujuan ini, kami mendiskritasikan u’(0) = 0 dengan perbedaan terpusat,


Persamaan4.78.jpg

Oleh karena itu, u-1 = u1, dan kita dapat menggunakan relasi ini untuk menghilangkan u1 di persamaan 4.77

Persamaan4.79.png

Dengan U0 = U0 dan u1 dihitung dari persamaan (4.78), kita dapat menghitung u2, u3, dan seterusnya dari persamaan (4.76). Latihan 4.19 meminta Anda untuk mengeksplorasi bagaimana langkah-langkah di atas diubah seandainya kita memiliki kondisi awal bukan nol u’ (0) = V0

Kita dapat memperkirakan kondisi awal U’(0) dengan menggunakan Forward difference

Persamaan4.80.jpg

Mengarah pada u1 = u0 . lalu kita dapat menggunakan persamaan (4.76) untuk langkah selanjutnya . Walaupun forward difference memiliki kesalahan proporsional ke ∆t . dimana centered difference yang kita gunakan memiliki error proporsional ke ∆t2. Yang dimana kompatibel dengan akurasi (erro yang enunjukan ∆t2) yang digunakan dalam diskritisasi persamaan diferensial. Metode untuk ODE orde kedua yang dijelaskan di atas berjalan di bawah nama metode Störmer atau integrasi Verlet 7. Ternyata metode ini secara matematis setara dengan skema Euler-Cromer Atau lebih tepatnya, rumus umum (4,76) setara dengan rumus Euler-Cromer.


4.3.13 Metode finite diference; damping linier

Sebuah isu kunci adalah bagaimana untuk mengkonferensi skema dari daerah 4.3.12 ke persamaan diferensial dengan lebih banyak istilah. Kita mulai dengan kasus linear penempatan f (u') = bu', kemungkinan gaya per nonlinear s(u), dan sebuah gaya excitation F(t):

4.79.png


Kita harus cari perkiraan perbedaan yang tepat untuk u' di dalam bu'. Sebuah pilihan yang baik adalah perbedaan berpusat


4.80.png

Sampling persamaan pada titik tn,

4.80a.png

Dan memasukkan perkiraan perbedaan terhingga pada u" dan u' hasil dalam


4.81.png


Dimana F" adalah notasi pendek untuk F(t). Persamaan (4.81) adalah linear dalam u^(n+1) tak diketahui kita dapat dengan mudah memecahkan untuk kuantitas ini:


4.82.png


Pada kasus tanpa redaman, kita membutuhkan formula khusus untuk u1. kondisi awal U`(0) = 0 menyatakan bahwa u-1 = u1, dan dengan persamaan (4.82) untuk n = 0, kita mendapatkan.

4.8.3casees.JPG

Pada kasus yang lebih unun dengan sebuah bentuk redaman nonlinier f(u`),


4.8.3.2.1.2.JPG


Kita mendapatkan

4.8.3.2.1.2.1.JPG

Dimana sebauh persamaan ajabar non linier untuk un+1 bahwa harus diseleseikan dengan metode numerik. Skema lebih bagus diperoleh dari penggunaan "backward difference" untuk u`,

4.8.3.2.1.2.1.2.JPG

Karena pada bagian redaman akan lebih diketahui, yang hanya melibatkan un dan un-1, dan kita dapat dengan mudah menyelesaikan untuk un+1. Kelemahan dari backward difference dibandingkan dengan centered difference (4.80) adalah ini mengurangi urutan akurasi dalam skema keseluruhan dari ∆t2 ke ∆t. Pada kenyataanya, skema Euler-Cromer mengevaluasi istilah redaman nonlinear sebagai f(vn), saat menghitung vn+1, dan ini setara dengan menggunakan backward difference di atas. Akibatnya, kenyamanan skema Euler-Cromer untuk redaman nonlinier datang dengan konsekuensi menurunkan akurasi keseluruhan skema dari urutan kedua ke urutan pertama pada ∆t. Menggunakan trik yang sama dalam skema beda hingga {finite difference} untuk persamaan diferensial orde kedua, yaitu, menggunakan backward difference dalam f(u’), membuat skema ini sama bagus dan akuratnya seperti skema Euler-Cromer pada kasus nonlinier umum mu”+f(u’)+s(u) = F.


Hasil tugas kaloborasi Oscillating 1-D Dynamic Artikel 1

Pembagian Tema

4.3.1 Oldy Fahlovvi, Muchalis Zikramansyah Masuku, Ahmad Zikri, Muhammad Irfan Dzaky

4.3.2 Andhika, faturahman

4.3.4 iqbal & Alghi & Adam, aji suryadi

4.3.5 Shabrina & Edo, jerry, raihan

4.3.6 ronald & Desy & yophie, ardi

4.3.7 Kania & Chandra, evi & Dieter

4.3.9 Wafirul & fajri & keni, maha

4.3.10 Daniel Meino Soedira & Paskal Rachman, Joko Triwardono & Aghnia Ilmiah Nurhudan

4.3.11 Bambang ali.

4.3.12 Bagus, maheka, adzanna, Adinda

4.3.13 Harry, wisnu, Ichwan, fadli

Artikel Kolaborasi : 1-D OSCILLATING SYSTEM arranged by Oldy Fahlovvi, Muchalis Zikramansyah Masuku, Ahmad Zikri, Muhammad Irfan Dzaky

Berikut ini kami lampirkan tugas kolaborasi tentang 1-D OSCILLATING SYSTEM dalam bentuk slideshow.

Artikel 1 Hasil diskusi : OSCILLATING 1-D DYNAMIC SYSTEM with 1 MASS, 3 SPRING and 1 DAMPING

Tema 4.3.13 Linear Damping Oscillation

Artikel Kolaborasi - Sistem Osilasi Satu Dimensi

Berikut adalah tugas Artikel Kolaborasi kelompok kami mengenai Penggunaan Perangkat Lunak Python untuk Menyelesaikan ODEs pada Sistem Mekanik Berosilasi dengan anggota sebagai berikut:

1. Ardy Lefran Lololau

2. I Gusti Agung Ayu Desy Wulandari

3. Ronald Akbar

4. Yophie Dikaimana

Selanjutnya berikut ini adalah hasil diskusi kelompok kami mengenai sistem osilasi satu dimensi dengan menggunakan metode Artificial Neural Network (ANN)

Tugas kolaborasi 4.3.4. ARTIKEL OSCILLATING ORDE 1 DUMPING PADA SISTEM SPRING MOBIL

Tugas kolaborasi 4.3.5 sistem osilasi satu dimensi runge - kutta : studi kasus shock breaker motor

Collab13 (5).jpg Collab13 (6).jpg Collab13 (7).jpg Collab13 (8).jpg Collab13 (9).jpg Collab13 (10).jpg Collab13 (11).jpg Collab13 (12).jpg Collab13 (13).jpg Collab13 (14).jpg Collab13 (15).jpg Collab13 (16).jpg Collab13 (1).jpg Collab13 (2).jpg Collab13 (3).jpg Collab13 (4).jpg


Artikel Tugas 4.3.9 Illustrasi redaman linear Kania Dyah Nastiti, Mohamad wafirul Hadi, Maha Hidayatullah Akbar, Fajri Octadiansyah

Halaman 1 artikell.png Halaman 2 artikell.png

Halaman 3 artikell.png

Halaman 4 artikell.png

Halaman 5 artikell.png

Halaman 6 artikell.png

Halaman 7 artikell.png

Halaman 8 artikell.png

Halaman 9 artikell.png

Halaman 10 artikell.png

Halaman 11 artikell.png

Halaman 12 artikell.png

Halaman 13 artikell.png


Artikel 4.3.10. Osilasi Pegas dengan Peredam

Tugas Kolaborasi Artikel Komputasi Teknik

Anggota Kelompok :

1. Aghnia Ilmiah Nurhudan

2. Daniel Meino Soedira

3. Joko Triwardono

4. Paskal Rachman

Berikut merupakan hasil dari Tugas Kolaborasi Artikel Komputasi Teknik:


Artikel 4.3.12. Metode finite diference; Undamped, Linear Case

Smoosilasi1.jpg

Smoosilasi2.jpg

Smoosilasi3.jpg

Smoosilasi4.jpg

Smoosilasi5.jpg


Artikel: Penggunaan ANN untuk menyelesaikan sistem osilasi 1D

Oleh: Adhika Satyadharma 1906323956


1. Pendahuluan

Dalam kasus ini, persamaan osilasi 1D akan dicoba diselesaikan dengan menggunakan Artificial Neural Network (ANN). Untuk lebih jelasnya, persamaan yang akan diselesaikan adalah d2x/dt2 + w2x = 0. Adapun untuk menyelesaikan persamaan ini, setidaknya diperlukan 3 buah input yaitu 2 boundary condition (x(t=0) dan dx/dt (t=0)), dan nilai w. Untuk outputnya sendiri, solusi yang ditargetkan dalam kasus ini adalah solusi persamaan defrential tersebut dalam rentang dt.


2. Detail Teknis

2.1 Initial Thinking

Karena output yang diinginkan adalah solusi diskrit dalam rentang dt, hal ini memunculkan pertanyaan berapakah nilai dt yang tepat digunakan dan berapa banyak interval yang akan digunakan? Mempertimbangkan agar data dari setiap contoh kasus yang digunakan untuk training dapat dinyatakan dengan format dan ukuran yang sama serta bahwa tidak ada nilai dt yang general, maka nilai dt dibiarkan bebas, jumlah interval diset dalam range 40 hinnga 50 data dan total waktu diasumsikan sebesar 5s.


Untuk lebih jelasnya berikut variasi dari input yang digunakan:

0<= x(t=0) <=1

-x(t=0)*w <= dx/dt(t=0) <= x(t=0)*w

0 <= w <= 10

T = 5

40 <= n <= 50

dt = T/n


2.2 Training Data Generation

Data training yang akan digunakan merupakan solusi analitis dari persamaan defrential tersebut yaitu x = x(t=0)*cos(wt+c), c=arcsin(x(t=0)/w/(dx/dt(t=0))). Data solusi analitis ini dihasilkan dengan excel dimana input-input yang diperlukan diset secara random dalam range masing-masing. Secara total 600 contoh kasus dihasilkan dari rumus ini, 500 contoh kasus akan digunakan untuk training ANN dan 100 untuk testing.


2.3. Kode ANN

Kode ANN yang digunakan dalam kasus ini mengikuti panduan dari https://iamtrask.github.io/2015/07/12/basic-python-network/. Codingan ini dilakukan dalam bahasa phython, yang dieksekusi oleh aplikasi Spyder. Mengenai detail ANN yang digunakan, terdapat 10 hidden neuron digunakan dalam ANN ini. Hanya saja karena dalam kodingan ini tidak terdapat nilai bias dalam setiap hidden neuron, untuk mengkompensasinya digunakanlah input ke-5 yang selalu bernilai 1.

Coding ANN yang dilakukan terdapat 2 jenis yaitu training dan testing. Codingan untuk yang training mempunyai format seperti berikut:

  1. Define Sigmoid Function
  2. Read Data (Input, and Output)
  3. Normalize
  4. Define Weights
  5. Set Iteration
  6. Forward Propagation
  7. Backwards Propagation
  8. Error Calculation
  9. Backup Data
  10. Set Break Condition
  11. Export Weights


Adapun untuk coding mengenai masalah testing, format yang digunakan adalah sebagai berikut:

  1. Read Data (Input, Analytical Output, Weights)
  2. Normalize
  3. Calculate ANN
  4. Unormalize
  5. Export Results
  6. Calculate Error + Print Error

(Kode tidak dilampirkan karena cukup panjang dan terdiri dari berberapa file)


3. Hasil

Sebelum masuk ke hasil, perlu diketahui bahwa training ANN ini belum mencapi hasil yang memuaskan. Error yang dihasilkan dari ANN ini masih cukup besar. Dari hasil proses iterasinya, menurut kami apabila dilakukan lebih banyak iterasi hasilnya dapat lebih baik, namun setelah kami perkirakan (secara kasar) hal ini dapat memakan waktu yang sangat-sangat lama. Karena kami sendiri kurang pengalaman dalam ANN, kami kurang mengetahui cara yang paling efektif untuk mempercepat proses iterasinya sehingga hasilnya akan konvergen.

Dari hasil testing sendiri, secara rata-rata terdapat perbedaan nilai 0.123 antara hasil analitis dan ANN. Adapun setelah berberapa contoh kasus data pada testing diplot, dapat terlihat bahwa ada kasus-kasus dimana hasil antara ANN dan analitis cukup mirip, dan ada juga yang berbeda jauh.


ANN-Osilasi1D 1312351 1 v1.png

Artikel 4.3.11. Sistem pegas-massa dengan gesekan luncur

File1.jpg File2.jpg File3.jpg File4.jpg


Artikel .... Hasil diskusi  : judul ...