Walker: Integrating the skew-normal SDE
This example runs Walker to integrate the skew-normal SDE (see DiffEq/SkewNormal.h) using constant coefficients.
Control file
title "Example problem" walker term 10.0 # Max time dt 0.001 # Time step size npar 10000 # Number of particles ttyi 1000 # TTY output interval rngs mkl_r250 seed 1 end end skew-normal depvar m init zero coeff const ncomp 2 T 1.0 3.5 end sigmasq 0.04 0.25 end lambda 100.0 -50.0 end rng mkl_r250 end statistics interval 2 <m1m1> <m2m2> end pdfs interval 10 filetype txt policy overwrite centering elem format scientific precision 4 p1( M1 : 1.0e-2 ) p2( M2 : 1.0e-2 ) end end
Example run on 4 CPUs
./charmrun +p4 Main/walker -v -c ../../tmp/test.q -u 0.9
Output
Running on 4 processors: Main/walker -v -c ../../tmp/skew.q -u 0.9 charmrun> /usr/bin/setarch x86_64 -R mpirun -np 4 Main/walker -v -c ../../tmp/skew.q -u 0.9 Charm++> Running on MPI version: 3.0 Charm++> level of thread support used: MPI_THREAD_SINGLE (desired: MPI_THREAD_SINGLE) Charm++> Running in non-SMP mode: numPes 4 Converse/Charm++ Commit ID: b8b2735 CharmLB> Load balancer assumes all CPUs are same. Charm++> Running on 1 unique compute nodes (4-way SMP). Charm++> cpu topology info is gathered in 0.000 seconds. ,::,` `. .;;;'';;;: ;;# ;;;@+ +;;; ;;;;;, ;;;;. ;;;;;, ;;;; ;;;; `;;;;;;: ;;; :;;@` :;;' .;;;@, ,;@, ,;;;@: .;;;' .;+;. ;;;@#:';;; ;;;;' ;;;# ;;;: ;;;' ;: ;;;' ;;;;; ;# ;;;@ ;;; ;+;;' .;;+ ;;;# ;;;' ;: ;;;' ;#;;;` ;# ;;@ `;;+ .;#;;;. ;;;# :;;' ;;;' ;: ;;;' ;# ;;; ;# ;;;@ ;;; ;# ;;;+ ;;;# .;;; ;;;' ;: ;;;' ;# ,;;; ;# ;;;# ;;;: ;@ ;;; ;;;# .;;' ;;;' ;: ;;;' ;# ;;;; ;# ;;;' ;;;+ ;', ;;;@ ;;;+ ,;;+ ;;;' ;: ;;;' ;# ;;;' ;# ;;;' ;;;' ;':::;;;; `;;; ;;;@ ;;;' ;: ;;;' ;# ;;;';# ;;;@ ;;;:,;+++++;;;' ;;;; ;;;@ ;;;# .;. ;;;' ;# ;;;;# `;;+ ;;# ;# ;;;' .;;; :;;@ ,;;+ ;+ ;;;' ;# ;;;# ;;; ;;;@ ;@ ;;;. ';;; ;;;@, ;;;;``.;;@ ;;;' ;+ .;;# ;;; :;;@ ;;; ;;;+ :;;;;;;;+@` ';;;;;'@ ;;;;;, ;;;; ;;+ +;;;;;;#@ ;;;;. .;;;;;; .;;#@' `#@@@: ;::::; ;:::: ;@ '@@@+ ;:::; ;:::::: :;;;;;;. __ __ .__ __ .;@+@';;;;;;' / \ / \_____ | | | | __ ___________ ` '#''@` \ \/\/ /\__ \ | | | |/ // __ \_ __ \ \ / / __ \| |_| <\ ___/| | \/ \__/\ / (____ /____/__|_ \\___ >__| \/ \/ \/ \/ < ENVIRONMENT > ------ o ------ * Build environment: -------------------- Hostname : sprout Executable : walker Version : 0.1 Release : LA-CC-XX-XXX Revision : e26d8f8514a11ade687ba460f42dfae5af53d4d6 CMake build type : DEBUG Asserts : on (turn off: CMAKE_BUILD_TYPE=RELEASE) Exception trace : on (turn off: CMAKE_BUILD_TYPE=RELEASE) MPI C++ wrapper : /opt/openmpi/1.8/clang/system/bin/mpicxx Underlying C++ compiler : /usr/bin/clang++-3.5 Build date : Fri Feb 6 06:39:01 MST 2015 * Run-time environment: ----------------------- Date, time : Sat Feb 7 07:41:35 2015 Work directory : /home/jbakosi/code/quinoa/build/clang Executable (rel. to work dir) : Main/walker Command line arguments : '-v -c ../../tmp/skew.q -u 0.9' Output : verbose (quiet: omit -v) Control file : ../../tmp/skew.q Parsed control file : success < FACTORY > ---- o ---- * Particle properties data layout policy (CMake: LAYOUT): --------------------------------------------------------- particle-major * Registered differential equations: ------------------------------------ Unique equation types : 8 With all policy combinations : 18 Legend: equation name : supported policies Policy codes: * i: initialization policy: R-raw, Z-zero * c: coefficients policy: C-const, J-jrrj Beta : i:RZ, c:CJ Diagonal Ornstein-Uhlenbeck : i:RZ, c:C Dirichlet : i:RZ, c:C Gamma : i:RZ, c:C Generalized Dirichlet : i:RZ, c:C Ornstein-Uhlenbeck : i:RZ, c:C Skew-Normal : i:RZ, c:C Wright-Fisher : i:RZ, c:C < PROBLEM > ---- o ---- * Title: Example problem ------------------------ * Differential equations integrated (1): ---------------------------------------- < Skew-Normal > kind : stochastic dependent variable : m initialization policy : Z coefficients policy : C start offset in particle array : 0 number of components : 2 random number generator : MKL R250 coeff T [2] : { 1 3.5 } coeff sigmasq [2] : { 0.04 0.25 } coeff lambda [2] : { 100 -50 } * Output filenames: ------------------- Statistics : stat.txt PDF : pdf * Discretization parameters: ---------------------------- Number of time steps : 18446744073709551615 Terminate time : 10 Initial time step size : 0.001 * Output intervals: ------------------- TTY : 1000 Statistics : 2 PDF : 10 * Statistical moments and distributions: ---------------------------------------- Estimated statistical moments : <M1> <M2> <m1m1> <m2m2> PDFs : p1(M1:0.01) p2(M2:0.01) PDF output file type : txt PDF output file policy : overwrite PDF output file centering : elem Text floating-point format : scientific Text precision in digits : 4 * Load distribution: -------------------- Virtualization [0.0...1.0] : 0.9 Load (number of particles) : 10000 Number of processing elements : 4 Number of work units : 40 (39*250+250) * Time integration: Differential equations testbed -------------------------------------------------- Legend: it - iteration count t - time dt - time step size ETE - estimated time elapsed (h:m:s) ETA - estimated time for accomplishment (h:m:s) out - output-saved flags (S: statistics, P: PDFs) it t dt ETE ETA out --------------------------------------------------------------- 1000 1.000000e+00 1.000000e-03 000:00:02 000:00:26 SP 2000 2.000000e+00 1.000000e-03 000:00:05 000:00:23 SP 3000 3.000000e+00 1.000000e-03 000:00:08 000:00:20 SP 4000 4.000000e+00 1.000000e-03 000:00:11 000:00:17 SP 5000 5.000000e+00 1.000000e-03 000:00:14 000:00:14 SP 6000 6.000000e+00 1.000000e-03 000:00:18 000:00:12 SP 7000 7.000000e+00 1.000000e-03 000:00:21 000:00:09 SP 8000 8.000000e+00 1.000000e-03 000:00:24 000:00:06 SP 9000 9.000000e+00 1.000000e-03 000:00:27 000:00:03 SP 10000 1.000000e+01 1.000000e-03 000:00:30 000:00:00 SP Normal finish, maximum time reached: 10.000000 * Timers (h:m:s): ----------------- Initial conditions : 0:0:0 Migration of global-scope data : 0:0:0 Total runtime : 0:0:30 [Partition 0][Node 0] End of program
Results
The left figure shows the first two moments indicating convergence to a statistically stationary state. The right one shows the estimated PDFs with their analytical solution (see DiffEq/SkewNormal.h).
Gnuplot commands to reproduce the above plots:
plot "stat.txt" u 2:3 w l t "<M1>", "stat.txt" u 2:4 w l t "<M2>", "stat.txt" u 2:5 w l t "<m1m1>", "stat.txt" u 2:6 w l t "<m2m2>" plot "pdf_p1.txt" w p, "pdf_p2.txt" w p, exp(-x*x/2/0.2/0.2)*(1+erf(100.0*x/sqrt(2)))/0.2/sqrt(2*pi) lt 1, exp(-x*x/2/0.5/0.5)*(1+erf(-50.0*x/sqrt(2)))/0.5/sqrt(2*pi) lt 2