Walker: Integrating the mass-fraction beta SDE
This example runs Walker to integrate the mass-fraction beta SDE (see DiffEq/MassFractionBeta.h) using constant coefficients. The mass-fraction beta SDE is based on the beta SDE, additionally computing two more stochastic variables that are functions of the beta variables integrated. For more detail on the beta SDE, see https:/
Control file
title "Test Ray's closure ideas for <y^2> and <rho v>" walker #nstep 1 # Max number of time steps term 25.0 # Max time dt 0.002 # Time step size npar 10000 # Number of particles ttyi 500 # TTY output interval rngs mkl_r250 end end massfracbeta # Select the mass-fraction beta SDE depvar Y # Dependent variable: Y, denoting mass fractions ncomp 15 # = 3x5 = 5 systems each computing 3 variables: # Y - mass fraction, # R - density, # V - specific volume init zero coeff const # alpha = Sb/kappa, beta = (1-S)b/kappa kappa 2.0 0.76923 0.5 0.15873 0.1 end b 0.4 1.0 1.0 1.0 100.0 end S 0.5 0.53846 0.5 0.39683 0.5 end rng mkl_r250 rho2 1.0 1.0 1.0 1.0 1.0 end #r 0.0101 0.0101 0.0101 0.0101 0.0101 end # low-A r 9.0 9.0 9.0 9.0 9.0 end # high-A end statistics #precision 5 #format default # <Y>, mass fraction means <Y1> # 3 <Y2> # 4 <Y3> # 5 <Y4> # 6 <Y5> # 7 # <rho>, mean densities <Y6> # 8 <Y7> # 9 <Y8> # 10 <Y9> # 11 <Y10> # 12 # <V>, mean specific volumes <Y11> # 13 <Y12> # 14 <Y13> # 15 <Y14> # 16 <Y15> # 17 # <y^2>, mass fraction variances <y1y1> # 23 <y2y2> # 24 <y3y3> # 25 <y4y4> # 26 <y5y5> # 27 # <rho^2>, density variances <y6y6> # 28 <y7y7> # 30 <y8y8> # 32 <y9y9> # 34 <y10y10> # 36 # <v^2>, specific volume variances <y11y11> # 38 <y12y12> # 39 <y13y13> # 40 <y14y14> # 41 <y15y15> # 42 # <rho v>, density-specific-volume covariances <y6y11> # 29 <y7y12> # 31 <y8y13> # 33 <y9y14> # 35 <y10y15> # 37 # <rho v^2> <Y6y11y11> # 18 <Y7y12y12> # 19 <Y8y13y13> # 20 <Y9y14y14> # 21 <Y10y15y15> # 22 end pdfs interval 100 filetype txt policy overwrite centering elem format scientific precision 4 p1( Y1 : 1.0e-2 ) p2( Y2 : 1.0e-2 ) p3( Y3 : 1.0e-2 ) p4( Y4 : 1.0e-2 ) p5( Y5 : 1.0e-2 ) end end
Example run on 8 CPUs
./charmrun +p8 Main/walker -v -c ../../tmp/massfracbeta.q
Output
Running on 8 processors: Main/walker -v -c massfracbeta.q charmrun> /usr/bin/setarch x86_64 -R mpirun -np 8 Main/walker -v -c massfracbeta.q 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 8 Converse/Charm++ Commit ID: 63927de CharmLB> Load balancer assumes all CPUs are same. Charm++> Running on 1 unique compute nodes (8-way SMP). Charm++> cpu topology info is gathered in 0.002 seconds. ,::,` `. .;;;'';;;: ;;# ;;;@+ +;;; ;;;;;, ;;;;. ;;;;;, ;;;; ;;;; `;;;;;;: ;;; :;;@` :;;' .;;;@, ,;@, ,;;;@: .;;;' .;+;. ;;;@#:';;; ;;;;' ;;;# ;;;: ;;;' ;: ;;;' ;;;;; ;# ;;;@ ;;; ;+;;' .;;+ ;;;# ;;;' ;: ;;;' ;#;;;` ;# ;;@ `;;+ .;#;;;. ;;;# :;;' ;;;' ;: ;;;' ;# ;;; ;# ;;;@ ;;; ;# ;;;+ ;;;# .;;; ;;;' ;: ;;;' ;# ,;;; ;# ;;;# ;;;: ;@ ;;; ;;;# .;;' ;;;' ;: ;;;' ;# ;;;; ;# ;;;' ;;;+ ;', ;;;@ ;;;+ ,;;+ ;;;' ;: ;;;' ;# ;;;' ;# ;;;' ;;;' ;':::;;;; `;;; ;;;@ ;;;' ;: ;;;' ;# ;;;';# ;;;@ ;;;:,;+++++;;;' ;;;; ;;;@ ;;;# .;. ;;;' ;# ;;;;# `;;+ ;;# ;# ;;;' .;;; :;;@ ,;;+ ;+ ;;;' ;# ;;;# ;;; ;;;@ ;@ ;;;. ';;; ;;;@, ;;;;``.;;@ ;;;' ;+ .;;# ;;; :;;@ ;;; ;;;+ :;;;;;;;+@` ';;;;;'@ ;;;;;, ;;;; ;;+ +;;;;;;#@ ;;;;. .;;;;;; .;;#@' `#@@@: ;::::; ;:::: ;@ '@@@+ ;:::; ;:::::: :;;;;;;. __ __ .__ __ .;@+@';;;;;;' / \ / \_____ | | | | __ ___________ ` '#''@` \ \/\/ /\__ \ | | | |/ // __ \_ __ \ \ / / __ \| |_| <\ ___/| | \/ \__/\ / (____ /____/__|_ \\___ >__| \/ \/ \/ \/ < ENVIRONMENT > ------ o ------ * Build environment: -------------------- Hostname : karman Executable : walker Version : 0.1 Release : LA-CC-XX-XXX Revision : ab18ecb7bd1b27706d963643c4c1ae144a21bd22 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 : Mon Aug 3 15:17:35 MDT 2015 * Run-time environment: ----------------------- Date, time : Tue Aug 4 08:01:36 2015 Work directory : /home/jbakosi/code/quinoa/build/clang Executable (rel. to work dir) : Main/walker Command line arguments : '-v -c massfracbeta.q' Output : verbose (quiet: omit -v) Control file : massfracbeta.q Parsed control file : success < FACTORY > ---- o ---- * Particle properties data layout policy (CMake: LAYOUT): --------------------------------------------------------- particle-major * Registered differential equations: ------------------------------------ Unique equation types : 12 With all policy combinations : 56 Legend: equation name : supported policies Policy codes: * i: initialization policy: R - raw Z - zero D - delta B - beta * c: coefficients policy: C - const D - decay H - homogeneous decay M - Monte Carlo homogeneous decay Beta : i:BDRZ, c:C Diagonal Ornstein-Uhlenbeck : i:BDRZ, c:C Dirichlet : i:BDRZ, c:C Gamma : i:BDRZ, c:C Generalized Dirichlet : i:BDRZ, c:C Mass-fraction beta : i:BDRZ, c:C Mix mass-fraction beta : i:BDRZ, c:DHM Mix number-fraction beta : i:BDRZ, c:D Number-fraction beta : i:BDRZ, c:C Ornstein-Uhlenbeck : i:BDRZ, c:C Skew-Normal : i:BDRZ, c:C Wright-Fisher : i:BDRZ, c:C < PROBLEM > ---- o ---- * Title: Test Ray's closure ideas for <y^2> and <rho v> ------------------------------------------------------- * Differential equations integrated (1): ---------------------------------------- < Mass-fraction beta > kind : stochastic dependent variable : Y initialization policy : Z coefficients policy : C start offset in particle array : 0 number of components : 5 random number generator : MKL R250 coeff b [5] : { 0.4 1 1 1 100 } coeff S [5] : { 0.5 0.53846 0.5 0.39683 0.5 } coeff kappa [5] : { 2 0.76923 0.5 0.15873 0.1 } coeff rho2 [5] : { 1 1 1 1 1 } coeff r [5] : { 9 9 9 9 9 } * Output filenames: ------------------- Statistics : stat.txt PDF : pdf * Discretization parameters: ---------------------------- Number of time steps : 18446744073709551615 Terminate time : 25 Initial time step size : 0.002 * Output intervals: ------------------- TTY : 500 Statistics : 1 PDF : 100 * Statistical moments and distributions: ---------------------------------------- Estimated statistical moments : <Y1> <Y2> <Y3> <Y4> <Y5> <Y6> <Y6y11y11> <Y7> <Y7y12y12> <Y8> <Y8y13y13> <Y9> <Y9y14y14> <Y10> <Y10y15y15> <Y11> <Y12> <Y13> <Y14> <Y15> <y1y1> <y2y2> <y3y3> <y4y4> <y5y5> <y6y6> <y6y11> <y7y7> <y7y12> <y8y8> <y8y13> <y9y9> <y9y14> <y10y10> <y10y15> <y11y11> <y12y12> <y13y13> <y14y14> <y15y15> Stats floating-point format : default Stats text precision, digits : 6 Estimated PDFs : p1(Y1:0.01) p2(Y2:0.01) p3(Y3:0.01) p4(Y4:0.01) p5(Y5:0.01) PDF output file type : txt PDF output file policy : overwrite PDF output file centering : elem PDF text floating-point format : scientific PDF text precision, digits : 4 * Load distribution: -------------------- Virtualization [0.0...1.0] : 0 Number of processing elements : 8 Number of work units : 8 User load (# of particles) : 10000 Chunksize (load per work unit) : 1250 Actual load (# of particles) : 10000 (=8*1250) * 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 --------------------------------------------------------------- 500 1.000000e+00 2.000000e-03 000:00:02 000:01:02 S 1000 2.000000e+00 2.000000e-03 000:00:04 000:00:54 S 1500 3.000000e+00 2.000000e-03 000:00:06 000:00:50 S 2000 4.000000e+00 2.000000e-03 000:00:09 000:00:47 S 2500 5.000000e+00 2.000000e-03 000:00:11 000:00:46 S 3000 6.000000e+00 2.000000e-03 000:00:13 000:00:44 S 3500 7.000000e+00 2.000000e-03 000:00:16 000:00:41 S 4000 8.000000e+00 2.000000e-03 000:00:18 000:00:39 S 4500 9.000000e+00 2.000000e-03 000:00:20 000:00:36 S 5000 1.000000e+01 2.000000e-03 000:00:23 000:00:34 S 5500 1.100000e+01 2.000000e-03 000:00:25 000:00:32 S 6000 1.200000e+01 2.000000e-03 000:00:28 000:00:30 S 6500 1.300000e+01 2.000000e-03 000:00:30 000:00:27 S 7000 1.400000e+01 2.000000e-03 000:00:32 000:00:25 S 7500 1.500000e+01 2.000000e-03 000:00:35 000:00:23 S 8000 1.600000e+01 2.000000e-03 000:00:37 000:00:21 S 8500 1.700000e+01 2.000000e-03 000:00:40 000:00:18 S 9000 1.800000e+01 2.000000e-03 000:00:42 000:00:16 S 9500 1.900000e+01 2.000000e-03 000:00:44 000:00:14 S 10000 2.000000e+01 2.000000e-03 000:00:46 000:00:11 S 10500 2.100000e+01 2.000000e-03 000:00:49 000:00:09 S 11000 2.200000e+01 2.000000e-03 000:00:51 000:00:07 S 11500 2.300000e+01 2.000000e-03 000:00:53 000:00:04 S 12000 2.400000e+01 2.000000e-03 000:00:56 000:00:02 S 12500 2.500000e+01 2.000000e-03 000:00:58 000:00:00 S Normal finish, maximum time reached: 25.000000 * Timers (h:m:s): ----------------- Migrate global-scope data : 0:0:0 Initial conditions : 0:0:0 Total runtime : 0:0:58 [Partition 0][Node 0] End of program
Results
The rationale for these runs is to integrate the system in time, extract the time evolution of various statistics, and then test several different closure hypotheses among the statistics. Example gnuplot commands to plot and test some closure ideas are:
# vim: filetype=gnuplot: # Nomenclature: # ------------- # V - instantaneous specific volume # R - instantaneous density # v = V - <V>, fluctuating specific volume # r = R - <R>, fluctuating density # b = -<rv>, density-specific-volume covariance # B - as postfix means the Bousinessq approximation # nm - no-mix value # <y^2> / <y^2>nm = <r^2> / <r^2>nm = b / bnm, low At plot "stat.txt" u 2:($23/($3*(1.0-$3))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($28/((0.0101**2)/((1.0+0.0101)**3)*$3*(1.0-$3))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$29/($3*(1.0+$3)*0.0101*0.0101/(1.0+0.0101))) w l t "b/bnm" # no-mix plot "stat.txt" u 2:($24/($4*(1.0-$4))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($30/((0.0101**2)/((1.0+0.0101)**3)*$4*(1.0-$4))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$31/($4*(1.0+$4)*0.0101*0.0101/(1.0+0.0101))) w l t "b/bnm" plot "stat.txt" u 2:($25/($5*(1.0-$5))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($32/((0.0101**2)/((1.0+0.0101)**3)*$5*(1.0-$5))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$33/($5*(1.0+$5)*0.0101*0.0101/(1.0+0.0101))) w l t "b/bnm" # uniform plot "stat.txt" u 2:($26/($6*(1.0-$6))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($34/((0.0101**2)/((1.0+0.0101)**3)*$6*(1.0-$6))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$35/($6*(1.0+$6)*0.0101*0.0101/(1.0+0.0101))) w l t "b/bnm" plot "stat.txt" u 2:($27/($7*(1.0-$7))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($36/((0.0101**2)/((1.0+0.0101)**3)*$7*(1.0-$7))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$37/($7*(1.0+$7)*0.0101*0.0101/(1.0+0.0101))) w l t "b/bnm" # almost mixed # <y^2> / <y^2>nm = <r^2> / <r^2>nm = b / bnm, high At plot "stat.txt" u 2:($23/($3*(1.0-$3))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($28/((9.0**2)/((1.0+9.0)**3)*$3*(1.0-$3))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$29/($3*(1.0+$3)*9.0*9.0/(1.0+9.0))) w l t "b/bnm" # no-mix plot "stat.txt" u 2:($24/($4*(1.0-$4))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($30/((9.0**2)/((1.0+9.0)**3)*$4*(1.0-$4))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$31/($4*(1.0+$4)*9.0*9.0/(1.0+9.0))) w l t "b/bnm" plot "stat.txt" u 2:($25/($5*(1.0-$5))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($32/((9.0**2)/((1.0+9.0)**3)*$5*(1.0-$5))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$33/($5*(1.0+$5)*9.0*9.0/(1.0+9.0))) w l t "b/bnm" # uniform plot "stat.txt" u 2:($26/($6*(1.0-$6))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($34/((9.0**2)/((1.0+9.0)**3)*$6*(1.0-$6))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$35/($6*(1.0+$6)*9.0*9.0/(1.0+9.0))) w l t "b/bnm" plot "stat.txt" u 2:($27/($7*(1.0-$7))) w l t "<y^2>/<y^2>nm", "stat.txt" u 2:($36/((9.0**2)/((1.0+9.0)**3)*$7*(1.0-$7))) w l t "<r^2>/<r^2>nm", "stat.txt" u 2:(-$37/($7*(1.0+$7)*9.0*9.0/(1.0+9.0))) w l t "b/bnm" # almost mixed # consistency test: <Rv^2> = b<V> = b(1+b)/<R>, both At plot "stat.txt" u 2:18 w l t "<Rv^2>", "stat.txt" u 2:(-$29*$13) w p ps 1.0 pt 7 t "b<V>", "stat.txt" u 2:(-$29*(1.0-$29)/$8) w p ps 1.0 pt 7 t "b(1+b)/<R>" # no-mix plot "stat.txt" u 2:19 w l t "<Rv^2>", "stat.txt" u 2:(-$31*$14) w p ps 1.0 pt 7 t "b<V>", "stat.txt" u 2:(-$31*(1.0-$31)/$9) w p ps 1.0 pt 7 t "b(1+b)/<R>" plot "stat.txt" u 2:20 w l t "<Rv^2>", "stat.txt" u 2:(-$33*$15) w p ps 1.0 pt 7 t "b<V>", "stat.txt" u 2:(-$33*(1.0-$33)/$10) w p ps 1.0 pt 7 t "b(1+b)/<R>" # uniform plot "stat.txt" u 2:21 w l t "<Rv^2>", "stat.txt" u 2:(-$35*$16) w p ps 1.0 pt 7 t "b<V>", "stat.txt" u 2:(-$35*(1.0-$35)/$11) w p ps 1.0 pt 7 t "b(1+b)/<R>" plot "stat.txt" u 2:22 w l t "<Rv^2>", "stat.txt" u 2:(-$37*$17) w p ps 1.0 pt 7 t "b<V>", "stat.txt" u 2:(-$37*(1.0-$37)/$12) w p ps 1.0 pt 7 t "b(1+b)/<R>" # almost mixed # b = <R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho1rho2/<R>^2 - 1) ] where <v^2>nm = (r/rho_2)^2 <y^2>nm, <y^2>nm = <Y>(1-<Y>), low At plot "stat.txt" u 2:(-$29) w l t "b", "stat.txt" u 2:($8*$8*$38*(1.0+$38/0.0101/0.0101/$3/(1.0+$3)*(0.99/$8/$8-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho_2/<R>^2 - 1) ]" # no-mix plot "stat.txt" u 2:(-$31) w l t "b", "stat.txt" u 2:($9*$9*$39*(1.0+$39/0.0101/0.0101/$4/(1.0+$4)*(0.99/$9/$9-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho1_rho2/<R>^2 - 1) ]" plot "stat.txt" u 2:(-$33) w l t "b", "stat.txt" u 2:($10*$10*$40*(1.0+$40/0.0101/0.0101/$5/(1.0+$5)*(0.99/$10/$10-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho2/<R>^2 - 1) ]" # uniform plot "stat.txt" u 2:(-$35) w l t "b", "stat.txt" u 2:($11*$11*$41*(1.0+$41/0.0101/0.0101/$6/(1.0+$6)*(0.99/$11/$11-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho_2/<R>^2 - 1) ]" plot "stat.txt" u 2:(-$37) w l t "b", "stat.txt" u 2:($12*$12*$42*(1.0+$42/0.0101/0.0101/$7/(1.0+$7)*(0.99/$12/$12-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho2/<R>^2 - 1) ]" # almost mixed # b = <R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho1rho2/<R>^2 - 1) ] where <v^2>nm = (r/rho_2)^2 <y^2>nm, <y^2>nm = <Y>(1-<Y>), high At plot "stat.txt" u 2:(-$29) w l t "b", "stat.txt" u 2:($8*$8*$38*(1.0+$38/9.0/9.0/$3/(1.0+$3)*(0.1/$8/$8-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho_2/<R>^2 - 1) ]" # no-mix plot "stat.txt" u 2:(-$31) w l t "b", "stat.txt" u 2:($9*$9*$39*(1.0+$39/9.0/9.0/$4/(1.0+$4)*(0.1/$9/$9-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho_2/<R>^2 - 1) ]" plot "stat.txt" u 2:(-$33) w l t "b", "stat.txt" u 2:($10*$10*$40*(1.0+$40/9.0/9.0/$5/(1.0+$5)*(0.1/$10/$10-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho_2/<R>^2 - 1) ]" # uniform plot "stat.txt" u 2:(-$35) w l t "b", "stat.txt" u 2:($11*$11*$41*(1.0+$41/9.0/9.0/$6/(1.0+$6)*(0.1/$11/$11-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho_2/<R>^2 - 1) ]" plot "stat.txt" u 2:(-$37) w l t "b", "stat.txt" u 2:($12*$12*$42*(1.0+$42/9.0/9.0/$7/(1.0+$7)*(0.1/$12/$12-1.0))) w p ps 1.0 pt 7 t "<R>^2<v^2> [1 + <v^2>/<v^2>nm * (rho_1rho_2/<R>^2 - 1) ]" # almost mixed # <r^2> = <R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho1^3rho2/<R>^4 - 1) ] where <v^2>nm = (r/rho_2)^2 <y^2>nm, <y^2>nm = <Y>(1-<Y>), low At plot "stat.txt" u 2:(-$29) w l t "b", "stat.txt" u 2:($8**4.0*$38*(1.0+$38/0.0101/0.0101/$3/(1.0+$3)*(0.99**2.0/$8**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^3rho_2/<R>^4 - 1) ]" # no-mix plot "stat.txt" u 2:(-$31) w l t "b", "stat.txt" u 2:($9**4.0*$39*(1.0+$39/0.0101/0.0101/$4/(1.0+$4)*(0.99**2.0/$9**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^3rho_2/<R>^4 - 1) ]" plot "stat.txt" u 2:(-$33) w l t "b", "stat.txt" u 2:($10**4.0*$40*(1.0+$40/0.0101/0.0101/$5/(1.0+$5)*(0.99**2.0/$10**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^3rho_2/<R>^4 - 1) ]" # uniform plot "stat.txt" u 2:(-$35) w l t "b", "stat.txt" u 2:($11**4.0*$41*(1.0+$41/0.0101/0.0101/$6/(1.0+$6)*(0.99**2.0/$11**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^3rho_2/<R>^4 - 1) ]" plot "stat.txt" u 2:(-$37) w l t "b", "stat.txt" u 2:($12**4.0*$42*(1.0+$42/0.0101/0.0101/$7/(1.0+$7)*(0.99**2.0/$12**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^3rho_2/<R>^4 - 1) ]" # almost mixed # <r^2> = <R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho1^3rho2/<R>^4 - 1) ] where <v^2>nm = (r/rho_2)^2 <y^2>nm, <y^2>nm = <Y>(1-<Y>), high At plot "stat.txt" u 2:(-$29) w l t "b", "stat.txt" u 2:($8**4.0*$38*(1.0+$38/9.0/9.0/$3/(1.0+$3)*(0.1**2.0/$8**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^2rho_2^2/<R>^4 - 1) ]" # no-mix plot "stat.txt" u 2:(-$31) w l t "b", "stat.txt" u 2:($9**4.0*$39*(1.0+$39/9.0/9.0/$4/(1.0+$4)*(0.1**2.0/$9**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^2rho_2^2/<R>^4 - 1) ]" plot "stat.txt" u 2:(-$33) w l t "b", "stat.txt" u 2:($10**4.0*$40*(1.0+$40/9.0/9.0/$5/(1.0+$5)*(0.1**2.0/$10**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^2rho_2^2/<R>^4 - 1) ]" # uniform plot "stat.txt" u 2:(-$35) w l t "b", "stat.txt" u 2:($11**4.0*$41*(1.0+$41/9.0/9.0/$6/(1.0+$6)*(0.1**2.0/$11**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^2rho_2^2/<R>^4 - 1) ]" plot "stat.txt" u 2:(-$37) w l t "b", "stat.txt" u 2:($12**4.0*$42*(1.0+$42/9.0/9.0/$7/(1.0+$7)*(0.1**2.0/$12**4.0-1.0))) w p ps 1.0 pt 7 t "<R>^4<v^2> [1 + <v^2>/<v^2>nm * (rho_1^2rho_2^2/<R>^4 - 1) ]" # almost mixed