[lammps-users] Implementing a Biasing Potential

Hi Steve and Lammps users,

I’m an MSci physics student at King’s College London in the UK, and I’m using Lammps to run simulations for my project.

I have implemented a harmonic biasing potential into Lammps; more specifically in the ‘fix_langevin’ routine. The purpose of it is to restrain two groups of atoms at a certain angle relative to each other, in order to run MD simulations for Umbrella Sampling. Once the biasing potential is found, the force components are calculated and added to the force array ( F[][] ) in ‘fix_langevin’.

My question is as follows: would I need to convert my force units to account for the fact that Lammps uses different internal units? Lammps seems to perform unit conversions in some routines, like using “ftm2v” in ‘fix_langevin’, but in others it does not. Furthermore, I’d appreciate your opinions on whether I’m approaching this in the right way. For example, is there a routine more suited to modifying the force?

My biasing potential is of the form V = k( cos( theta( x ) ) - cos( theta0 ) )^2, where cos( theta( x ) ) depends on the atomic positions, taken from Lammps’ memory.

Thanking you in advance,

Pepe Falahat

Hi Steve and Lammps users,

hi pepe,

I'm an MSci physics student at King's College London in the UK, and
I'm using Lammps to run simulations for my project.

excellent!

I have implemented a harmonic biasing potential into Lammps; more
specifically in the 'fix_langevin' routine. The purpose of it is to

ouch! that is a *bad* idea. you should implement a biasing potential
as a separate fix. have a look at fix spring for example.

restrain two groups of atoms at a certain angle relative to each
other, in order to run MD simulations for Umbrella Sampling. Once the
biasing potential is found, the force components are calculated and
added to the force array ( F[][] ) in 'fix_langevin'.

the way to do this is to compute and add the forces int a "post_force()"
method of a new fix. see the documentation about extending lammps and
particularly adding new fixes.

My question is as follows: would I need to convert my force units to
account for the fact that Lammps uses different internal units? Lammps
seems to perform unit conversions in some routines, like using “ftm2v”

that one is only required for the time integration. for computing
forces this is not really needed, since the force constant that
you have to provide as a user has to be in the proper units.
have a look at the "compute()" method in angle_harmonic.cpp

in 'fix_langevin', but in others it does not. Furthermore, I'd
appreciate your opinions on whether I'm approaching this in the right
way. For example, is there a routine more suited to modifying the
force?

see above. the design of LAMMPS is highly modular and allows
very flexible modifications of your system through the "fix"
infrastructure. during the processing of the input file, LAMMPS
collects lists of classes and executes the appropriate methods
during whenever suitable. this takes advantage of the abstraction
through derived classes and has lots of benefits for writing
add-on modules as for maintaining the code as a whole.

hope this helps,
    axel.

Hi Axel,

Thanks for your swift reply.

the way to do this is to compute and add the forces int a “post_force()”
method of a new fix. see the documentation about extending lammps and
particularly adding new fixes.

Right, I will give this a go. This is probably a stupid question, but out of interest, why can’t I simply add the forces in the “post_force()” method an existing fix command like ‘fix_langevin’?

that one is only required for the time integration. for computing
forces this is not really needed, since the force constant that
you have to provide as a user has to be in the proper units.
have a look at the “compute()” method in angle_harmonic.cpp

I see. That makes things clearer. By ‘proper units’ I assume you mean the units chosen by the user in the input file (e.g. “real”), not Lammps’ internal units?

Pepe

Hi Axel,

Thanks for your swift reply.

> the way to do this is to compute and add the forces int a
"post_force()"
> method of a new fix. see the documentation about extending lammps
and
> particularly adding new fixes.

Right, I will give this a go. This is probably a stupid question, but
out of interest, why can't I simply add the forces in the
"post_force()" method an existing fix command like 'fix_langevin'?

in principle, you can do whatever you want, but it
is *very* bad program design. not only within LAMMPS
but overall. you want to mix up separate functionality.
that makes programs unmaintainable.

the whole point of the overall design of LAMMPS is to make
it easy to add new features in a modular way and without
having to worry too much about the rest of the code.
this way all additions have the maximum benefit for everybody.

just assume you would like to use your angle restraint
potential with a different thermostat.

> that one is only required for the time integration. for computing
> forces this is not really needed, since the force constant that
> you have to provide as a user has to be in the proper units.
> have a look at the "compute()" method in angle_harmonic.cpp

I see. That makes things clearer. By 'proper units' I assume you mean
the units chosen by the user in the input file (e.g. "real"), not
Lammps' internal units?

"external" units are the "internal" units. what changes
are certain constants that are needed internally, e.g.
to compute the potential or kinetic energy.

cheers,
   axel.