Merge branch 'develop' into always-exceptions

# Conflicts:
#	src/library.cpp
This commit is contained in:
Axel Kohlmeyer
2023-08-03 16:46:32 -04:00
626 changed files with 8068 additions and 4460 deletions

View File

@ -1,23 +1,23 @@
# SRD diffusion demo - ellipsoids
units lj
atom_style ellipsoid
atom_modify first big
dimension 2
units lj
atom_style ellipsoid
atom_modify first big
dimension 2
# create big ellipsoidal particles
lattice sq 0.14
region box block 0 10 0 10 -0.5 0.5
create_box 2 box
create_atoms 1 region box
lattice sq 0.14
region box block 0 10 0 10 -0.5 0.5
create_box 2 box
create_atoms 1 region box
set type 1 mass 1.0
set type 1 shape 3.0 1.0 1.0
group big type 1
set group big quat/random 29898
set type 1 mass 1.0
set type 1 shape 3.0 1.0 1.0
group big type 1
set group big quat/random 29898
velocity big create 1.44 87287 loop geom
velocity big create 1.44 87287 loop geom
# equilibrate big particles
@ -26,57 +26,59 @@ pair_coeff 1 1 1.0 1.0 1 1 1 1 1 1
pair_coeff 1 2 1.0 1.0 1 1 1 1 1 1 0.0
pair_coeff 2 2 1.0 1.0 1 1 1 1 1 1 0.0
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
fix 1 big nve/asphere
fix 2 all enforce2d
fix 1 big nve/asphere
fix 2 all enforce2d
compute rot big temp/asphere
compute rot big temp/asphere
compute 0 all property/atom quatw quati quatj quatk shapex shapey shapez
#dump 1 all custom 10 dump.ellipsoid.equil id type x y z &
# quatw quati quatj quatk
#dump 1 all custom 10 dump.ellipsoid.equil id type x y z c_0[*]
#dump_modify 1 colname c_0[1] quatw colname c_0[2] quati colname c_0[3] quatj colname c_0[4] quatk &
# colname c_0[5] shapex colname c_0[6] shapey colname c_0[7] shapez
thermo_style custom step temp c_rot epair etotal press
thermo 100
thermo_style custom step temp c_rot epair etotal press
thermo 100
run 1000
run 1000
#undump 1
unfix 1
unfix 2
#undump 1
unfix 1
unfix 2
# add small particles as hi density lattice
region plane block INF INF INF INF -0.001 0.001 units box
lattice sq 120.0
create_atoms 2 region plane
region plane block INF INF INF INF -0.001 0.001 units box
lattice sq 120.0
create_atoms 2 region plane
set type 2 mass 0.01
group small type 2
velocity small create 1.0 593849 loop geom
set type 2 mass 0.01
group small type 2
velocity small create 1.0 593849 loop geom
# delete overlaps
# must set 1-2 cutoff to non-zero value
pair_style lj/cut 2.5
pair_coeff 1 1 1.0 1.0
pair_coeff 2 2 0.0 1.0 0.0
pair_coeff 1 2 0.0 1.0 2.0
pair_style lj/cut 2.5
pair_coeff 1 1 1.0 1.0
pair_coeff 2 2 0.0 1.0 0.0
pair_coeff 1 2 0.0 1.0 2.0
neigh_modify one 10000
neigh_modify one 10000
delete_atoms overlap 1.6 small big
delete_atoms overlap 1.6 small big
# SRD run
reset_timestep 0
reset_timestep 0
neighbor 0.3 multi
neigh_modify delay 0 every 1 check yes
neighbor 0.3 multi
neigh_modify delay 0 every 1 check yes
comm_modify mode multi group big vel yes
neigh_modify include big
comm_modify mode multi group big vel yes
neigh_modify include big
# no pairwise interactions with small particles
@ -88,29 +90,30 @@ pair_coeff 2 2 1.0 1.0 1 1 1 1 1 1 0.0
# use fix SRD to push small particles out from inside big ones
# if comment out, big particles won't see SRD particles
timestep 0.0005
timestep 0.0005
fix 1 big nve/asphere
fix 2 small srd 20 big 1.0 0.25 49894 shift yes 54979 &
collision noslip search 0.2 inside ignore exact no bounce 50
fix 3 all enforce2d
fix 1 big nve/asphere
fix 2 small srd 20 big 1.0 0.25 49894 shift yes 54979 &
collision noslip search 0.2 inside ignore exact no bounce 50
fix 3 all enforce2d
# diagnostics
compute tbig big temp/asphere
variable pebig equal pe*atoms/count(big)
variable ebig equal etotal*atoms/count(big)
thermo_style custom step temp c_rot f_2[9] etotal v_pebig v_ebig press &
f_2[1] f_2[2] f_2[3] f_2[4] f_2[5] &
f_2[6] f_2[7] f_2[8] f_2[9] f_2[10] f_2[11] f_2[12]
compute tbig big temp/asphere
variable pebig equal pe*atoms/count(big)
variable ebig equal etotal*atoms/count(big)
thermo_style custom step temp c_rot f_2[9] etotal v_pebig v_ebig press &
f_2[1] f_2[2] f_2[3] f_2[4] f_2[5] &
f_2[6] f_2[7] f_2[8] f_2[9] f_2[10] f_2[11] f_2[12]
thermo_modify temp tbig
thermo 1000
thermo_modify temp tbig
thermo 1000
#dump 1 all custom 1000 dump.ellipsoid id type x y z &
# quatw quati quatj quatk
#dump 1 all custom 1000 dump.ellipsoid id type x y z c_0[*]
#dump_modify 1 colname c_0[1] quatw colname c_0[2] quati colname c_0[3] quatj colname c_0[4] quatk &
# colname c_0[5] shapex colname c_0[6] shapey colname c_0[7] shapez
#dump 1 all image 1000 image.*.jpg type type zoom 1.6
#dump_modify 1 pad 6 adiam 1 1 adiam 2 0.2
#dump 2 all image 1000 image.*.jpg type type zoom 1.6
#dump_modify 2 pad 6 adiam 1 1 adiam 2 0.2
run 100000
run 100000

View File

@ -1,23 +1,23 @@
# SRD viscosity demo - ellipsoids
units lj
atom_style ellipsoid
atom_modify first big
dimension 2
units lj
atom_style ellipsoid
atom_modify first big
dimension 2
# create big ellipsoidal particles
lattice sq 0.14
region box block 0 10 0 10 -0.5 0.5
create_box 2 box
create_atoms 1 region box
lattice sq 0.14
region box block 0 10 0 10 -0.5 0.5
create_box 2 box
create_atoms 1 region box
set type 1 mass 1.0
set type 1 shape 3.0 1.0 1.0
group big type 1
set group big quat/random 29898
set type 1 mass 1.0
set type 1 shape 3.0 1.0 1.0
group big type 1
set group big quat/random 29898
velocity big create 1.44 87287 loop geom
velocity big create 1.44 87287 loop geom
# equilibrate big particles
@ -26,57 +26,59 @@ pair_coeff 1 1 1.0 1.0 1 1 1 1 1 1
pair_coeff 1 2 1.0 1.0 1 1 1 1 1 1 0.0
pair_coeff 2 2 1.0 1.0 1 1 1 1 1 1 0.0
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
fix 1 big nve/asphere
fix 2 all enforce2d
fix 1 big nve/asphere
fix 2 all enforce2d
compute rot big temp/asphere
compute rot big temp/asphere
compute 0 all property/atom quatw quati quatj quatk shapex shapey shapez
#dump 1 all custom 10 dump.ellipsoid.equil id type x y z &
# quatw quati quatj quatk
#dump 1 all custom 10 dump.ellipsoid.equil id type x y z c_0[*]
#dump_modify 1 colname c_0[1] quatw colname c_0[2] quati colname c_0[3] quatj colname c_0[4] quatk &
# colname c_0[5] shapex colname c_0[6] shapey colname c_0[7] shapez
thermo_style custom step temp c_rot epair etotal press
thermo 100
thermo_style custom step temp c_rot epair etotal press
thermo 100
run 1000
run 1000
#undump 1
unfix 1
unfix 2
#undump 1
unfix 1
unfix 2
# add small particles as hi density lattice
region plane block INF INF INF INF -0.001 0.001 units box
lattice sq 120.0
create_atoms 2 region plane
region plane block INF INF INF INF -0.001 0.001 units box
lattice sq 120.0
create_atoms 2 region plane
set type 2 mass 0.01
group small type 2
velocity small create 1.0 593849 loop geom
set type 2 mass 0.01
group small type 2
velocity small create 1.0 593849 loop geom
# delete overlaps
# must set 1-2 cutoff to non-zero value
pair_style lj/cut 2.5
pair_coeff 1 1 1.0 1.0
pair_coeff 2 2 0.0 1.0 0.0
pair_coeff 1 2 0.0 1.0 2.0
pair_style lj/cut 2.5
pair_coeff 1 1 1.0 1.0
pair_coeff 2 2 0.0 1.0 0.0
pair_coeff 1 2 0.0 1.0 2.0
neigh_modify one 10000
neigh_modify one 10000
delete_atoms overlap 1.6 small big
delete_atoms overlap 1.6 small big
# SRD run
reset_timestep 0
reset_timestep 0
neighbor 0.3 multi
neigh_modify delay 0 every 1 check yes
neighbor 0.3 multi
neigh_modify delay 0 every 1 check yes
comm_modify mode multi group big vel yes
neigh_modify include big
comm_modify mode multi group big vel yes
neigh_modify include big
# no pairwise interactions with small particles
@ -88,30 +90,31 @@ pair_coeff 2 2 1.0 1.0 1 1 1 1 1 1 0.0
# use fix SRD to push small particles out from inside big ones
# if comment out, big particles won't see SRD particles
timestep 0.0005
timestep 0.0005
fix 1 big nve/asphere
fix 2 small srd 20 big 1.0 0.25 49894 shift yes 54979 &
fix 1 big nve/asphere
fix 2 small srd 20 big 1.0 0.25 49894 shift yes 54979 &
collision noslip search 0.2 inside ignore exact no bounce 50
fix 3 small viscosity 20 x y 50
fix 4 all enforce2d
fix 3 small viscosity 20 x y 50
fix 4 all enforce2d
# diagnostics
compute tbig big temp/asphere
variable pebig equal pe*atoms/count(big)
variable ebig equal etotal*atoms/count(big)
thermo_style custom step temp c_rot f_2[9] etotal v_pebig v_ebig press &
f_2[1] f_2[2] f_2[3] f_2[4] f_2[5] &
f_2[6] f_2[7] f_2[8] f_2[9] f_2[10] f_2[11] f_2[12]
compute tbig big temp/asphere
variable pebig equal pe*atoms/count(big)
variable ebig equal etotal*atoms/count(big)
thermo_style custom step temp c_rot f_2[9] etotal v_pebig v_ebig press &
f_2[1] f_2[2] f_2[3] f_2[4] f_2[5] &
f_2[6] f_2[7] f_2[8] f_2[9] f_2[10] f_2[11] f_2[12]
thermo_modify temp tbig
thermo 1000
thermo_modify temp tbig
thermo 1000
#dump 1 all custom 500 dump.ellipsoid.mp id type x y z &
# quatw quati quatj quatk
#dump 1 all custom 500 dump.ellipsoid id type x y z c_0[*]
#dump_modify 1 colname c_0[1] quatw colname c_0[2] quati colname c_0[3] quatj colname c_0[4] quatk &
# colname c_0[5] shapex colname c_0[6] shapey colname c_0[7] shapez
#dump 1 all image 500 image.*.jpg type type zoom 1.6
#dump_modify 1 pad 6 adiam 1 1 adiam 2 0.2
#dump 1 all image 500 image.*.jpg type type zoom 1.6
#dump_modify 1 pad 6 adiam 1 1 adiam 2 0.2
run 50000
run 50000

View File

