Berikut adalah source code simulasi Monte Carlo pada kasus peluruhan :
Source code dalam bentuk .cpp dapat didownload disini.
//********************************************************************************
#include <iostream>
#include <ctime>
#include <cstdlib>
#include <cmath>
#include <fstream>
using namespace std;
int main()
{
double N0 = 10000000; //inti induk awal;
double N0_an =N0; //untuk analitik
double lambda=0.2; //hitung konstanta peluruhan (konversi)
double p = 1-(exp(-lambda)); //peluang meluruh
double hl_step=0;
double rat_nN0=0; //peluang meluruh
double sum_hl=0; //half life wadah kosong
double hl_an=log(2)/lambda; //analitical half-life
int t=1; //untuk t=0 tanpa iterasi
srand( time(0) ); //avoid same generated number
int n = 0; // counter untuk menghitung inti yang meluruh
ofstream file_; //buat file
file_.open ("decay_out10000000.txt"); //buka, lalu tulis file
cout <<"#t(sec)"<<"\t"<<"#Undecayed (MC)"<<"\t"<<"#Undecayed (analitic)"<< "\n"; //saat t masih 0
file_<<t-1<<"\t"<<N0<< "\t"<< N0_an<< "\n"; //saat t masih 0
cout <<t-1<<"\t"<<N0<< "\t"<< N0_an<< "\n"; //saat t masih 0
Source code dalam bentuk .cpp dapat didownload disini.
while (N0>0) // hanaya run ketika : kondisi berhenti (according to both N0 or time)
{
for (int i = 0; i < N0; ++i)
{
double random=static_cast <double> (rand()) /( static_cast <double> (RAND_MAX)); //random keluaran
if (random<p){ //marking decayed isotop ("1" dice)
n = ++n; }
}
if (n>0){ //jika tak ada yang meluruh (n=0) random tidak dicount
rat_nN0=(n/N0); //rasio n/N0 untuk menghitung half life (estimasi MC)
hl_step=(log(0.5)/(log(1-(rat_nN0))));
sum_hl=sum_hl+hl_step;
N0=N0-n; //sisa inti induk MC
double N0_a=N0_an*(exp(-lambda*t)); //sisa inti induk (analitik)
cout <<t<<"\t"<<N0<<"\t"<<N0_a<<"\t"<<"\t"<<hl_step<<"\t"<<hl_an<<"\n";
file_<<t<<"\t"<<N0<<"\t"<<N0_a<<"\t"<<"\t"<<hl_step<<"\t"<<hl_an<<"\n";
if (N0<1){
break;
}
n=0; // set ulang counter
t=++t;
}
}
file_.close ();
double hl_MC=sum_hl/t;
cout<<"\n \n" ;
cout<<"--------------------------------\n";
cout<<"half-life (analitic)"<<"\t"<<"="<<hl_an<<"(s)\n";
cout<<"half-life (MC)"<<"\t \t"<<"="<<hl_MC<<"(s)\n";
file_.open ("hl_out10000000.txt");
file_<<hl_MC<<"\t"<<hl_an<<"\n";
file_.close ();
return 0;
}
//*******************************************************************************
Source code dalam bentuk .cpp dapat didownload disini.
Contoh hasilnya adalah :
Penjelasan lengkap dalam bentuk slide tersedia pada link youtube di sini.
Membutuhkan jasa pembuatan simulasi Matlab, C++ , dll ? Klik di sini
Simulasi Peluruhan Menggunakan Metode Acak Monte Carlo C++
May 12, 2017
3
kak, saya sangat tertarik dgn koding ini. tapi ketika dibuka di mcnp kok tidak mau terlihat inputnya ya?
ReplyDeletekak, saya sangat tertarik dgn koding ini. tapi ketika dibuka di mcnp kok tidak mau terlihat inputnya ya?
ReplyDeleteini C++ biasa mba, bukan code MCNP. Thanks
ReplyDelete