Dear all,
I am trying to do a polymer simulation (20000 atoms) where I have an additional set of forces in a text (eg: fx_lammps.txt which contains 20000 different numerical force values each corresponding to the force acting on each atom) file. I would like to add these extra forces into the polymer beads as an additional force. I used the “variable file” command to read the text file followed by the “variable atom” to store forces on a per-atom basis. As a starting point I tried to use fix setforce so that all other forces will be nullified. I used “dump fx” to output the forces acting on each atom.
Problem was that the force that is acting on each atom in every time step is the same and it is the first value of force in the file. It is not assigning force values in the file to its corresponding atom. Only on the next time step the second line is read from the file and the second value is assigned to all the atoms for the second run (similarly 3rd value in 3rd run and so on). It would be really helpful if I can know whether in each time step, I can use the 20000 force values (each force value corresponding to each atom) using setforce/addforce command. I am attaching the code below. Am I missing here something really fundamental?
############# Variable Declaration ############################
variable fx file fx_lammps.txt #fx_lammps is the file name which contains 20000 force values
variable fy file fy_lammps.txt
variable fz file fz_lammps.txt
variable fxval atom next(fx) # Here I tried to store the value read from file to an atom-style variable
variable fyval atom next(fy)
variable fzval atom next(fz)
############ Initialization ###################################
units real
boundary p p p
atom_style molecular
log log.{simuname}.txt
read_data {dataname}
############ LJ Potential Information #########################
neighbor 0.02 bin
bond_style harmonic
bond_coeff * $k 0.387
pair_style lj/cut 1.0
pair_coeff * * 0.005 0.387
############ Equilibration Stage 1 (NVT Dynamics) #############
velocity all create 300 67233 dist gaussian
fix nvt_run all nvt temp 300.0 300.0 100.0 # I tried without fix nvt as well, but it gave the same result
fix extfor all setforce v_fxval v_fyval v_fzval
########### Output Data ####################################
dump coord all custom 10 coord.txt id x y z
dump force all custom 1 fatom.txt id fx fy fz
dump_modify coord sort id
dump_modify force sort id
run 10
Vaidyanathan M S
PhD Scholar
Department of Chemical Engineering
The University of Texas at Austin