@ -239,7 +239,7 @@ struct _liblammpsplugin {
void (*free)(void *);
void (*is_running)(void *);
int (*is_running)(void *);
void (*force_timeout)(void *);
int (*has_error)(void *);

View File

@ -30,7 +30,7 @@ PROGRAM f_driver
REAL (kind=8), ALLOCATABLE :: x(:)
REAL (kind=8), PARAMETER :: epsilon=0.1
CHARACTER (len=64) :: arg
CHARACTER (len=1024) :: arg
CHARACTER (len=1024) :: line
! setup MPI and various communicators

View File

@ -1,11 +1,11 @@
# This test evaluates Hardy fields for a small block of LJ material
# that's undergoing translation in the x-direction through the
# that's undergoing translation in the x-direction through the
# periodic boundary.
echo both
log log.flying_cube
log log.flying_cube
units real
atom_style atomic
variable lattice_constant equal 5.256227487
variable lattice_constant equal 5.256227487
variable c equal 6
variable L equal (${lattice_constant}*$c)
variable V equal $L*$L*$L
@ -23,15 +23,15 @@ create_box 1 box
create_atoms 1 region box
mass 1 $m
group all region box
pair_style lj/cut 13.5
pair_coeff 1 1 0.238 3.405
pair_style lj/cut 13.5
pair_coeff 1 1 0.238 3.405
# define region
neighbor 1.0 bin
neigh_modify delay 0 every 200 check no
# neigh_modify delay 0 every $n
reset_timestep 0
thermo $s
compute_modify thermo_temp extra 0
compute_modify thermo_temp extra/dof 0
thermo_style custom step temp pxx pyy pzz etotal lx ly lz
timestep ${dt}
min_modify line quadratic
@ -46,10 +46,10 @@ variable zHI equal zhi
# region BOX block ${xLO} ${xHI} ${yLO} ${yHI} ${zLO} ${zHI} units box
region BOX block ${xLO} ${xHI} -0.5 31.7 -0.5 31.7 units box
fix ATC all atc field
fix_modify ATC mesh create 1 1 1 box p p p
fix_modify ATC mesh create 1 1 1 box p p p
fix_modify ATC fields add stress velocity displacement
fix_modify ATC fields add temperature kinetic_temperature
fix_modify ATC fields add internal_energy energy
fix_modify ATC fields add internal_energy energy
fix_modify ATC output flying_cubeFE $s text binary tensor_components
# fix_modify ATC atomic_output flying_cubeMD $m text binary tensor_components
dump dumpfc all atom $s flying_cube.dmp

View File

@ -1,34 +1,34 @@
# This test compares Hardy and Cauchy-Born metrics of energy density and stress
# for an Ar system undergoing dynamics at a finite temperature.
echo both
log ftcb_constV.log
log ftcb_constV.log
units real
atom_style atomic
variable lattice_constant equal 5.256227487
variable lattice_constant equal 5.256227487
variable c equal 6
variable L equal (${lattice_constant}*$c)
variable V equal $L*$L*$L
print "Volume : $V"
variable Ti equal 0.0001
variable Ti equal 0.0001
variable dT equal 0.2
variable n equal 10000
variable m equal 200
variable dt equal 0.1
# create system
lattice fcc ${lattice_constant}
lattice fcc ${lattice_constant}
region box block 0 $c 0 $c 0 $c
boundary p p p
pair_style lj/cut 13.5
pair_style lj/cut 13.5
read_data ftcb_constV_setup.init
mass 1 39.95
pair_coeff 1 1 0.238 3.405
pair_coeff 1 1 0.238 3.405
# define region
neighbor 1.0 bin
neigh_modify delay 0 every 2000 check no
reset_timestep 0
thermo $m
compute_modify thermo_temp extra 0
variable nrepeat equal $n/$m
compute_modify thermo_temp extra/dof 0
variable nrepeat equal $n/$m
variable sxx equal -pxx
variable intenergy equal etotal
fix TDAVE all ave/time $m ${nrepeat} $n c_thermo_temp v_sxx v_intenergy ave one file ftcb_constV.profile
@ -37,7 +37,7 @@ timestep ${dt}
# filtered
fix ATCFILT all atc field Ar_CauchyBorn.mat
#fix_modify ATCFILT reset_atomic_reference_positions ftcb_constV_setup.data
fix_modify ATCFILT mesh create 1 1 1 box p p p
fix_modify ATCFILT mesh create 1 1 1 box p p p
fix_modify ATCFILT fields none
fix_modify ATCFILT fields add internal_energy stress cauchy_born_energy cauchy_born_stress
fix_modify ATCFILT fields add temperature displacement
@ -51,7 +51,7 @@ fix_modify ATCFILT filter type step
# not filtered
fix ATC all atc field Ar_CauchyBorn.mat
#fix_modify ATC reset_atomic_reference_positions ftcb_constV_setup.data
fix_modify ATC mesh create 1 1 1 box p p p
fix_modify ATC mesh create 1 1 1 box p p p
fix_modify ATC fields none
fix_modify ATC fields add internal_energy stress cauchy_born_energy cauchy_born_stress
fix_modify ATC fields add temperature displacement
@ -59,7 +59,7 @@ fix_modify ATC gradients add displacement
fix_modify ATC output ftcb_constVFE $m text binary tensor_components
velocity all create ${Ti} 102486 mom yes rot yes dist gaussian
# step iu T,V space
variable i loop 1
variable i loop 1
label loop_i
print ">>> step $i, T: ${Ti}, V: $V"
variable Tf equal ${Ti}+${dT}
@ -69,6 +69,6 @@ label loop_i
fix NVT all nvt temp ${Tf} ${Tf} 20 drag 0.5 tchain 1
run $n
unfix NVT
variable Ti equal ${Tf}
variable Ti equal ${Tf}
next i
jump in.ftcb_constV loop_i

View File

@ -5,50 +5,50 @@
# out # to the FEM on the right. Insufficient time is captured to reach the
# linear # steady state, but heat crossing both boundaries should be observed.
log bar1d_fluids.log
units real
atom_style atomic
units real
atom_style atomic
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
read_data temp.init
#read_restart post_eq.dat
#read_restart post_eq.dat
lattice fcc 5.405 origin 0.25 0.25 0.25
# create atoms, NOTE commented out for restart
#region mdRegion block -7 7 -3 3 -3 3
#boundary f p p
#create_box 2 mdRegion
#create_atoms 1 region mdRegion
#mass * 39.95
#region mdRegion block -7 7 -3 3 -3 3
#boundary f p p
#create_box 2 mdRegion
#create_atoms 1 region mdRegion
#mass * 39.95
# specify interal/ghost atoms
#region mdInternal block -6 6 -3 3 -3 3
#region mdGhost block -6 6 -3 3 -3 3 side out
#region mdInternal block -6 6 -3 3 -3 3
#region mdGhost block -6 6 -3 3 -3 3 side out
#create_atoms 1 region mdGhost
#lattice fcc 6.5 origin 0.25 0.25 0.25
#create_atoms 2 region mdInternal
group internal type 2
group ghost type 1
#create_atoms 2 region mdInternal
group internal type 2
group ghost type 1
# velocities have Vcm = 0, NOTE next four lines commented out for restart
#velocity internal create 100. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff * * .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#write_restart tinit.dat
#velocity internal create 100. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff * * .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#write_restart tinit.dat
# ID group atc PhysicsType ParameterFile
lattice fcc 5.405 origin 0.25 0.25 0.25
region atcRegion block -6.25 6.25 -3 3 -3 3
region atcRegion block -6.25 6.25 -3 3 -3 3
fix AtC internal atc thermal Ar_thermal.mat
# ID part keywords nx ny nz region periodicity
fix_modify AtC mesh create 6 1 1 atcRegion f p p
fix_modify AtC time_integration fractional_step
fix_modify AtC atom_element_map eulerian 10 # NOTE this introduces a minor amount of time integration error tracking the atomic temperature
fix_modify AtC atom_element_map eulerian 10 # NOTE this introduces a minor amount of time integration error tracking the atomic temperature
fix_modify AtC internal_quadrature off
#fix_modify AtC mass_matrix fe
#fix_modify AtC mass_matrix fe
# fix a temperature
fix_modify AtC fix temperature all 100.
fix_modify AtC fix temperature all 100.
# add nodesets
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
fix_modify AtC mesh create_nodeset lbc -6.3 -6.2 -INF INF -INF INF
@ -56,43 +56,43 @@ fix_modify AtC mesh create_nodeset rbc 6.2 6.3 -INF INF -INF INF
#fix_modify AtC fix temperature lbc 120.
#fix_modify AtC fix temperature rbc 100.
#used for restarting
#fix_modify AtC consistent_fe_initialization on
#fix_modify AtC consistent_fe_initialization on
#fix_modify AtC initial temperature all 100.
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
# equilibrate MD field
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
#fix_modify AtC output bar1d_fluids_eqFE 200 text binary
#dump D1 all atom 200 dump_eq.bar1d
timestep 2
thermo 200
run 2000
#write_restart post_eq.dat
timestep 2
thermo 200
run 2000
#write_restart post_eq.dat
# change thermostat to hoover
reset_timestep 0
fix_modify AtC reset_time
fix_modify AtC reset_time
fix_modify AtC unfix temperature all
variable deltaT equal 2.*2000.
fix_modify AtC fix temperature lbc temporal_ramp 100. 120. ${deltaT}
#fix_modify AtC fix temperature lbc 100.
#fix_modify AtC fix temperature lbc 100.
fix_modify AtC fix temperature rbc 100.
fix_modify AtC control thermal flux no_boundary
fix_modify AtC control tolerance 1.e-14
fix_modify AtC control localized_lambda on
fix_modify AtC control tolerance 1.e-14
fix_modify AtC control localized_lambda on
fix_modify AtC filter type exponential
fix_modify AtC filter scale 1000.0
fix_modify AtC filter on
# output commands
fix_modify AtC output bar1d_fluidsFE 100 text
#undump D1
#dump D2 all atom 200 dump.bar1d
#dump D2 all custom 200 dump.bar1d id type xs ys zs vx vy vz
#dump D2 all atom 200 dump.bar1d
#dump D2 all custom 200 dump.bar1d id type xs ys zs vx vy vz
# set-up non-equilibrium IC
thermo 100
run 2000
thermo 100
run 2000
# run non-equilibrium condition
fix_modify AtC fix temperature lbc 120.
run 4000
fix_modify AtC fix temperature lbc 120.
run 4000

View File

@ -1,10 +1,10 @@
echo both
units real
units real
atom_style full # charge
### NOTE p p p is required for both lammps & atc since periodic images are employed in ExtrinsicModelElectrostatic::correct_electrostatic_forces
########### BEGIN PARAMETERS ####################################
variable T equal 300
variable a equal 1.0 # 5.719025032
variable T equal 300
variable a equal 1.0 # 5.719025032
variable i equal 10 # thermo
variable f equal 50 # neighbor & conc interval
variable s equal 100 # 10 # 100 # output
@ -12,18 +12,18 @@ variable n equal 800 # 20 # 200 # 300 # duration
variable x equal 4 # 40 # 2 # 4 # exchanges
variable e equal 100. # 1. # 100. # 10. # energy
variable h equal 5 # nelems
variable dt equal 4. #1. # 0 # 4.0
variable dt equal 4. #1. # 0 # 4.0
############## END PARAMETERS #################################
dimension 3
boundary p p p
pair_style lj/cut/coul/cut 13.0
lattice sc $a
dimension 3
boundary p p p
pair_style lj/cut/coul/cut 13.0
lattice sc $a
read_data concentration_init.data
atom_modify sort 0 1
mass * 39.948
pair_coeff * * 0.2381 3.405
pair_coeff 1 * 0.4 3.405
dielectric 80.0
dielectric 80.0
variable xlo equal xlo
variable xhi equal xhi
variable xmid equal 0.5*(${xhi}+${xlo})
@ -37,9 +37,9 @@ region FLUID block ${xlo} ${xhi} ${ylo} ${yhi} ${zlo} ${zhi} units box
region R block ${xmid} ${xhi} ${ylo} ${yhi} ${zlo} ${zhi} units box
group SOLID type 1
variable xdof equal 3*count(SOLID)
compute_modify thermo_temp extra ${xdof}
#set group SOLID charge 0
group NEUTRAL type 2
compute_modify thermo_temp extra/dof ${xdof}
#set group SOLID charge 0
group NEUTRAL type 2
group FLUID type 2 3 4
group NION type 3
group PION type 4
@ -72,16 +72,16 @@ compute PMIN PION reduce min x
compute PMAX PION reduce max x
compute NMIN NION reduce min x
compute NMAX NION reduce max x
timestep ${dt}
timestep ${dt}
#neighbor 13 bin
neigh_modify every $i check no
# coulombic interactions
fix ATC FLUID atc species_electrostatic Ar_electrostatic.mat
fix ATC FLUID atc species_electrostatic Ar_electrostatic.mat
#fix_modify ATC parallel_consistency off
fix_modify ATC extrinsic short_range off
#fix_modify ATC boundary SOLID
fix_modify ATC atom_element_map eulerian $i
fix_modify ATC internal_quadrature off
fix_modify ATC internal_quadrature off
fix_modify ATC consistent_fe_initialization on
fix_modify ATC filter type exponential
fix_modify ATC filter scale $i
@ -97,31 +97,31 @@ fix_modify ATC mesh create_elementset BOX BOX
fix_modify ATC mass_matrix fe
fix_modify ATC include atomic_charge
fix_modify ATC add_species N type 3
fix_modify ATC add_species P type 4
fix_modify ATC add_species P type 4
fix_modify ATC add_species Nt type 5
fix_modify ATC add_species Pt type 6
fix_modify ATC add_species Pt type 6
## CC
fix_modify ATC control concentration N R 0.001 Nt # deletions R
fix_modify ATC control concentration N frequency $f
fix_modify ATC control concentration N max_energy $e
fix_modify ATC control concentration N frequency $f
fix_modify ATC control concentration N max_energy $e
fix_modify ATC control concentration N max_attempts 100
fix_modify ATC control concentration N max_exchanges $x
#-
fix_modify ATC control concentration P R 0.002 Pt # addtions R
fix_modify ATC control concentration P frequency $f
fix_modify ATC control concentration P frequency $f
fix_modify ATC control concentration P max_energy $e
fix_modify ATC control concentration P max_attempts 100
fix_modify ATC control concentration P max_attempts 100
fix_modify ATC control concentration P max_exchanges $x
fix_modify ATC output volume_integral all mass_density
fix_modify ATC output volume_integral all charge_density
# data reduction
fix PP all atc field
#fix_modify PP add_species ions type 3 4
#fix_modify PP add_species ions type 3 4
#fix_modify PP add_species IONS type 3 4 1 2
fix_modify PP add_species n type 3
fix_modify PP add_species p type 4
fix_modify PP add_species s type 1
fix_modify PP add_species o type 2
fix_modify PP add_species n type 3
fix_modify PP add_species p type 4
fix_modify PP add_species s type 1
fix_modify PP add_species o type 2
fix_modify PP fields add species_concentration mass_density charge_density temperature
fix_modify PP fields add charge_flux species_flux
fix_modify PP output volume_integral all mass_density
@ -129,7 +129,7 @@ fix_modify PP output volume_integral all charge_density
fix_modify PP mesh create $h 1 1 BOX f p p
fix_modify PP atom_element_map eulerian $i
# output
thermo $i
thermo $i
variable dN equal f_ATC[1]
variable dP equal f_ATC[2]
variable Nc equal f_ATC[3]
@ -141,19 +141,19 @@ variable Nu equal f_ATC[7]
variable Pu equal f_ATC[8]
variable Ng equal f_ATC[9]
variable Pg equal f_ATC[10]
thermo_style custom step temp etotal &
thermo_style custom step temp etotal &
atoms v_F v_P v_N v_dP v_dN v_Pc v_Nc v_Pr v_Nr v_Pt v_Nt &
c_PMIN c_NMIN c_PMAX c_NMAX c_Q c_M c_Qf c_Mf v_Pu v_Nu v_Pg v_Ng
# NOTE this doesn't seem to work
#thermo_modify format 4 %5d format 5 %5d format 6 %5d format 7 %5d format 8 %5d
#thermo_modify format 4 %5d
#thermo_modify format 5 %5d
#thermo_modify format 6 %5d
#thermo_modify format 7 %5d
#thermo_modify format 8 %5d
#thermo_modify format 9 %5d
#thermo_modify format 10 %5d
#thermo_modify format 11 %5d
# NOTE this doesn't seem to work
#thermo_modify format 4 %5d format 5 %5d format 6 %5d format 7 %5d format 8 %5d
#thermo_modify format 4 %5d
#thermo_modify format 5 %5d
#thermo_modify format 6 %5d
#thermo_modify format 7 %5d
#thermo_modify format 8 %5d
#thermo_modify format 9 %5d
#thermo_modify format 10 %5d
#thermo_modify format 11 %5d
log concentration.log
fix_modify ATC output concentrationFE $s text binary
fix_modify PP output concentrationPP $s text

View File

@ -2,42 +2,42 @@
# This benchmark tests momentum flux as a BC to the MD region to generate shear flow.
# Currently heat will be generated as it will have no where to go.
log shear_flow.log
units real
atom_style atomic
units real
atom_style atomic
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
#read_data temp.init
read_data post_eq.init
read_data post_eq.init
lattice fcc 5.405 origin 0.25 0.25 0.25
# create atoms, NOTE commented out for restart
#region mdRegion block -8 8 -3 3 -3 3
#boundary f p p
#create_box 2 mdRegion
#create_atoms 1 region mdRegion
#mass * 39.95
#region mdRegion block -8 8 -3 3 -3 3
#boundary f p p
#create_box 2 mdRegion
#create_atoms 1 region mdRegion
#mass * 39.95
# specify interal/ghost atoms
#region mdInternal block -6 6 -3 3 -3 3
#region mdGhost block -6 6 -3 3 -3 3 side out
#region mdInternal block -6 6 -3 3 -3 3
#region mdGhost block -6 6 -3 3 -3 3 side out
#create_atoms 1 region mdGhost
#lattice fcc 6.5 origin 0.25 0.25 0.25
#create_atoms 2 region mdInternal
group internal type 2
group ghost type 1
#create_atoms 2 region mdInternal
group internal type 2
group ghost type 1
# velocities have Vcm = 0, NOTE next four lines commented out for restart
#velocity internal create 100. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff * * .238 3.405 13.5
#velocity internal create 100. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff * * .238 3.405 13.5
#write_restart tinit.dat
neighbor 5. bin
neigh_modify every 10 delay 0 check no
neighbor 5. bin
neigh_modify every 10 delay 0 check no
# equilibrate MD field
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe temp
thermo_modify format 1 %6i format 2 %7.2g
#fix NVT internal nvt temp 100 100 10 drag 0.2
@ -49,28 +49,28 @@ timestep 2
# ID group atc PhysicsType ParameterFile
lattice fcc 5.405 origin 0.25 0.25 0.25
region atcRegion block -6.25 6.25 -3 3 -3 3
region atcRegion block -6.25 6.25 -3 3 -3 3
fix AtC internal atc shear Ar_viscosity.mat
# ID part keywords nx ny nz region periodicity
fix_modify AtC mesh create 6 1 1 atcRegion f p p
#fix_modify AtC time_integration fractional_step
##fix_modify AtC atom_element_map eulerian 10 # NOTE this introduces a minor amount of time integration error tracking the atomic temperature
##fix_modify AtC atom_element_map eulerian 10 # NOTE this introduces a minor amount of time integration error tracking the atomic temperature
fix_modify AtC internal_quadrature off
#fix_modify AtC mass_matrix fe
#fix_modify AtC mass_matrix fe
# add nodesets
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
fix_modify AtC mesh create_elementset rbc 4.12 6.3 -INF INF -INF INF
#used for restarting
fix_modify AtC consistent_fe_initialization on
fix_modify AtC consistent_fe_initialization on
#fix_modify AtC initial velocity x all 0.
# turn on kinetostat
reset_timestep 0
fix_modify AtC reset_time
fix_modify AtC source velocity y rbc .0000001
fix_modify AtC reset_time
fix_modify AtC source velocity y rbc .0000001
fix_modify AtC control momentum flux no_boundary
fix_modify AtC control tolerance 1.e-14
fix_modify AtC control localized_lambda on
fix_modify AtC control tolerance 1.e-14
fix_modify AtC control localized_lambda on
#fix_modify AtC filter type exponential
#fix_modify AtC filter scale 1000.0
#fix_modify AtC filter on
@ -78,7 +78,7 @@ fix_modify AtC control localized_lambda on
# output commands
fix_modify AtC output shear_flowFE 100 text #binary
#undump D1
#dump D1 all custom 100 shear_flow.dmp id type xs ys zs vx vy vz
#dump D1 all custom 100 shear_flow.dmp id type xs ys zs vx vy vz
# set-up non-equilibrium IC
thermo 100
run 2000
thermo 100
run 2000

View File

@ -1,4 +1,4 @@
# This test calculates the J-integral for a sequence of different loads
# This test calculates the J-integral for a sequence of different loads
# for a system with a crack tip.
#
# loops around crack and one around undefected material
@ -7,17 +7,17 @@
# loop 2 - medium contour around left crack tip
# loop 3 - largest contour around left crack tip
# loop 4 - medium contour around right crack tip
#for Lennard-Jones potential for gold:
# surface energy (gamma) of 0.1599 eV/Å^2 = 2.56191 J/m2.
#for Lennard-Jones potential for gold:
# surface energy (gamma) of 0.1599 eV/Å^2 = 2.56191 J/m2.
# 2*gamma = 0.3198 eV/Å^2 = 5.12382 J/m2
# using the internal nktv2p conversion factor,
# using the internal nktv2p conversion factor,
# the units of eshelby_stress are eV/Å^3, and
# the resulting boundary integral is eV/Å since we also include thickness (a
# boundary surface area rather than a periphery). If we multiply the 2*gamma
# quantity times the thickness t = 3*4.08 = 12.24 Å, we obtain 2*gamma*t =
# 3.914352 eV/Å, which agrees with the mesh calculation.
# 3.914352 eV/Å, which agrees with the mesh calculation.
echo both
log eshelby_static.log
log eshelby_static.log
units metal
dimension 3
boundary p p p
@ -29,12 +29,12 @@ lattice fcc 4.08 origin 0.25 0.25 0.25
region SYSTEM block -10 10 -10 10 0 3
region UPPER block -10 10 0 10 0 3
region LOWER block -10 0 -10 10 0 3
create_box 1 SYSTEM
create_box 1 SYSTEM
create_atoms 1 region SYSTEM
mass 1 63.55
group internal region SYSTEM
group upper region UPPER
group lower region LOWER
group internal region SYSTEM
group upper region UPPER
group lower region LOWER
#pair_style eam
#pair_coeff * * ../../../../potentials/Cu_u3.eam
#pair_style lj/smooth 8.0 10.0
@ -42,27 +42,27 @@ group lower region LOWER
pair_style lj/smooth/linear 5.456108274435118
pair_coeff * * 0.7242785984051078 2.598146797350056
# define region
region abovecrack block -2.01 2.01 0 8 INF INF units lattice
region abovecrack block -2.01 2.01 0 8 INF INF units lattice
group abovecrack region abovecrack
region belowcrack block -2.01 2.01 -8 0 INF INF units lattice
region belowcrack block -2.01 2.01 -8 0 INF INF units lattice
group belowcrack region belowcrack
neighbor 1.0 bin
neigh_modify delay 1000000
neigh_modify exclude group abovecrack belowcrack
neigh_modify exclude group abovecrack belowcrack
thermo 10
compute_modify thermo_temp extra 0
compute_modify thermo_temp extra/dof 0
region FORCECHECK block -6 0 6 10 0 3
group FORCECHECK region FORCECHECK
compute fxsum FORCECHECK reduce sum fx
compute fysum FORCECHECK reduce sum fy
compute fzsum FORCECHECK reduce sum fz
thermo_style custom step ke pe press c_fxsum c_fysum c_fzsum
thermo_modify format 1 %4i format 2 %3.1g format 3 %20.16g
thermo_modify format 1 %4i format 2 %3.1g format 3 %20.16g
timestep 0.0
# (1) minimize the system to get a relaxed configuration
min_modify line quadratic
variable tol equal 1.e-11
minimize 0.0 ${tol} 100000 1000000
minimize 0.0 ${tol} 100000 1000000
write_restart eshelby.restart
# (2) pull crack apart
reset_timestep 0
@ -70,13 +70,13 @@ fix PK1 internal atc field
# NOTE tune this relative to lattice size
fix_modify PK1 mesh create 10 10 1 SYSTEM p f p
fix_modify PK1 fields none
fix_modify PK1 fields add mass_density internal_energy temperature stress displacement
fix_modify PK1 fields add eshelby_stress transformed_stress
fix_modify PK1 gradients add displacement
fix_modify PK1 fields add mass_density internal_energy temperature stress displacement
fix_modify PK1 fields add eshelby_stress transformed_stress
fix_modify PK1 gradients add displacement
fix_modify PK1 set reference_potential_energy
fix_modify PK1 output eshelby_staticFE 1 text binary tensor_components
#fix_modify PK1 on_the_fly
#-- make concentric loops around one/both of the crack tips
#fix_modify PK1 on_the_fly
#-- make concentric loops around one/both of the crack tips
#-- & another around undefected material (loop0)
fix_modify PK1 mesh create_faceset loop0 box -6 0 6 10 -INF INF outward
fix_modify PK1 output boundary_integral eshelby_stress faceset loop0

View File

@ -5,12 +5,12 @@ atom_style atomic
echo both
dimension 3
boundary p p p
lattice fcc 5.376 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
region box block 0 5 0 5 0 5
read_data nvt.init
lattice fcc 5.376 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
region box block 0 5 0 5 0 5
read_data nvt.init
pair_style lj/cut 13.0
pair_coeff 1 1 0.2381 3.405
group internal region box
pair_coeff 1 1 0.2381 3.405
group internal region box
fix ATC internal atc field
fix_modify ATC mesh create 1 1 1 box p p p
@ -27,12 +27,12 @@ fix_modify PP atom_element_map eulerian 100
fix_modify PP fields add mass_density energy stress temperature kinetic_temperature
fix_modify PP fields add velocity thermal_energy kinetic_energy
fix_modify PP output nvtPP 100 text
thermo 100
compute_modify thermo_temp extra 0
thermo_style custom step temp etotal ke pe press vol
timestep 4.0
reset_timestep 0
thermo 100
compute_modify thermo_temp extra/dof 0
thermo_style custom step temp etotal ke pe press vol
timestep 4.0
reset_timestep 0
# NOTE fixes performing time integration (init_integrate/final_integrate) should be defined after atc
fix NVT all nvt temp 30 30 10.0 drag 0.2 tchain 1
log nvt.log
run 1000
fix NVT all nvt temp 30 30 10.0 drag 0.2 tchain 1
log nvt.log
run 1000

View File

@ -1,38 +1,38 @@
#AtC Thermal Coupling
# This benchmark tests heat conducting into and out of the MD region. The
# temperature is initially 20 everywhere and the left boundary BC is fixed at 40.
# The result should show heat diffusing through the FEM to the MD and back out
# to the FEM on the right. Insufficient time is captured to reach the linear
# The result should show heat diffusing through the FEM to the MD and back out
# to the FEM on the right. Insufficient time is captured to reach the linear
# steady state, but heat crossing both boundaries should be observed.
units real
atom_style atomic
units real
atom_style atomic
echo both
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
read_data temp.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region simRegion block -12 12 -3 3 -3 3
region simRegion block -12 12 -3 3 -3 3
# create atoms, NOTE commented out for restart
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
# specify interal/ghost atoms
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
# velocities have Vcm = 0, NOTE next four lines commented out for restart
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
# ID group atc PhysicsType ParameterFile
fix AtC internal atc thermal Ar_thermal.mat
fix_modify AtC boundary ghost
@ -40,26 +40,26 @@ fix_modify AtC boundary ghost
fix_modify AtC mesh create 12 1 1 simRegion f p p
fix_modify AtC internal_quadrature off
# fix a temperature
fix_modify AtC fix temperature all 20.
fix_modify AtC fix temperature all 20.
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
# equilibrate MD field
#dump D1 all atom 100 dump.bar1d
#fix_modify AtC output bar1dfe 100
timestep 5
timestep 5
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
thermo 100
run 400
thermo 100
run 400
# change thermostat to hoover
fix_modify AtC unfix temperature all
fix_modify AtC control thermal hoover
#fix_modify AtC localized_lambda on
fix_modify AtC filter type exponential
fix_modify AtC filter scale 10000.0
fix_modify AtC filter on
fix_modify AtC unfix temperature all
fix_modify AtC control thermal hoover
#fix_modify AtC localized_lambda on
fix_modify AtC filter type exponential
fix_modify AtC filter scale 10000.0
fix_modify AtC filter on
# add nodesets
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
fix_modify AtC mesh create_nodeset lbc -12.1 -11.9 -INF INF -INF INF
@ -67,13 +67,13 @@ fix_modify AtC mesh create_nodeset rbc 11.9 12.1 -INF INF -INF INF
fix_modify AtC fix temperature lbc 40.
fix_modify AtC fix temperature rbc 20.
# initialize filter
thermo 1
run 100
thermo 1
run 100
# set up output, should be before a "run"
fix_modify AtC output bar1dFE 100 text binary
# output command
#dump D1 all atom 1000 dump.bar1d
# run with FE
reset_timestep 0
thermo 100
run 10000
thermo 100
run 10000

View File

@ -4,33 +4,33 @@
# 40.# The result should show heat diffusing through the FEM to the MD and back
# out # to the FEM on the right. Insufficient time is captured to reach the
# linear # steady state, but heat crossing both boundaries should be observed.
units real
atom_style atomic
units real
atom_style atomic
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
read_data all_atoms.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region simRegion block -7 7 -3 3 -3 3
region simRegion block -7 7 -3 3 -3 3
# create atoms, NOTE commented out for restart
#boundary f p p
#create_box 1 simRegion
#create_atoms 1 region simRegion
#mass 1 39.95
#boundary f p p
#create_box 1 simRegion
#create_atoms 1 region simRegion
#mass 1 39.95
# specify interal/ghost atoms
region mdInternal block -6 6 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
region mdInternal block -6 6 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
# velocities have Vcm = 0, NOTE next four lines commented out for restart
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
#write_restart all_atoms.init
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
#write_restart all_atoms.init
neighbor 5. bin
neigh_modify every 10 delay 0 check no
# ID group atc PhysicsType ParameterFile
fix AtC internal atc thermal Ar_thermal.mat
# ID part keywords nx ny nz region
@ -39,34 +39,34 @@ fix_modify AtC mesh create_faceset obndy box -6.0 6.0 -INF INF -INF INF out
fix_modify AtC mesh create_faceset lbndy box -6.0 16.0 -INF INF -INF INF outward
fix_modify AtC mesh create_faceset rbndy box -16.0 6.0 -INF INF -INF INF outward
# fix a temperature
fix_modify AtC fix temperature all 20.
fix_modify AtC fix temperature all 20.
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
fix_modify AtC internal_quadrature off
# equilibrate MD field
#dump D1 all atom 100 dump.bar1d
#fix_modify AtC output bar1dfe 100
#fix_modify AtC time_integration frac_step
timestep 5
#fix_modify AtC time_integration frac_step
timestep 5
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
thermo 100
run 400
thermo 100
run 400
# change thermostat to hoover
fix_modify AtC control thermal hoover
#fix_modify AtC thermal control flux interpolate
#fix_modify AtC thermal control flux faceset obndy
fix_modify AtC filter scale 1000.0
#fix_modify AtC filter equilibrate
#fix_modify AtC filter on
#fix_modify AtC localized_lambda on
fix_modify AtC control thermal hoover
#fix_modify AtC thermal control flux interpolate
#fix_modify AtC thermal control flux faceset obndy
fix_modify AtC filter scale 1000.0
#fix_modify AtC filter equilibrate
#fix_modify AtC filter on
#fix_modify AtC localized_lambda on
fix_modify AtC boundary_faceset is obndy
# initialize filter
#fix_modify AtC output bar1d_all_atoms_preqFE 1
reset_timestep 0
thermo 1
thermo 1
#run 5000
# add nodesets and ramp temperature
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
@ -76,26 +76,26 @@ thermo 1
#fix_modify AtC fix temperature rbc 20.
# initialize filter
#fix_modify AtC output bar1d_all_atoms_eqFE 200
reset_timestep 0
fix_modify AtC reset_time
reset_timestep 0
fix_modify AtC reset_time
#run 5000
# release and let run
#fix_modify AtC filter on
fix_modify AtC unfix temperature all
#fix_modify AtC filter on
fix_modify AtC unfix temperature all
#fix_modify AtC localized_lambda on
fix_modify AtC control lumped_lambda_solve on
fix_modify AtC control lumped_lambda_solve on
# dirichlet boundaries
#fix_modify AtC fix temperature lbc 40.
#fix_modify AtC fix temperature rbc 20.
# neumann boundaries
fix_modify AtC control thermal flux no_boundary
fix_modify AtC fix_flux temperature lbndy 0.0000000001
fix_modify AtC fix_flux temperature lbndy 0.0000000001
fix_modify AtC fix_flux temperature rbndy -0.0000000001
# set up output, should be before a "run"
fix_modify AtC output bar1d_all_atomsFE 200 text binary
# output command
#dump D1 all atom 1000 dump.bar1d
# run with FE
reset_timestep 0
thermo 200
run 5000
reset_timestep 0
thermo 200
run 5000

View File

@ -6,62 +6,62 @@
# linear # steady state, but heat crossing the boundaries should be observed,
# as should the interaction of the two themostat types.
#echo both
units real
atom_style atomic
units real
atom_style atomic
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
boundary f p p
pair_style lj/cut 13.5
read_data temp.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region simRegion block -16 4 -3 3 -3 3
region simRegion block -16 4 -3 3 -3 3
# create atoms, NOTE commented out for restart
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
# specify interal/ghost atoms
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
region mdGhost block -5 -4 -3 3 -3 3
group ghost region mdGhost
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
region mdGhost block -5 -4 -3 3 -3 3
group ghost region mdGhost
# velocities have Vcm = 0, NOTE next four lines commented out for restart
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#write_restart tinit.dat
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#write_restart tinit.dat
# ID group atc PhysicsType ParameterFile
fix AtC internal atc thermal Ar_thermal.mat
fix_modify AtC boundary ghost
fix_modify AtC time_integration fractional_step
fix_modify AtC internal_quadrature off
fix_modify AtC internal_quadrature off
# ID part keywords nx ny nz region
fix_modify AtC mesh create 10 1 1 simRegion f p p
fix_modify AtC mesh create_faceset ibndy box -4.0 40.0 -INF INF -INF INF in
fix_modify AtC mesh create_faceset obndy box -4.0 40.0 -INF INF -INF INF outward
# fix a temperature
fix_modify AtC fix temperature all 20.
fix_modify AtC fix temperature all 20.
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
# equilibrate MD field
#fix_modify AtC output bar1dfe_combined_init 100
timestep 5
thermo 10
run 400
timestep 5
thermo 10
run 400
# change thermostat to flux
fix_modify AtC unfix temperature all
fix_modify AtC control thermal flux faceset obndy
#fix_modify AtC control thermal hoover
fix_modify AtC control localized_lambda on
fix_modify AtC unfix temperature all
fix_modify AtC control thermal flux faceset obndy
#fix_modify AtC control thermal hoover
fix_modify AtC control localized_lambda on
#fix_modify AtC control tolerance 1.e-16 # this tolerance seems necessary to prevent noticeable drift
# add nodesets
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
@ -74,8 +74,8 @@ fix_modify AtC fix temperature rbc 20.
fix_modify AtC filter type exponential
fix_modify AtC filter scale 1000.0
#fix_modify AtC filter on
thermo 100
#run 100
thermo 100
#run 100
# set up output, should be before a "run"
fix_modify AtC output bar1d_combinedFE 100 text
@ -83,5 +83,5 @@ fix_modify AtC output bar1d_combinedFE 100 text
#dump D1 all atom 100 dump.bar1d_combined
# run with FE
reset_timestep 0
thermo 100
run 5000
thermo 100
run 5000

View File

@ -5,34 +5,34 @@
# out # to the FEM on the right. Insufficient time is captured to reach the
# linear # steady state, but heat crossing both boundaries should be observed.
#echo both
units real
atom_style atomic
units real
atom_style atomic
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
boundary f p p
pair_style lj/cut 13.5
read_data temp.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region simRegion block -12 12 -3 3 -3 3
region simRegion block -12 12 -3 3 -3 3
# create atoms, NOTE commented out for restart
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
# specify interal/ghost atoms
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
# velocities have Vcm = 0, NOTE next four lines commented out for restart
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#write_restart tinit.dat
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#write_restart tinit.dat
# ID group atc PhysicsType ParameterFile
fix AtC internal atc thermal Ar_thermal.mat
fix_modify AtC boundary ghost
@ -42,25 +42,25 @@ fix_modify AtC mesh create 12 1 1 simRegion f p p
fix_modify AtC mesh create_faceset ibndy box -4.0 4.0 -INF INF -INF INF in
fix_modify AtC mesh create_faceset obndy box -4.0 4.0 -INF INF -INF INF outward
# fix a temperature
fix_modify AtC fix temperature all 20.
#fix_modify AtC fix temperature all linear 0 0 0 -0.154 0 0 30.
fix_modify AtC fix temperature all 20.
#fix_modify AtC fix temperature all linear 0 0 0 -0.154 0 0 30.
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
# equilibrate MD field
#fix_modify AtC output bar1dfe_flux_init 100
timestep 5
thermo 10
run 400
timestep 5
thermo 10
run 400
# change thermostat to hoover
fix_modify AtC unfix temperature all
fix_modify AtC control thermal flux faceset obndy
fix_modify AtC unfix temperature all
fix_modify AtC control thermal flux faceset obndy
fix_modify AtC filter type exponential
fix_modify AtC filter scale 10000.0
#fix_modify AtC filter on
fix_modify AtC filter scale 10000.0
#fix_modify AtC filter on
# add nodesets
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
fix_modify AtC mesh create_nodeset lbc -12.1 -11.9 -INF INF -INF INF
@ -68,12 +68,12 @@ fix_modify AtC mesh create_nodeset rbc 11.9 12.1 -INF INF -INF INF
fix_modify AtC fix temperature lbc 40.
fix_modify AtC fix temperature rbc 20.
# initialize filter
run 100
run 100
# set up output, should be before a "run"
fix_modify AtC output bar1d_fluxFE 100 text binary
# output command
#dump D1 all atom 1000 dump.bar1d
# run with FE
reset_timestep 0
thermo 100
run 10000
thermo 100
run 10000

View File

@ -4,33 +4,33 @@
# 40.# The result should show heat diffusing through the FEM to the MD and back
# out # to the FEM on the right. Insufficient time is captured to reach the
# linear # steady state, but heat crossing both boundaries should be observed.
units real
atom_style atomic
units real
atom_style atomic
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
read_data temp.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region simRegion block -12 12 -3 3 -3 3
region simRegion block -12 12 -3 3 -3 3
# create atoms, NOTE commented out for restart
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
# specify interal/ghost atoms
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
# velocities have Vcm = 0, NOTE next four lines commented out for restart
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
# ID group atc PhysicsType ParameterFile
fix AtC internal atc thermal Ar_thermal.mat
# ID part keywords nx ny nz region
@ -40,7 +40,7 @@ fix_modify AtC internal_atom_integrate off
fix iNVE internal nve
fix_modify AtC internal_quadrature off
# fix a temperature
fix_modify AtC fix temperature all 20.
fix_modify AtC fix temperature all 20.
# add nodesets
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
fix_modify AtC mesh create_nodeset lbc -4.1 -3.9 -INF INF -INF INF
@ -48,25 +48,25 @@ fix_modify AtC mesh create_nodeset rbc 3.9 4.1 -INF INF -INF INF
fix_modify AtC fix temperature lbc 40.
fix_modify AtC fix temperature rbc 20.
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
# equilibrate MD field
#dump D1 all atom 100 dump.bar1d
#fix_modify AtC output bar1dfe 100 text
timestep 5
timestep 5
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
thermo 100
run 400
thermo 100
run 400
# change thermostat to hoover
fix_modify AtC unfix temperature all
fix_modify AtC fix temperature lbc 40.
fix_modify AtC fix temperature rbc 20.
fix_modify AtC control thermal flux no_boundary
fix_modify AtC control tolerance 1.e-9 # this is to remove looser benchmark tolerance because the lambda solution was reacting to numerical noise, but it results in a slight drift in the constraints
#fix_modify AtC thermal control none
fix_modify AtC control localized_lambda on
fix_modify AtC control tolerance 1.e-9 # this is to remove looser benchmark tolerance because the lambda solution was reacting to numerical noise, but it results in a slight drift in the constraints
#fix_modify AtC thermal control none
fix_modify AtC control localized_lambda on
fix_modify AtC filter type exponential
fix_modify AtC filter scale 1000.0
fix_modify AtC filter on
@ -76,5 +76,5 @@ fix_modify AtC output bar1d_frac_stepFE 200 text
#dump D1 all atom 1000 dump.bar1d
# run with FE
reset_timestep 0
thermo 100
run 10000
thermo 100
run 10000

View File

@ -4,35 +4,35 @@
# 40.# The result should show heat diffusing through the FEM to the MD and back
# out # to the FEM on the right. Insufficient time is captured to reach the
# linear # steady state, but heat crossing both boundaries should be observed.
units real
atom_style atomic
units real
atom_style atomic
log bar1d_hoover.log
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
read_data all_atoms.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region simRegion block -7 7 -3 3 -3 3
region simRegion block -7 7 -3 3 -3 3
# create atoms, NOTE commented out for restart
#boundary f p p
#create_box 1 simRegion
#create_atoms 1 region simRegion
#mass 1 39.95
#atom_modify sort 0 0
#boundary f p p
#create_box 1 simRegion
#create_atoms 1 region simRegion
#mass 1 39.95
#atom_modify sort 0 0
# specify interal/ghost atoms
region mdInternal block -6 6 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
region mdInternal block -6 6 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
# velocities have Vcm = 0, NOTE next four lines commented out for restart
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
#write_restart all_atoms.init
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
##pair_coeff 1 1 0.010323166 3.405 13.5
#pair_coeff 1 1 .238 3.405 13.5
#write_restart all_atoms.init
neighbor 5. bin
neigh_modify every 10 delay 0 check no
# ID group atc PhysicsType ParameterFile
fix AtC internal atc thermal Ar_thermal.mat
# ID part keywords nx ny nz region
@ -40,9 +40,9 @@ fix AtC internal atc thermal Ar_thermal.mat
# switch depending on application of BCs or coupling
#variable nEls equal 6 # for boundary conditions
variable nEls equal 12 # for coupling
fix_modify AtC boundary ghost # for coupling
fix_modify AtC boundary ghost # for coupling
region atcRegion block -${nEls} ${nEls} -3 3 -3 3
region atcRegion block -${nEls} ${nEls} -3 3 -3 3
fix_modify AtC mesh create ${nEls} 1 1 atcRegion f p p
fix_modify AtC mesh create_faceset obndy box -6.0 6.0 -INF INF -INF INF outward
fix_modify AtC mesh create_faceset lbndy box -6.0 16.0 -INF INF -INF INF outward
@ -50,39 +50,39 @@ fix_modify AtC mesh create_faceset rbndy box -16.0 6.0 -INF INF -INF INF ou
variable tol equal 0.1
variable uBnd equal ${nEls}+${tol}
variable lBnd equal ${nEls}-${tol}
fix_modify AtC mesh create_nodeset lbc -${uBnd} -${lBnd} -INF INF -INF INF
fix_modify AtC mesh create_nodeset lbc -${uBnd} -${lBnd} -INF INF -INF INF
fix_modify AtC mesh create_nodeset rbc ${lBnd} ${uBnd} -INF INF -INF INF
# fix a temperature
fix_modify AtC fix temperature all 20.
#fix_modify AtC initial temperature all 20.
#fix_modify AtC consistent_fe_initialization on
fix_modify AtC fix temperature all 20.
#fix_modify AtC initial temperature all 20.
#fix_modify AtC consistent_fe_initialization on
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
fix_modify AtC internal_quadrature off
# equilibrate MD field
#dump D1 all atom 100 dump.bar1d_hoover_init
#fix_modify AtC output bar1d_hoover_init 100 text binary
fix_modify AtC time_integration fractional_step
fix_modify AtC temperature_definition total
timestep 5
fix_modify AtC time_integration fractional_step
fix_modify AtC temperature_definition total
timestep 5
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
thermo 100
run 400
thermo 100
run 400
fix_modify AtC unfix temperature all
fix_modify AtC unfix temperature all
fix_modify AtC filter type exponential
fix_modify AtC filter type exponential
fix_modify AtC filter scale 1000.0
fix_modify AtC filter on
# boundary conditions------------------------
# dirichlet with ramp (unstable when filtered due to instantaneous change in rate)
#fix_modify AtC fix temperature lbc temporal_ramp 20. 40. 25000.
#fix_modify AtC fix temperature rbc 20.
#fix_modify AtC fix temperature lbc temporal_ramp 20. 40. 25000.
#fix_modify AtC fix temperature rbc 20.
#fix_modify AtC localized_lambda on
#fix_modify AtC control thermal flux no_boundary
#fix_modify AtC output test_rampFE 200 text binary
@ -94,7 +94,7 @@ fix_modify AtC filter on
#fix_modify AtC fix temperature lbc 40.
# dirichlet with rescaling
#fix_modify AtC filter off
#fix_modify AtC filter off
#fix_modify AtC fix temperature all 20.
#fix_modify AtC fix temperature lbc 40.
#fix_modify AtC output test_rampFE 200 text binary
@ -126,13 +126,13 @@ fix_modify AtC fix temperature rbc 20.
#fix_modify AtC lumped_lambda_solve on
# hoover
fix_modify AtC control thermal hoover
fix_modify AtC control thermal hoover
fix_modify AtC control localized_lambda on
fix_modify AtC control tolerance 1.e-14
fix_modify AtC control tolerance 1.e-14
fix_modify AtC output bar1d_hooverFE 200 text binary
#dump D1 all atom 200 dump.bar1d_hoover
reset_timestep 0
fix_modify AtC reset_time
thermo 100
thermo 100
run 5000

View File

@ -4,65 +4,65 @@
# 40.# The result should show heat diffusing through the FEM to the MD and back
# out # to the FEM on the right. Insufficient time is captured to reach the
# linear # steady state, but heat crossing both boundaries should be observed.
units real
atom_style atomic
units real
atom_style atomic
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
read_data temp.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region simRegion block -12 12 -3 3 -3 3
region simRegion block -12 12 -3 3 -3 3
# create atoms, NOTE commented out for restart
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
# specify interal/ghost atoms
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
# velocities have Vcm = 0, NOTE next three lines commented out for restart
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
# ID group atc PhysicsType ParameterFile
fix AtC internal atc thermal Ar_thermal.mat
fix_modify AtC boundary ghost
# numerical parameters
fix_modify AtC time_integration fractional_step
fix_modify AtC time_integration fractional_step
# ID part keywords nx ny nz region
fix_modify AtC mesh create 12 1 1 simRegion f p p
fix_modify AtC mesh create_faceset ibndy box -4.0 4.0 -INF INF -INF INF in
fix_modify AtC mesh create_faceset obndy box -4.0 4.0 -INF INF -INF INF outward
# fix a temperature
fix_modify AtC fix temperature all 20.
#fix_modify AtC fix temperature all linear 0 0 0 -0.154 0 0 30.
#fix_modify AtC fix temperature all quadratic 0 0 0 -0.154 0 0 30. 0.01 0.0 0.0 0 0 0
fix_modify AtC fix temperature all 20.
#fix_modify AtC fix temperature all linear 0 0 0 -0.154 0 0 30.
#fix_modify AtC fix temperature all quadratic 0 0 0 -0.154 0 0 30. 0.01 0.0 0.0 0 0 0
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
# equilibrate MD field
#fix_modify AtC output bar1dfe_flux_init 100
timestep 5
timestep 5
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
thermo 100
run 400
thermo 100
run 400
# change thermostat to hoover
fix_modify AtC unfix temperature all
fix_modify AtC control thermal none
fix_modify AtC control thermal flux interpolate
fix_modify AtC control localized_lambda on
#fix_modify AtC thermal control flux faceset obndy
fix_modify AtC unfix temperature all
fix_modify AtC control thermal none
fix_modify AtC control thermal flux interpolate
fix_modify AtC control localized_lambda on
#fix_modify AtC thermal control flux faceset obndy
fix_modify AtC filter type exponential
fix_modify AtC filter scale 10000.0
fix_modify AtC filter on
fix_modify AtC filter scale 10000.0
fix_modify AtC filter on
# add nodesets
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
fix_modify AtC mesh create_nodeset lbc -12.1 -11.9 -INF INF -INF INF
@ -70,11 +70,11 @@ fix_modify AtC mesh create_nodeset rbc 11.9 12.1 -INF INF -INF INF
fix_modify AtC fix temperature lbc 40.
fix_modify AtC fix temperature rbc 20.
# initialize filter
run 100
run 100
# set up output, should be before a "run"
fix_modify AtC output bar1d_interpolateFE 100 text
# run with FE
reset_timestep 0
thermo 100
#run 10
run 10000
thermo 100
#run 10
run 10000

View File

@ -5,32 +5,32 @@
# out # to the FEM on the right. Insufficient time is captured to reach the
# linear # steady state, but heat crossing both boundaries should be observed.
echo both
units real
atom_style atomic
units real
atom_style atomic
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
# NOTE following 3 lines added for restart
boundary f p p
pair_style lj/cut 13.5
read_data temp.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region simRegion block -12 12 -3 3 -3 3
region simRegion block -12 12 -3 3 -3 3
# create atoms, NOTE commented out for restart
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
#region mdRegion block -5 5 -3 3 -3 3
#boundary f p p
#create_box 1 mdRegion
#create_atoms 1 region mdRegion
#mass 1 39.95
# specify interal/ghost atoms
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
region mdInternal block -4 4 -3 3 -3 3
group internal region mdInternal
group ghost subtract all internal
# velocities have Vcm = 0, NOTE next three lines commented out for restart
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
#velocity internal create 40. 87287 mom yes loop geom
#pair_style lj/cut 13.5
#pair_coeff 1 1 .238 3.405 13.5
neighbor 5. bin
neigh_modify every 10 delay 0 check no
# ID group atc PhysicsType ParameterFile
fix AtC internal atc thermal Ar_thermal.mat
fix_modify AtC boundary ghost
@ -42,30 +42,30 @@ fix_modify AtC mesh create 12 1 1 simRegion f p p
fix_modify AtC mesh create_faceset ibndy box -4.0 4.0 -INF INF -INF INF in
fix_modify AtC mesh create_faceset obndy box -4.0 4.0 -INF INF -INF INF outward
# fix a temperature
fix_modify AtC fix temperature all 20.
#fix_modify AtC fix temperature all linear 0 0 0 -0.154 0 0 30.
#fix_modify AtC fix temperature all quadratic 0 0 0 -0.154 0 0 30. 0.01 0.0 0.0 0 0 0
fix_modify AtC fix temperature all 20.
#fix_modify AtC fix temperature all linear 0 0 0 -0.154 0 0 30.
#fix_modify AtC fix temperature all quadratic 0 0 0 -0.154 0 0 30. 0.01 0.0 0.0 0 0 0
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
# output
thermo 100
thermo 100
variable xdof equal 3*count(ghost)
compute_modify thermo_temp extra ${xdof}
compute_modify thermo_temp extra/dof ${xdof}
thermo_style custom step cpu etotal pe f_AtC[1] temp f_AtC[2]
thermo_modify format 1 %6i format 2 %7.2g
# equilibrate MD field
#fix_modify AtC output bar1d_lumped_initFE 1 text
timestep 5
run 400
#run 100
timestep 5
run 400
#run 100
# change thermostat to hoover
fix_modify AtC unfix temperature all
#fix_modify AtC thermal control flux interpolate
fix_modify AtC control thermal flux faceset obndy
fix_modify AtC control localized_lambda on
fix_modify AtC unfix temperature all
#fix_modify AtC thermal control flux interpolate
fix_modify AtC control thermal flux faceset obndy
fix_modify AtC control localized_lambda on
fix_modify AtC filter type exponential
fix_modify AtC filter scale 10000.0
fix_modify AtC filter on
fix_modify AtC filter scale 10000.0
fix_modify AtC filter on
# add nodesets
# ID mesh create_nodeset tag xmin xmax ymin ymax zmin zmax
fix_modify AtC mesh create_nodeset lbc -12.1 -11.9 -INF INF -INF INF
@ -73,13 +73,13 @@ fix_modify AtC mesh create_nodeset rbc 11.9 12.1 -INF INF -INF INF
fix_modify AtC fix temperature lbc 80.
fix_modify AtC fix temperature rbc 80.
# initialize
thermo 10
run 100
thermo 10
run 100
# relax
fix_modify AtC unfix temperature all
# output
fix_modify AtC output bar1d_lumpedFE 100 text
# run with FE
reset_timestep 0
thermo 100
run 10000
thermo 100
run 10000

View File

@ -1,41 +1,41 @@
# needs description
#AtC Thermal Coupling
echo both
units real
atom_style atomic
boundary f p p
units real
atom_style atomic
boundary f p p
# create domain
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
pair_style lj/cut 13.5
read_data uniform_heating.init
#lattice type reduced density rho* = 4*(sigma/a)^3, where N=4 for fcc, s = 3.405 A (Wagner) and a = 5.25 A (Ashcroft & Mermin, p. 70)
pair_style lj/cut 13.5
read_data uniform_heating.init
lattice fcc 5.405 origin 0.25 0.25 0.25
region mdRegion block -8 8 -3 3 -3 3
region mdInternal block -6 6 -3 3 -3 3
group internal region mdInternal
neighbor 5. bin
neigh_modify every 10 delay 0 check no
region mdRegion block -8 8 -3 3 -3 3
region mdInternal block -6 6 -3 3 -3 3
group internal region mdInternal
neighbor 5. bin
neigh_modify every 10 delay 0 check no
# ID group atc PhysicsType ParameterFile
fix AtC internal atc two_temperature Ar_ttm.mat
# ID part keywords nx ny nz region
fix_modify AtC mesh create 6 1 1 mdInternal f p p
fix_modify AtC internal_quadrature on
# fix a temperature
fix_modify AtC initial temperature all 20.
fix_modify AtC fix temperature all 20.
fix_modify AtC initial electron_temperature all 20.
fix_modify AtC fix electron_temperature all 20.
fix_modify AtC initial temperature all 20.
fix_modify AtC fix temperature all 20.
fix_modify AtC initial electron_temperature all 20.
fix_modify AtC fix electron_temperature all 20.
# turn on thermostat
fix_modify AtC control thermal rescale 10
fix_modify AtC control thermal rescale 10
# output
thermo_style custom step cpu pe etotal temp f_AtC[1] f_AtC[2]
thermo 100
thermo 100
log uniform_heating.log
# make thermo output the correct temperature by removing ghost dof
# variable xdof equal 3*count(ghost)
# compute_modify thermo_temp extra ${xdof}
# compute_modify thermo_temp extra/dof ${xdof}
# equilibrate MD field
timestep 5
run 400
timestep 5
run 400
# change thermostat
fix_modify AtC output uniform_heatingFE 100 text
fix_modify AtC unfix temperature all
@ -50,4 +50,4 @@ fix_modify AtC fix temperature rbc 20.
fix_modify AtC source temperature all 6.3363e-10
# run with FE
reset_timestep 0
run 1000
run 1000

View File

@ -25,7 +25,7 @@ pair_modify shift yes
######### VERLET INTEGRATION WITH LANGEVIN THERMOSTAT ###########
fix fnve all nve
compute dtemp all temp
compute_modify dtemp dynamic yes
compute_modify dtemp dynamic/dof yes
fix fT all langevin 1.0 1.0 1.0 123
fix_modify fT temp dtemp

View File

@ -6,12 +6,12 @@ neighbor 10.0 bin
read_data data.chreg-acid-real
#real units
variable sigma equal 7.2 # particle diameter 0.72 nm
variable temperature equal 298 # temperature 298 K
variable sigma equal 7.2 # particle diameter 0.72 nm
variable temperature equal 298 # temperature 298 K
variable kb index 0.0019872067 # kB in Kcal/mol/K
variable epsilon equal ${kb}*${temperature}
variable tunit equal 2000 # time unit is 2000 fs
variable timestep equal 0.005*${tunit}
variable epsilon equal ${kb}*${temperature}
variable tunit equal 2000 # time unit is 2000 fs
variable timestep equal 0.005*${tunit}
variable cut_long equal 12.5*${sigma}
variable nevery equal 100
@ -27,13 +27,13 @@ velocity all create ${temperature} 8008 loop geom
pair_style lj/cut/coul/long ${cut_lj} ${cut_long}
pair_coeff * * ${epsilon} ${sigma}
kspace_style pppm 1.0e-3
dielectric 78
dielectric 78
pair_modify shift yes
######### VERLET INTEGRATION WITH LANGEVIN THERMOSTAT ###########
fix fnve all nve
compute dtemp all temp
compute_modify dtemp dynamic yes
compute_modify dtemp dynamic/dof yes
fix fT all langevin $(v_temperature) $(v_temperature) $(v_tunit) 123
fix_modify fT temp dtemp

View File

@ -17,7 +17,7 @@ dielectric 1.0
######### VERLET INTEGRATION WITH LANGEVIN THERMOSTAT ###########
fix fnve all nve
compute dtemp all temp
compute_modify dtemp dynamic yes
compute_modify dtemp dynamic/dof yes
fix fT all langevin 1.0 1.0 1.0 123
fix_modify fT temp dtemp

View File

@ -153,7 +153,7 @@ def mdi_engine(other_options):
# engine --> driver
mdi.MDI_Register_Command("@DEFAULT","<PE")
mdi.MDI_Register_Command("@DEFAULT","<ENERGY")
mdi.MDI_Register_Command("@DEFAULT","<FORCES")
mdi.MDI_Register_Command("@DEFAULT",">LATTICE_FORCES")
mdi.MDI_Register_Command("@DEFAULT","<STRESS")
@ -232,10 +232,10 @@ def execute_command(command,mdicomm,object_ptr):
# MDI commands which retreive quantum results
# each may also trigger the quantum calculation
elif command == "<PE":
elif command == "<ENERGY":
evaluate()
ierr = mdi.MDI_Send(qm_pe,1,mdi.MDI_DOUBLE,mdicomm)
if ierr: error("MDI: <PE data")
if ierr: error("MDI: <ENERGY data")
elif command == "<FORCES":
evaluate()

View File

@ -159,7 +159,7 @@ def mdi_engine(other_options):
# engine --> driver
mdi.MDI_Register_Command("@DEFAULT","<PE")
mdi.MDI_Register_Command("@DEFAULT","<ENERGY")
mdi.MDI_Register_Command("@DEFAULT","<FORCES")
mdi.MDI_Register_Command("@DEFAULT",">LATTICE_FORCES")
#mdi.MDI_Register_Command("@DEFAULT","<STRESS")
@ -236,10 +236,10 @@ def execute_command(command,mdicomm,object_ptr):
# MDI commands which retreive quantum results
# each may also trigger the quantum calculation
elif command == "<PE":
elif command == "<ENERGY":
evaluate()
ierr = mdi.MDI_Send(qm_pe,1,mdi.MDI_DOUBLE,mdicomm)
if ierr: error("MDI: <PE data")
if ierr: error("MDI: <ENERGY data")
elif command == "<FORCES":
evaluate()

View File

@ -189,7 +189,7 @@ def mdi_engine(other_options):
# engine --> driver
mdi.MDI_Register_Command("@DEFAULT","<PE")
mdi.MDI_Register_Command("@DEFAULT","<ENERGY")
mdi.MDI_Register_Command("@DEFAULT","<FORCES")
mdi.MDI_Register_Command("@DEFAULT","<LATTICE_FORCES")
#mdi.MDI_Register_Command("@DEFAULT","<STRESS")
@ -262,10 +262,10 @@ def execute_command(command,mdicomm,object_ptr):
# MDI commands which retreive quantum results
# each may also trigger the quantum calculation
elif command == "<PE":
elif command == "<ENERGY":
evaluate()
ierr = mdi.MDI_Send(qm_pe,1,mdi.MDI_DOUBLE,mdicomm)
if ierr: error("MDI: <PE data")
if ierr: error("MDI: <ENERGY data")
elif command == "<FORCES":
evaluate()

View File

@ -2,34 +2,34 @@
variable b index 0
variable x index 50
variable y index 20
variable d index 20
variable v index 5
variable w index 2
units lj
variable x index 50
variable y index 20
variable d index 20
variable v index 5
variable w index 2
units lj
dimension 2
atom_style atomic
atom_style atomic
boundary f f p
lattice hex 0.85
region box block 0 $x 0 $y -0.5 0.5
create_box 1 box
region circle sphere $(v_d/2+1) $(v_d/2/sqrt(3.0)+1) 0.0 $(v_d/2)
create_atoms 1 region circle
mass 1 1.0
lattice hex 0.85
region box block 0 $x 0 $y -0.5 0.5
create_box 1 box
region circle sphere $(v_d/2+1) $(v_d/2/sqrt(3.0)+1) 0.0 $(v_d/2)
create_atoms 1 region circle
mass 1 1.0
velocity all create 0.5 87287 loop geom
velocity all create 0.5 87287 loop geom
velocity all set $v $w 0 sum yes
pair_style lj/cut 2.5
pair_coeff 1 1 10.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff 1 1 10.0 1.0 2.5
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
fix 1 all nve
fix 1 all nve
fix 2 all wall/lj93 xlo 0.0 1 1 2.5 xhi $x 1 1 2.5
fix 3 all wall/lj93 ylo 0.0 1 1 2.5 yhi $y 1 1 2.5
@ -49,6 +49,6 @@ fix 10 all balance 50 0.9 rcb
#dump_modify 3 pad 5 amap 0 10 sa 1 10 ${colors}
thermo_style custom step temp epair press f_10[3] f_10
thermo 100
thermo 100
run 10000
run 10000

View File

@ -2,29 +2,29 @@
variable b index 0
variable x index 50
variable y index 20
variable d index 20
variable v index 5
variable w index 2
units lj
variable x index 50
variable y index 20
variable d index 20
variable v index 5
variable w index 2
units lj
dimension 2
atom_style bond
atom_style bond
boundary f f p
lattice hex 0.85
region box block 0 $x 0 $y -0.5 0.5
create_box 1 box bond/types 1 extra/bond/per/atom 6
region circle sphere $(v_d/2+1) $(v_d/2/sqrt(3.0)+1) 0.0 $(v_d/2)
create_atoms 1 region circle
mass 1 1.0
lattice hex 0.85
region box block 0 $x 0 $y -0.5 0.5
create_box 1 box bond/types 1 extra/bond/per/atom 6
region circle sphere $(v_d/2+1) $(v_d/2/sqrt(3.0)+1) 0.0 $(v_d/2)
create_atoms 1 region circle
mass 1 1.0
velocity all create 0.5 87287 loop geom
velocity all create 0.5 87287 loop geom
velocity all set $v $w 0 sum yes
pair_style lj/cut 2.5
pair_coeff 1 1 10.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff 1 1 10.0 1.0 2.5
bond_style harmonic
bond_coeff 1 10.0 1.2
@ -34,10 +34,10 @@ bond_coeff 1 10.0 1.2
special_bonds lj/coul 0 1 1
create_bonds many all all 1 1.0 1.5
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
fix 1 all nve
fix 1 all nve
fix 2 all wall/lj93 xlo 0.0 1 1 2.5 xhi $x 1 1 2.5
fix 3 all wall/lj93 ylo 0.0 1 1 2.5 yhi $y 1 1 2.5
@ -58,6 +58,6 @@ fix 10 all balance 50 0.9 rcb
#dump_modify 3 pad 5 amap 0 10 sa 1 10 ${colors}
thermo_style custom step temp epair press f_10[3] f_10
thermo 100
thermo 100
run 10000
run 10000

View File

@ -2,42 +2,42 @@
variable b index 0
variable x index 50
variable y index 20
variable d index 20
variable x index 50
variable y index 20
variable d index 20
# careful not to slam into wall too hard
variable v index 0.3
variable w index 0.08
units lj
variable v index 0.3
variable w index 0.08
units lj
dimension 2
atom_style bond
atom_style bond
boundary f f p
lattice hex 0.85
region box block 0 $x 0 $y -0.5 0.5
create_box 1 box bond/types 1 extra/bond/per/atom 6
region circle sphere $(v_d/2+1) $(v_d/2/sqrt(3.0)+1) 0.0 $(v_d/2)
create_atoms 1 region circle
mass 1 1.0
lattice hex 0.85
region box block 0 $x 0 $y -0.5 0.5
create_box 1 box bond/types 1 extra/bond/per/atom 6
region circle sphere $(v_d/2+1) $(v_d/2/sqrt(3.0)+1) 0.0 $(v_d/2)
create_atoms 1 region circle
mass 1 1.0
velocity all create 0.5 87287 loop geom
velocity all create 0.5 87287 loop geom
velocity all set $v $w 0 sum yes
pair_style lj/cut 2.5
pair_coeff 1 1 10.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff 1 1 10.0 1.0 2.5
bond_style harmonic
bond_coeff 1 10.0 1.2
create_bonds many all all 1 1.0 1.5
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
fix 1 all nve
fix 1 all nve
fix 2 all wall/lj93 xlo 0.0 1 1 2.5 xhi $x 1 1 2.5
fix 3 all wall/lj93 ylo 0.0 1 1 2.5 yhi $y 1 1 2.5
@ -57,6 +57,6 @@ fix 10 all balance 50 0.9 rcb
#dump_modify 3 pad 5 amap 0 10 sa 1 10 ${colors}
thermo_style custom step temp epair press f_10[3] f_10
thermo 100
thermo 100
run 40000
run 40000

View File

@ -1,54 +1,54 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
units lj
atom_style atomic
processors * 1 1
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
fix p all property/atom d_WEIGHT
compute p all property/atom d_WEIGHT
fix 0 all balance 50 1.0 shift x 10 1.0 &
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
fix p all property/atom d_WEIGHT
compute p all property/atom d_WEIGHT
fix 0 all balance 50 1.0 shift x 10 1.0 &
weight time 1.0 weight store WEIGHT
variable maximb equal f_0[1]
variable iter equal f_0[2]
variable prev equal f_0[3]
variable final equal f_0
variable maximb equal f_0[1]
variable iter equal f_0[2]
variable prev equal f_0[3]
variable final equal f_0
#fix 3 all print 50 "${iter} ${prev} ${final} ${maximb}"
#fix 3 all print 50 "${iter} ${prev} ${final} ${maximb}"
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all custom 50 dump.lammpstrj id type x y z c_p
#dump id all atom 50 dump.melt
#dump id all custom 50 dump.lammpstrj id type x y z c_p
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
thermo 50
run 500
run 500
fix 0 all balance 50 1.0 shift x 5 1.0 &
thermo 50
run 500
run 500
fix 0 all balance 50 1.0 shift x 5 1.0 &
weight neigh 0.5 weight time 0.66 weight store WEIGHT
run 500
run 500

View File

@ -1,48 +1,48 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
units lj
atom_style atomic
processors * 1 1
variable factor index 1.0
variable factor index 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
group fast type 1
group slow type 2
balance 1.0 shift x 5 1.1 weight time 1.0 # out unweighted.txt
group fast type 1
group slow type 2
balance 1.0 shift x 5 1.1 weight time 1.0 # out unweighted.txt
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all atom 50 dump.melt
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
thermo 50
run 250
balance 1.0 shift x 5 1.1 weight time 1.0 # out unweighted.txt
run 250
balance 1.0 shift x 5 1.1 weight time 1.0 # out unweighted.txt
run 250
thermo 50
run 250
balance 1.0 shift x 5 1.1 weight time 1.0 # out unweighted.txt
run 250
balance 1.0 shift x 5 1.1 weight time 1.0 # out unweighted.txt
run 250

View File

@ -1,47 +1,47 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
units lj
atom_style atomic
processors * 1 1
variable factor index 1.0
variable factor index 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
group fast type 1
group slow type 2
balance 1.0 shift x 5 1.1 &
group fast type 1
group slow type 2
balance 1.0 shift x 5 1.1 &
weight group 2 fast 1.0 slow ${factor} # out weighted.txt
fix 0 all balance 10 1.0 shift x 5 1.1 &
fix 0 all balance 10 1.0 shift x 5 1.1 &
weight group 2 fast 1.0 slow ${factor}
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all atom 50 dump.melt
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
thermo 50
run 250
thermo 50
run 250

View File

@ -1,54 +1,54 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
units lj
atom_style atomic
processors * 1 1
variable factor index 1.0
variable factor index 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
balance 1.0 shift x 5 1.1 # out unweighted.txt
balance 1.0 shift x 5 1.1 # out unweighted.txt
balance 1.0 x uniform
balance 1.0 x uniform
variable weight atom (type==1)*1.0+(type==2)*v_factor
balance 1.0 shift x 5 1.1 weight var weight # out weighted_var.txt
balance 1.0 shift x 5 1.1 weight var weight # out weighted_var.txt
balance 1.0 x uniform
balance 1.0 x uniform
group fast type 1
group slow type 2
balance 1.0 shift x 5 1.1 &
group fast type 1
group slow type 2
balance 1.0 shift x 5 1.1 &
weight group 2 fast 1.0 slow ${factor} # out weighted_group.txt
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all atom 50 dump.melt
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
thermo 50
run 250
thermo 50
run 250

View File

@ -1,47 +1,47 @@
# 3d Lennard-Jones melt
units lj
#atom_style charge
units lj
#atom_style charge
processors * 1 1
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
#set type 1:2 charge 0.0
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/long/coul/long long off 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/long/coul/long long off 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
kspace_style pppm/disp 1.0e-4
kspace_modify gewald/disp 0.1
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
group fast type 1
group slow type 2
fix 0 all balance 20 1.0 shift x 5 1.0 &
group fast type 1
group slow type 2
fix 0 all balance 20 1.0 shift x 5 1.0 &
weight group 2 fast 1.0 slow 2.0 weight time 0.66
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all atom 50 dump.melt
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
thermo 50
run 500
thermo 50
run 500

View File

@ -1,53 +1,53 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
units lj
atom_style atomic
processors * 1 1
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
fix p all property/atom d_WEIGHT
fix 0 all balance 50 1.0 shift x 5 1.0 &
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
fix p all property/atom d_WEIGHT
fix 0 all balance 50 1.0 shift x 5 1.0 &
weight neigh 0.8 weight store WEIGHT
compute p all property/atom d_WEIGHT
variable maximb equal f_0[1]
variable iter equal f_0[2]
variable prev equal f_0[3]
variable final equal f_0
compute p all property/atom d_WEIGHT
variable maximb equal f_0[1]
variable iter equal f_0[2]
variable prev equal f_0[3]
variable final equal f_0
#fix 3 all print 50 "${iter} ${prev} ${final} ${maximb}"
#fix 3 all print 50 "${iter} ${prev} ${final} ${maximb}"
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all custom 50 dump.lammpstrj id type x y z c_p
#dump id all atom 50 dump.melt
#dump id all custom 50 dump.lammpstrj id type x y z c_p
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mp4 c_p type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3 amap 0.0 2.0 cf 0.1 3 min blue 0.5 green max red
#dump 3 all movie 25 movie.mp4 c_p type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3 amap 0.0 2.0 cf 0.1 3 min blue 0.5 green max red
thermo 50
run 500
run 500
thermo 50
run 500
run 500
run 500
run 500

View File

@ -1,53 +1,53 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
units lj
atom_style atomic
processors * 1 1
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
comm_style tiled
comm_style tiled
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
fix p all property/atom d_WEIGHT
fix 0 all balance 50 1.0 rcb weight neigh 0.8 weight store WEIGHT
compute p all property/atom d_WEIGHT
variable maximb equal f_0[1]
variable iter equal f_0[2]
variable prev equal f_0[3]
variable final equal f_0
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
fix p all property/atom d_WEIGHT
fix 0 all balance 50 1.0 rcb weight neigh 0.8 weight store WEIGHT
compute p all property/atom d_WEIGHT
variable maximb equal f_0[1]
variable iter equal f_0[2]
variable prev equal f_0[3]
variable final equal f_0
#fix 3 all print 50 "${iter} ${prev} ${final} ${maximb}"
#fix 3 all print 50 "${iter} ${prev} ${final} ${maximb}"
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all custom 50 dump.lammpstrj id type x y z c_p
#dump id all atom 50 dump.melt
#dump id all custom 50 dump.lammpstrj id type x y z c_p
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
thermo 50
thermo 50
run 250
run 250
run 250

View File

@ -1,51 +1,51 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
units lj
atom_style atomic
processors * 1 1
variable factor index 1.0
variable factor index 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
balance 1.0 shift x 10 1.0 weight neigh 0.8 # out weighted_var.txt
balance 1.0 shift x 10 1.0 weight neigh 0.8 # out weighted_var.txt
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all atom 50 dump.melt
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
thermo 50
run 250 post no
balance 1.0 shift x 10 1.0 weight neigh 0.8
thermo 50
run 250 post no
balance 1.0 shift x 10 1.0 weight neigh 0.8
balance 1.0 shift x 10 1.0 weight neigh 0.8
run 250 post no
balance 1.0 shift x 10 1.0 weight neigh 0.8 weight time 0.6
balance 1.0 shift x 10 1.0 weight neigh 0.8
run 250 post no
balance 1.0 shift x 10 1.0 weight neigh 0.8 weight time 0.6
run 250
balance 1.0 shift x 10 1.0 weight neigh 0.8 weight time 0.6
balance 1.0 shift x 10 1.0 weight neigh 0.8 weight time 0.6
run 250

View File

@ -1,65 +1,65 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
units lj
atom_style atomic
processors * 1 1
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 3 box
create_atoms 1 box
mass * 1.0
region long block 3 6 0 10 0 10
region long block 3 6 0 10 0 10
set region long type 2
velocity all create 1.0 87287
velocity all create 1.0 87287
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_style lj/cut 2.5
pair_coeff * * 1.0 1.0 2.5
pair_coeff * 2 1.0 1.0 5.0
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
fix p all property/atom d_WEIGHT
compute p all property/atom d_WEIGHT
neighbor 0.3 bin
neigh_modify every 2 delay 4 check yes
fix p all property/atom d_WEIGHT
compute p all property/atom d_WEIGHT
group fast type 1
group slow type 2
balance 1.0 shift x 10 1.0 &
group fast type 1
group slow type 2
balance 1.0 shift x 10 1.0 &
weight group 2 fast 0.8 slow 2.5 weight store WEIGHT
variable lastweight atom c_p
variable lastweight atom c_p
fix 0 all balance 50 1.0 shift x 10 1.0 &
fix 0 all balance 50 1.0 shift x 10 1.0 &
weight var lastweight weight time 0.5 weight store WEIGHT
variable maximb equal f_0[1]
variable iter equal f_0[2]
variable prev equal f_0[3]
variable final equal f_0
variable maximb equal f_0[1]
variable iter equal f_0[2]
variable prev equal f_0[3]
variable final equal f_0
#fix 3 all print 50 "${iter} ${prev} ${final} ${maximb}"
#fix 3 all print 50 "${iter} ${prev} ${final} ${maximb}"
fix 1 all nve
fix 1 all nve
#dump id all atom 50 dump.melt
#dump id all custom 50 dump.lammpstrj id type x y z c_p
#dump id all atom 50 dump.melt
#dump id all custom 50 dump.lammpstrj id type x y z c_p
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
thermo 50
run 500
run 500
thermo 50
run 500
run 500
balance 1.0 shift x 10 1.0 &
balance 1.0 shift x 10 1.0 &
weight group 2 fast 0.8 slow 2.5 weight store WEIGHT
fix 0 all balance 50 1.0 shift x 5 1.0 &
fix 0 all balance 50 1.0 shift x 5 1.0 &
weight var lastweight weight neigh 0.5 weight store WEIGHT
run 500

View File

@ -1,23 +1,23 @@
# GayBerne ellipsoids in LJ background fluid
units lj
units lj
atom_style ellipsoid
dimension 2
lattice sq 0.02
region box block 0 20 0 20 -0.5 0.5
lattice sq 0.02
region box block 0 20 0 20 -0.5 0.5
create_box 2 box
create_atoms 1 box
set group all type/fraction 2 0.1 95392
set type 1 mass 1.0
set type 2 mass 1.5
set type 1 shape 1 1 1
set type 2 shape 3 1 1
set group all quat/random 18238
set group all type/fraction 2 0.1 95392
set type 1 mass 1.0
set type 2 mass 1.5
set type 1 shape 1 1 1
set type 2 shape 3 1 1
set group all quat/random 18238
compute rot all temp/asphere
group spheroid type 1
compute rot all temp/asphere
group spheroid type 1
variable dof equal count(spheroid)+2
compute_modify rot extra/dof ${dof}
@ -31,36 +31,37 @@ pair_coeff 2 2 1.0 1.0 1 1 0.2 0 0 0
neighbor 0.8 bin
thermo_style custom step c_rot epair etotal press vol
thermo 100
thermo 100
timestep 0.002
compute q all property/atom quatw quati quatj quatk
compute 0 all property/atom quatw quati quatj quatk shapex shapey shapez
#dump 1 all custom 100 dump.ellipse.gayberne &
# id type x y z c_q[1] c_q[2] c_q[3] c_q[4]
dump 1 all custom 100 dump.ellipse.gayberne id type x y z c_0[*]
dump_modify 1 colname c_0[1] quatw colname c_0[2] quati colname c_0[3] quatj colname c_0[4] quatk &
colname c_0[5] shapex colname c_0[6] shapey colname c_0[7] shapez
#dump 2 all image 100 image.*.jpg type type &
# zoom 1.6 center d 0.5 0.5 0.5
#dump 2 all image 100 image.*.jpg type type &
# zoom 1.6 center d 0.5 0.5 0.5
#dump_modify 2 pad 4 adiam 1 1.0 adiam 2 2.0
#dump 3 all movie 100 movie.mpg type type &
# zoom 1.6 center d 0.5 0.5 0.5
#dump 3 all movie 100 movie.mpg type type &
# zoom 1.6 center d 0.5 0.5 0.5
#dump_modify 3 pad 4 adiam 1 1.0 adiam 2 2.0
fix 1 all npt/asphere temp 2.0 2.0 0.1 iso 0.0 1.0 1.0 &
mtk no pchain 0 tchain 1
fix 2 all enforce2d
fix 1 all npt/asphere temp 2.0 2.0 0.1 iso 0.0 1.0 1.0 &
mtk no pchain 0 tchain 1
fix 2 all enforce2d
compute_modify 1_temp extra/dof ${dof}
# equilibrate to shrink box around dilute system
run 2000
run 2000
# run dynamics on dense system
unfix 1
fix 1 all nve/asphere
unfix 1
fix 1 all nve/asphere
run 2000
run 2000

View File

@ -1,23 +1,23 @@
# RESquared ellipsoids in LJ background fluid
units lj
units lj
atom_style ellipsoid
dimension 2
lattice sq 0.02
region box block 0 20 0 20 -0.5 0.5
lattice sq 0.02
region box block 0 20 0 20 -0.5 0.5
create_box 2 box
create_atoms 1 box
set group all type/fraction 2 0.1 95392
set type 1 mass 1.0
set type 2 mass 1.5
set type 1 shape 1 1 1
set type 2 shape 3 1 1
set group all quat/random 18238
set group all type/fraction 2 0.1 95392
set type 1 mass 1.0
set type 2 mass 1.5
set type 1 shape 1 1 1
set type 2 shape 3 1 1
set group all quat/random 18238
compute rot all temp/asphere
group spheroid type 1
compute rot all temp/asphere
group spheroid type 1
variable dof equal count(spheroid)+2
compute_modify rot extra/dof ${dof}
@ -31,36 +31,40 @@ pair_coeff 2 2 1.0 1 1 1 0.2 0 0 0
neighbor 0.8 bin
thermo_style custom step c_rot epair etotal press vol
thermo 100
thermo 100
timestep 0.002
compute q all property/atom quatw quati quatj quatk
compute 0 all property/atom quatw quati quatj quatk shapex shapey shapez
#dump 1 all custom 100 dump.ellipse.resquared &
# id type x y z c_q[1] c_q[2] c_q[3] c_q[4]
dump 1 all custom 100 dump.ellipse.resquared id type x y z c_0[*]
dump_modify 1 colname c_0[1] quatw colname c_0[2] quati colname c_0[3] quatj colname c_0[4] quatk &
colname c_0[5] shapex colname c_0[6] shapey colname c_0[7] shapez
#dump 2 all image 100 image.*.jpg type type &
# zoom 1.6 center d 0.5 0.5 0.5
#dump 1 all custom 100 dump.ellipse.resquared &
# id type x y z c_q[1] c_q[2] c_q[3] c_q[4]
#dump 2 all image 100 image.*.jpg type type &
# zoom 1.6 center d 0.5 0.5 0.5
#dump_modify 2 pad 4 adiam 1 1.0 adiam 2 2.0
#dump 3 all movie 100 movie.mpg type type &
# zoom 1.6 center d 0.5 0.5 0.5
#dump 3 all movie 100 movie.mpg type type &
# zoom 1.6 center d 0.5 0.5 0.5
#dump_modify 3 pad 4 adiam 1 1.0 adiam 2 2.0
fix 1 all npt/asphere temp 2.0 2.0 0.1 iso 0.0 1.0 1.0 &
mtk no pchain 0 tchain 1
fix 2 all enforce2d
fix 1 all npt/asphere temp 2.0 2.0 0.1 iso 0.0 1.0 1.0 &
mtk no pchain 0 tchain 1
fix 2 all enforce2d
compute_modify 1_temp extra/dof ${dof}
# equilibrate to shrink box around dilute system
run 2000
run 2000
# run dynamics on dense system
unfix 1
fix 1 all nve/asphere
unfix 1
fix 1 all nve/asphere
run 2000
run 2000

View File

@ -3,35 +3,35 @@
# rho ~ 0.5
# p ~ 1.5
# mu_ex ~ 0.0
# comparable to Frenkel and Smit GCMC Case Study, Figure 5.8
# comparable to Frenkel and Smit GCMC Case Study, Figure 5.8
# variables modifiable using -var command line switch
variable mu index -1.25
variable temp index 2.0
variable disp index 1.0
variable disp index 1.0
variable lbox index 5.0
# global model settings
units lj
atom_style atomic
pair_style lj/cut 3.0
pair_modify tail no # turn of to avoid triggering full_energy
pair_style lj/cut 3.0
pair_modify tail no # turn of to avoid triggering full_energy
# box
region box block 0 ${lbox} 0 ${lbox} 0 ${lbox}
create_box 1 box
region box block 0 ${lbox} 0 ${lbox} 0 ${lbox}
create_box 1 box
# lj parameters
pair_coeff * * 1.0 1.0
mass * 1.0
pair_coeff * * 1.0 1.0
mass * 1.0
# we recommend setting up a dedicated group for gcmc
group gcmcgroup type 1
group gcmcgroup type 1
# gcmc
@ -39,29 +39,29 @@ fix mygcmc gcmcgroup gcmc 1 100 100 1 29494 ${temp} ${mu} ${disp}
# atom count
variable type1 atom "type==1"
group type1 dynamic gcmcgroup var type1
variable type1 atom "type==1"
group type1 dynamic gcmcgroup var type1
variable n1 equal count(type1)
# averaging
variable rho equal density
variable p equal press
variable nugget equal 1.0e-8
variable rho equal density
variable p equal press
variable nugget equal 1.0e-8
variable lambda equal 1.0
variable muex equal ${mu}-${temp}*ln(density*${lambda}+${nugget})
fix ave all ave/time 10 100 1000 v_rho v_p v_muex v_n1 ave one file rho_vs_p.dat
variable rhoav equal f_ave[1]
variable pav equal f_ave[2]
variable muexav equal f_ave[3]
variable n1av equal f_ave[4]
variable muex equal ${mu}-${temp}*ln(density*${lambda}+${nugget})
fix ave all ave/time 10 100 1000 v_rho v_p v_muex v_n1 ave one file rho_vs_p.dat
variable rhoav equal f_ave[1]
variable pav equal f_ave[2]
variable muexav equal f_ave[3]
variable n1av equal f_ave[4]
# output
variable tacc equal f_mygcmc[2]/(f_mygcmc[1]+${nugget})
variable iacc equal f_mygcmc[4]/(f_mygcmc[3]+${nugget})
variable dacc equal f_mygcmc[6]/(f_mygcmc[5]+${nugget})
compute_modify thermo_temp dynamic yes
variable tacc equal f_mygcmc[2]/(f_mygcmc[1]+${nugget})
variable iacc equal f_mygcmc[4]/(f_mygcmc[3]+${nugget})
variable dacc equal f_mygcmc[6]/(f_mygcmc[5]+${nugget})
compute_modify thermo_temp dynamic/dof yes
thermo_style custom step temp press pe ke density atoms v_iacc v_dacc v_tacc v_rhoav v_pav v_muexav v_n1av
thermo 1000

View File

@ -0,0 +1,87 @@
# Running JAX from LAMMPS
### Getting started
First make a Python environment with dependencies:
conda create --name jax python=3.10
conda activate jax
# Upgrade pip
python -m pip install --upgrade pip
# Install JAX:
python -m pip install --upgrade "jax[cpu]"
# Install other dependencies:
python -m pip install numpy scipy torch scikit-learn virtualenv psutil tabulate mpi4py Cython
Install LAMMPS:
cd /path/to/lammps
mkdir build-jax; cd build-jax
cmake ../cmake -DLAMMPS_EXCEPTIONS=yes \
-DBUILD_SHARED_LIBS=yes \
-DMLIAP_ENABLE_PYTHON=yes \
-DPKG_PYTHON=yes \
-DPKG_ML-SNAP=yes \
-DPKG_ML-IAP=yes \
-DPYTHON_EXECUTABLE:FILEPATH=`which python`
make -j4
make install-python
### Kokkos install
Use same Python dependencies as above, with some extra changes:
1. Make sure you install cupy properly! E.g.
python -m pip install cupy-cuda12x
2. Install JAX for GPU/CUDA:
python -m pip install --trusted-host storage.googleapis.com --upgrade "jax[cuda12_local]" -f https://storage.googleapis.com/jax-releases/jax_cuda_releases.html
3. Install cudNN: https://developer.nvidia.com/cudnn
Install LAMMPS. Take care to change `Kokkos_ARCH_*` flag:
cmake ../cmake -DLAMMPS_EXCEPTIONS=yes \
-DBUILD_SHARED_LIBS=yes \
-DPKG_PYTHON=yes \
-DPKG_ML-SNAP=yes \
-DPKG_ML-IAP=yes \
-DMLIAP_ENABLE_PYTHON=yes \
-DPKG_KOKKOS=yes \
-DKokkos_ARCH_TURING75=yes \
-DKokkos_ENABLE_CUDA=yes \
-DKokkos_ENABLE_OPENMP=yes \
-DCMAKE_CXX_COMPILER=${HOME}/lammps/lib/kokkos/bin/nvcc_wrapper \
-DPYTHON_EXECUTABLE:FILEPATH=`which python`
make -j
make install-python
Run example:
mpirun -np 1 lmp -k on g 1 -sf kk -pk kokkos newton on -in in.run
### Deploying JAX models on CPU
Use `deploy_script.py`, which will wrap model with `write_unified_jax`.
python deploy_script.py
This creates `.pkl` file to be loaded by LAMMPS ML-IAP Unified.
Run LAMMPS with the model:
mpirun -np P lmp -in in.run
### Deploying JAX models in Kokkos
Use `deploy_script_kokkos.py`, which will wrap model with `write_unified_jax_kokkos`.
python deploy_script_kokkos.py
This creates `.pkl` file to be loaded by LAMMPS ML-IAP Unified.
Run LAMMPS with the model:
mpirun -np 1 lmp -k on g 1 -sf kk -pk kokkos newton on -in in.run

View File

@ -0,0 +1,11 @@
import lammps
import lammps.mliap
#from lammps.mliap.mliap_unified_lj import MLIAPUnifiedLJ
from mliap_unified_jax import MLIAPUnifiedJAX
def create_pickle():
unified = MLIAPUnifiedJAX(["Ar"])
unified.pickle('mliap_unified_jax_Ar.pkl')
create_pickle()

37
examples/mliap/jax/in.run Normal file
View File

@ -0,0 +1,37 @@
# 3d Lennard-Jones melt
units lj
atom_style atomic
lattice fcc 0.8442
region box block 0 10 0 10 0 10
create_box 1 box
create_atoms 1 box
mass 1 1.0
velocity all create 3.0 87287 loop geom
pair_style mliap unified mliap_unified_jax_Ar.pkl 0
pair_coeff * * Ar
neighbor 0.3 bin
neigh_modify every 20 delay 0 check no
fix 1 all nve
#dump id all atom 50 dump.melt
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 1 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
#dump 4 all custom 1 forces.xyz fx fy fz
dump 1 all xyz 10 dump.xyz
thermo 1
run 250

Binary file not shown.

View File

@ -0,0 +1,61 @@
from lammps.mliap.mliap_unified_abc import MLIAPUnified
import numpy as np
import jax
import jax.numpy as jnp
from jax import jit
from functools import partial
import os
# Required else get `jaxlib.xla_extension.XlaRuntimeError: RESOURCE_EXHAUSTED: Out of memory`
os.environ["XLA_PYTHON_CLIENT_PREALLOCATE"]="false"
os.environ["XLA_PYTHON_CLIENT_MEM_FRACTION"]=".XX"
os.environ["XLA_PYTHON_CLIENT_ALLOCATOR"]="platform"
@jax.jit
def lj_potential(epsilon, sigma, rij):
def _tot_e(rij):
"""A differentiable fn for total energy."""
r2inv = 1.0 / jnp.sum(rij ** 2, axis=1)
r6inv = r2inv * r2inv * r2inv
lj1 = 4.0 * epsilon * sigma**12
lj2 = 4.0 * epsilon * sigma**6
eij = r6inv * (lj1 * r6inv - lj2)
return 0.5 * jnp.sum(eij), eij
# Compute _tot_e and its derivative.
(_, eij), fij = jax.value_and_grad(_tot_e, has_aux=True)(rij)
return eij, fij
class MLIAPUnifiedJAX(MLIAPUnified):
"""Test implementation for MLIAPUnified."""
epsilon: float
sigma: float
def __init__(self, element_types, epsilon=1.0, sigma=1.0, rcutfac=1.25):
# ARGS: interface, element_types, ndescriptors, nparams, rcutfac
super().__init__(None, element_types, 1, 3, rcutfac)
# Mimicking the LJ pair-style:
# pair_style lj/cut 2.5
# pair_coeff * * 1 1
self.epsilon = epsilon
self.sigma = sigma
def compute_gradients(self, data):
"""Test compute_gradients."""
def compute_descriptors(self, data):
"""Test compute_descriptors."""
def compute_forces(self, data):
"""Test compute_forces."""
# NOTE: Use data.rij_max with JAX.
rij = data.rij_max
eij, fij = lj_potential(self.epsilon, self.sigma, rij)
data.update_pair_energy(np.array(eij, dtype=np.float64))
data.update_pair_forces(np.array(fij, dtype=np.float64))

Binary file not shown.

View File

@ -0,0 +1,69 @@
from lammps.mliap.mliap_unified_abc import MLIAPUnified
import numpy as np
import jax
import jax.dlpack
import jax.numpy as jnp
from jax import jit
from functools import partial
import cupy
import os
# Required else get `jaxlib.xla_extension.XlaRuntimeError: RESOURCE_EXHAUSTED: Out of memory`
# Does not fix GPU problem with larger num. atoms.
#os.environ["XLA_PYTHON_CLIENT_PREALLOCATE"]="false"
#os.environ["XLA_PYTHON_CLIENT_MEM_FRACTION"]=".XX"
#os.environ["XLA_PYTHON_CLIENT_ALLOCATOR"]="platform"
@jax.jit
def lj_potential(epsilon, sigma, rij):
# A pure function we can differentiate:
def _tot_e(rij):
r2inv = 1.0 / jnp.sum(rij ** 2, axis=1)
r6inv = r2inv * r2inv * r2inv
lj1 = 4.0 * epsilon * sigma**12
lj2 = 4.0 * epsilon * sigma**6
eij = r6inv * (lj1 * r6inv - lj2)
return 0.5 * jnp.sum(eij), eij
# Construct a function computing _tot_e and its derivative
(_, eij), fij = jax.value_and_grad(_tot_e, has_aux=True)(rij)
return eij, fij
class MLIAPUnifiedJAXKokkos(MLIAPUnified):
"""JAX wrapper for MLIAPUnified."""
epsilon: float
sigma: float
def __init__(self, element_types, epsilon=1.0, sigma=1.0, rcutfac=1.25):
# ARGS: interface, element_types, ndescriptors, nparams, rcutfac
super().__init__(None, element_types, 1, 3, rcutfac)
# Mimicking the LJ pair-style:
# pair_style lj/cut 2.5
# pair_coeff * * 1 1
self.epsilon = epsilon
self.sigma = sigma
def compute_gradients(self, data):
"""Test compute_gradients."""
def compute_descriptors(self, data):
"""Test compute_descriptors."""
def compute_forces(self, data):
"""Test compute_forces."""
# NOTE: Use data.rij_max with JAX.
# dlpack requires cudnn:
rij = jax.dlpack.from_dlpack(data.rij_max.toDlpack())
eij, fij = lj_potential(self.epsilon, self.sigma, rij)
# Convert back to cupy.
eij = cupy.from_dlpack(jax.dlpack.to_dlpack(eij)).astype(np.float64)
fij = cupy.from_dlpack(jax.dlpack.to_dlpack(fij)).astype(np.float64)
# Send to LAMMPS.
data.update_pair_energy(eij)
data.update_pair_forces(fij)

View File

@ -0,0 +1,87 @@
"""
interface for creating LAMMPS MLIAP Unified models.
"""
import pickle
import numpy as np
from lammps.mliap.mliap_unified_abc import MLIAPUnified
#from deploy_script import MyModel
class MLIAPInterface(MLIAPUnified):
"""
Class for creating ML-IAP Unified model based on hippynn graphs.
"""
def __init__(self, model, element_types, cutoff=4.5, ndescriptors=1):
"""
:param model: class defining the model
:param element_types: list of atomic symbols corresponding to element types
:param ndescriptors: the number of descriptors to report to LAMMPS
:param model_device: the device to send torch data to (cpu or cuda)
"""
super().__init__()
self.model = model
self.element_types = element_types
self.ndescriptors = ndescriptors
#self.model_device = model_device
# Build the calculator
# TODO: Make this cutoff depend on model cutoff, ideally from deployed model itself but could
# be part of deploy step.
#rc = 4.5
self.rcutfac = 0.5*cutoff # Actual cutoff will be 2*rc
#print(self.model.nparams)
self.nparams = 10
#self.rcutfac, self.species_set, self.graph = setup_LAMMPS()
#self.nparams = sum(p.nelement() for p in self.graph.parameters())
#self.graph.to(torch.float64)
def compute_descriptors(self, data):
pass
def compute_gradients(self, data):
pass
def compute_forces(self, data):
#print(">>>>> hey!")
#elems = self.as_tensor(data.elems).type(torch.int64).reshape(1, data.ntotal)
"""
elems = self.as_tensor(data.elems).type(torch.int64) + 1
#z_vals = self.species_set[elems+1]
pair_i = self.as_tensor(data.pair_i).type(torch.int64)
pair_j = self.as_tensor(data.pair_j).type(torch.int64)
rij = self.as_tensor(data.rij).type(torch.float64).requires_grad_(True)
nlocal = self.as_tensor(data.nlistatoms)
"""
rij = data.rij
#(total_energy, fij) = self.network(rij, None, None, None, nlocal, elems, pair_i, pair_j, "cpu", dtype=torch.float64, mode="lammps")
test = self.model(rij)
#data.update_pair_forces(fij)
#data.energy = total_energy.item()
pass
def setup_LAMMPS(energy):
"""
:param energy: energy node for lammps interface
:return: graph for computing from lammps MLIAP unified inputs.
"""
model = TheModelClass(*args, **kwargs)
save_state_dict = torch.load("Ta_Pytorch.pt")
model.load_state_dict(save_state_dict["model_state_dict"])
#model.load_state_dict(torch.load(PATH))
model.eval()
#model.eval()
return model