Dear Lammps experts,
I am simulating impact of a particle on a similar type substrate with LAMMPS . my problem is about calculating the particle temperature as a moving object actually.
First the initial structure of the particle and its energy minimize at the room temperature(300K) with fix nve and fix temp/rescale and then the particle is thrown toward the substrate. I applied the compute temp/com in order to calculate the particle temperature as a moving object as per the Lammps document.
These are the problems:
1- Exactly upon throwing the particle ( after unfix temp/rescale) the particle temperature decreases about 150K which is not expected. I have no idea why temperature decreases TOO MUCH upon unfixing. How about minimizing the particle at the room temperature longer???
2- and then the temperature increases until the end of the simulation. However, it is expected that the temp. only increases at the impact moment and then decreases to a steady or constant temp at the end of the simulation.
So, any comment could be highly appreciated.
below you can find my input file:
#Phase 1 ------------------------------------------Simulation main setup----------------------------------
dimension 3
units real
atom_style charge
boundary p p p
variable radius equal 50 # makes 100 Ang or 10nm particle
variable n equal 50 # makes 81*3.786 length=306
variable a1 equal 3.786
variable a2 equal 3.786
variable a3 equal 9.514
variable boxx equal “v_nv_a1"
variable boxy equal "v_nv_a2”
variable diameter equal “2v_radius"
variable x0 equal "0.5v_boxx”
variable y0 equal “0.5v_boxy"
variable subs_thick equal "3v_a3” #v_a3 equal 9.514
variable subslow_thick equal “1v_a3"
variable z_gap equal "0.5v_radius”
variable distance equal “v_subs_thick+v_subslow_thick+v_z_gap+v_radius”
variable boxz equal “v_distance+v_diameter+v_radius”
#simulation box
region box block 0 {boxx} 0 {boxy} 0 ${boxz} units box
create_box 2 box #2 is number of atoms
it was very difficult to finally find these positions
lattice custom 1 a1 3.786 0.00000 0.00000 a2 0.0000 3.786 0.00000 a3 0 0 9.514 &
basis 0.0000 0.2500 0.3750 &
basis 0.0000 0.7500 0.6250 &
basis 0.5000 0.7500 0.8750 &
basis 0.5000 0.2500 0.1250 &
basis 0.0000 0.0000 0.1700 &
basis 0.0000 0.7500 0.4200 &
basis 0.5000 0.2500 0.3300 &
basis 0.5000 0.7500 0.0800 &
basis 0.500 0.5000 0.6700 &
basis 0.500 0.2500 0.9200 &
basis 0.000 0.7500 0.8300 &
basis 0.000 0.2500 0.5800
mass 1 47.86000
mass 2 15.99940
#particle
region particle sphere {x0} {y0} {distance} {radius} units box
create_atoms 2 region particle &
basis 1 1 &
basis 2 1 &
basis 3 1 &
basis 4 1 &
basis 5 2 &
basis 6 2 &
basis 7 2 &
basis 8 2 &
basis 9 2 &
basis 10 2 &
basis 11 2 &
basis 12 2
group particle region particle
set type 1 charge 2.196
set type 2 charge -1.098
#substrate
region substrate block 0 {boxx} 0 {boxy} {subslow_thick} {subs_thick} units box
create_atoms 2 region substrate &
basis 1 1 &
basis 2 1 &
basis 3 1 &
basis 4 1 &
basis 5 2 &
basis 6 2 &
basis 7 2 &
basis 8 2 &
basis 9 2 &
basis 10 2 &
basis 11 2 &
basis 12 2
group substrate region substrate
set type 1 charge 2.196
set type 2 charge -1.098
#group model union particle substrate
#lower_substrate
region lower_substrate block 0 {boxx} 0 {boxy} 0 ${subslow_thick} units box
create_atoms 2 region lower_substrate &
basis 1 1 &
basis 2 1 &
basis 3 1 &
basis 4 1 &
basis 5 2 &
basis 6 2 &
basis 7 2 &
basis 8 2 &
basis 9 2 &
basis 10 2 &
basis 11 2 &
basis 12 2
group lower_substrate region lower_substrate
set type 1 charge 2.196
set type 2 charge -1.098
group model union particle substrate
#–Phase 2----------------------------------------Buckingham Potential-----------------------------------------------
pair_style buck/coul/long 15
pair_coeff 1 1 717647.40 0.154 121.067
pair_coeff 1 2 391049.10 0.194 290.331
pair_coeff 2 2 271716.30 0.234 696.888
neighbor 2.0 bin # skin distance for real units is by default 2.0
neigh_modify every 1 delay 0 check yes
kspace_style pppm 0.0001
neigh_modify delay 5
#—Phase 3----------------------------------------Compute—how can I monitor Vparticle!!! ------------------------
compute temp_substrate all temp/region substrate
compute temp_particle_com particle temp/com
compute temp_particle all temp/region particle
#----Phase 4-------------------------------------Initial Equilibration at 300K ----------------------------------------
reset_timestep 0
timestep 1.0 # or 2
velocity all create 300 12345 mom yes rot no
fix 1 lower_substrate setforce 0.0 0.0 0.0
fix 2 particle nve
fix 3 substrate nve
fix 4 substrate temp/rescale 100 300.0 300.0 10.0 1.0
fix 5 particle temp/rescale 100 300.0 300.0 10.0 1.0
thermo 100
thermo_modify lost ignore flush yes
thermo_style custom step c_temp_particle c_temp_particle_com c_temp_substrate
dump 1 all xyz 100 dump1000.txt
run 7000
unfix 5
#----Phase 5---------------------------------------Particle Impact at the room temperature -------------------
velocity particle set 0 0 -0.010 units box
thermo 100
thermo_style custom step c_temp_particle c_temp_particle_com
run 50000
Best regards
Bahman Daneshian