Hello everyone,
I am trying to shear my microgel and I would like to calculate its gyration tensor components along the way as well as the intrinsic viscosity and I can’t seem to get it to work:
LAMMPS script for polymer gel under shear flow
Box and units (LJ units and periodic boundaries)
units lj
atom_style bond
boundary p p p
Neighbor list settings
neighbor 0.3 bin
neigh_modify every 1 delay 1 check yes
RESTART file generation
restart 100000 mgel.restart.a mgel.restart.b
READ initial configuration
read_data data_final_alfa0
reset_timestep 0
Convert box to triclinic for shear deformation
change_box all triclinic
Define groups
group all type 1
Dump configurations
dump mgelconf all custom 50000 mgel.dump.* id type x y z vx vy vz ix iy iz
dump_modify mgelconf format line “%d d .10f .10f .10f .10f .10f %.10f %d %d %d”
dump_modify mgelconf pbc yes
dump_modify mgelconf sort id
Pair interactions
pair_style lj/cut 1.122462048309
pair_modify shift yes
pair_coeff * * 1.0 1.0
Bond interactions (FENE)
bond_style fene
special_bonds fene
bond_coeff 1 30 1.5 1 1
Shear flow setup using SLLOD algorithm
Shear rate = 0.001 in LJ units
variable shear_rate equal 0.001
Fixes for SLLOD integration
First: deform the box at constant shear rate
fix deform all deform 1 xy erate ${shear_rate} remap v
Second: Use SLLOD equations with NVT thermostat
The “temp” thermostat acts on thermal velocities relative to shear flow
fix sllod all nvt/sllod temp 1.0 1.0 0.2 tchain 30
Compute for stress tensor components
compute stress all stress/atom NULL virial
compute p all pressure NULL temp c_stress
compute p all pressure temp c_stress
compute p all pressure NULL
Compute gyration tensor and related properties
compute gyrtensor all gyration tensor
compute rg all gyration
Extract individual components of gyration tensor
variable Rg2 equal c_rg
variable Rgxx equal c_gyrtensor[1]
variable Rgyy equal c_gyrtensor[2]
variable Rgzz equal c_gyrtensor[3]
variable Rgxy equal c_gyrtensor[4]
Compute intrinsic viscosity
σ_xy = -1/V * ∑(R_ix * F_iy)
Intrinsic viscosity = σ_xy / (shear_rate)
variable box_vol equal vol
variable sigma_xy equal -(c_p[4]) # c_p[4] is the xy component of stress tensor
variable intrinsic_viscosity equal v_sigma_xy/${shear_rate}
Additional computes for monitoring
compute ke all ke
compute pe all pe
compute temp_all all temp
Density profile calculation
compute densprof all chunk/atom bin/sphere 0 0 0 0.0 120.0 120
fix fordensprof all ave/chunk 2 20 10000 densprof density/number file densprofmgel
Recentering fix
fix recentermgel all recenter 0 0 0 shift all units box
Thermodynamic output with enhanced observables
thermo 1000
thermo_style custom step temp etotal pe ke press c_rg v_Rgxx v_Rgyy v_Rgzz v_Rgxy v_sigma_xy v_intrinsic_viscosity
Timestep
timestep 0.002
Minimization
minimize 1.0e-4 1.0e-6 1000 1000
Production run
run 50000000
Write final configuration
write_data data_final_sheared
I am happy for any feedback, thank you in advance!!
Leah
PS.: I also later on would like to shear a bulk of tighly packed microgels, do I simply remove the recentering command and lammps takes proper care of it (Lees-Edwards-boundary conditions, …) or is it not possible this way. Right now with one microgel I can recenter it and no boundaries are in the way but I am aware that this is not the case in bulk simulations.