access forces due to force-fields only

Dear Users,

I am running a simulation in which a BN-nanotube is subjected both to internal forces due to the adopted force-fields (tersoff and ILP/graphene/hBN) and to external forces that I apply using the fix commands:

#Add viscous forces to B and N atoms in the nanotube com frame of reference
fix visB Boron viscous {vissubB} fix visN Nitrogen viscous {vissubN}
fix fadsubB Boron addforce v_fdampBcmx 0.0 0.0
fix fadsubN Nitrogen addforce v_fdampNcmx 0.0 0.0

#Add external driving to the nanotube
fix fadtube nanotube addforce v_fext 0.0 0.0

I want to keep the x-velocity of the com ‘vcmtubex’ of the nanotube constant. A way of doing this is to compute the external driving (apart from coefficients) like this:

variable fext equal (${vtarget}-v_vcmtubex)-v/c_Ftotinternalx

Where ‘Ftotinternalx’ is the total force acting on the nanotube in the x-direction, due only to the internal forces.

In general, is there a way in LAMMPS to access the total force due to the force-fields only, before external forces are applied?

Thanks in advance,
Davide

Dear Users,

I am running a simulation in which a BN-nanotube is subjected both to internal forces due to the adopted force-fields (tersoff and ILP/graphene/hBN) and to external forces that I apply using the fix commands:

#Add viscous forces to B and N atoms in the nanotube com frame of reference
fix visB Boron viscous \{vissubB\} fix visN Nitrogen viscous {vissubN}
fix fadsubB Boron addforce v_fdampBcmx 0.0 0.0
fix fadsubN Nitrogen addforce v_fdampNcmx 0.0 0.0

#Add external driving to the nanotube
fix fadtube nanotube addforce v_fext 0.0 0.0

I want to keep the x-velocity of the com 'vcmtubex' of the nanotube constant. A way of doing this is to compute the external driving (apart from coefficients) like this:

as simpler way of keeping the velocity constant is to just set an
initial velocity and then set the force to zero. e.g. with:

velocity nanotube set 0.1 NULL NULL
fix fadtube nanotube setforce 0.0 NULL NULL

without force, the velocity does not change. what you are doing seems
needlessly complex. is there a reason you cannot use this simple
method?

axel.

Dear Axel,

Thanks for the quick response!

as simpler way of keeping the velocity constant is to just set an
initial velocity and then set the force to zero. e.g. with:

velocity nanotube set 0.1 NULL NULL
fix fadtube nanotube setforce 0.0 NULL NULL

without force, the velocity does not change. what you are doing seems
needlessly complex. is there a reason you cannot use this simple
method?

Yes, the quantity I am interested in is exactly that external force fext, which after the first time step equals the internal force contribution only.

Dear Axel,

Thanks for the quick response!

as simpler way of keeping the velocity constant is to just set an
initial velocity and then set the force to zero. e.g. with:

velocity nanotube set 0.1 NULL NULL
fix fadtube nanotube setforce 0.0 NULL NULL

without force, the velocity does not change. what you are doing seems
needlessly complex. is there a reason you cannot use this simple
method?

Yes, the quantity I am interested in is exactly that external force fext, which after the first time step equals the internal force contribution only.

sorry, but i am confused here. i don't think you are answering my
question, since there *is* no external force in what i am proposing.

axel.

The system is a nanotube sliding over a rigid substrate, I missed this information, sorry for that.
The total internal force acting on the nanotube is given by (1) the sum of all forces due to the boron-nitrogen interactions within the nanotube plus (2) the sum of all forces from the interactions of the nanotube with the substrate.
Contribution (1) is zero. Contribution (2) is what I am interested in monitoring, and I would like to do it at fixed center-of-mass velocity of the nanotube.
Now the picture should be more clear.

If I understand the script you are suggesting, the force acting along x on each atom of the nanotube would be set to zero.
The result of a compute command like “compute FFforcex nanotube reduce ave fx” would be “FFforcex=0”.
So I would still have no access to the internal force during simulation (though I could obtain it by postprocessing the trajectory).

Please correct me if I am wrong.

In any case, your first command “velocity nanotube set 0.1 NULL NULL” is setting the velocity along x of all atoms in the nanotube to a constant value.
This will not allow the nanotube to deform along x. Moreover, the nanotube may want to roll instead of just sliding over the substrate.

These are all motivations why I am trying to implement the protocol I reported in my first email.

Best,
Davide

The system is a nanotube sliding over a rigid substrate, I missed this information, sorry for that.
The total internal force acting on the nanotube is given by (1) the sum of all forces due to the boron-nitrogen interactions within the nanotube plus (2) the sum of all forces from the interactions of the nanotube with the substrate.
Contribution (1) is zero. Contribution (2) is what I am interested in monitoring, and I would like to do it at fixed center-of-mass velocity of the nanotube.
Now the picture should be more clear.

If I understand the script you are suggesting, the force acting along x on each atom of the nanotube would be set to zero.
The result of a compute command like "compute FFforcex nanotube reduce ave fx" would be "FFforcex=0".
So I would still have no access to the internal force during simulation (though I could obtain it by postprocessing the trajectory).
Please correct me if I am wrong.

In any case, your first command "velocity nanotube set 0.1 NULL NULL" is setting the velocity along x of all atoms in the nanotube to a constant value.
This will not allow the nanotube to deform along x. Moreover, the nanotube may want to roll instead of just sliding over the substrate.
These are all motivations why I am trying to implement the protocol I reported in my first email.

ok. i see. this sounds like a candidate for fix smd in constant
velocity mode. fix smd attaches a spring between the center of mass of
a group and a reference point and you can move that reference point at
constant velocity and implicitly get the force between the moving
object and its environment from the spring force, which you can
monitor.

axel.

Yes, that’s an option, of course the com velocity of the nanotube will not be constant in time using fix smd.
The point is that this fix introduces an external parameter via the spring constant, which needs to be tested/tuned for different systems, and I’d avoid that if possible.

I already did simulations implementing the fixed-com-velocity protocol in my own MD code, which is straightforward if you have access to the internal forces,

I wanted to check if there was a way to do it also in LAMMPS.

Thanks for the help,
Davide

Yes, that's an option, of course the com velocity of the nanotube will not be constant in time using fix smd.
The point is that this fix introduces an external parameter via the spring constant, which needs to be tested/tuned for different systems, and I'd avoid that if possible.

I already did simulations implementing the fixed-com-velocity protocol in my own MD code, which is straightforward if you have access to the internal forces,
I wanted to check if there was a way to do it also in LAMMPS.

since LAMMPS does velocity verlet timestepping and if you want to
constrain the velocity exactly, you would have to write a custom
integrator fix, so that you can enforce your constraint on both
half-updates of the velocity.

this is usually done in C++, but you could also try to use fix
python/move and then implement your protocol in python script code
invoked from there and access the internal force arrays from the
python wrapper for LAMMPS.

axel.

Thanks for the tip!
Davide

I’m not following these details, but this Q:

In general, is there a way in LAMMPS to access the total force due to the force-fields only, before external forces are applied?

makes me think of fix store/force which archives the forces
after pair/bond/etc forces are computed but before
external forces. If that fix is defined, any other part
of LAMMPS can potentially access them. E.g. to dump
the interatomic forces, w/out external forces added.

Steve

Thanks Steve.
I read the documentation on fix store/force, and it looks like it is what I need.
I will give it a try.
Davide

“fix store/force” solves the problem.
Thanks,
Davide