How to define a force toward a specific point

Hello everybody, I am having trouble in finding a way to define a force that would drag all my atoms towards a point with specific coordinates. This is my setting: lj style, 3D, a box size 50. I defined a region and I would like that each atom passing through said region feel a force dragging it toward the center of the region. Now, the easiest way I found is the command “fix ID group-ID drag x y z fmag delta”, so I tried it, but I found two invenkiences (one understandable, the other one less so): 1) it doesn’t work with a specific region (or at least I haven’t found a way around it), and 2) when I set the xyz coordinates… the point fall short of the real coordinates. It’s hard to explain, but from the OVITO representation, whenever I input “25 25 0” to point towards the center of the bottom face, the particles point towards something close, but not quite the exact point.
After the drag command, I checked “aveforce”, but then I am having difficult to set dynamic coordinates in order to define a vector that would point exactly towards the intended point: the “variable” command eludes me, unfortunately.

Is there any example/documentation outside the “official” one you could point me at? Thank you in advance!

The way to do this would be to use fix addforce with atom style variables for the force in x, y, and z.
For variables there is an rmask() function which is 1 or 0 when an atom is in the region or not. So you just compute the force based on the distance from the target and then multiply that value with rmask() so the added force is zero for all atoms outside the region.

No. There is no documentation beyond that. All available features are documented in the LAMMPS manual.

fix aveforce is not doing what you described, fix addforce is the ticket.

Why not have a fixed atom at the exact location which has an attractive interaction with all other atoms?
Although you will need a large cutoff and though increase the computation time!

Ok, thanks! I asked about more documentation because, to be honest, I find the LAMMPS documentation a bit hard to understand, especially the examples. As for the rmask(), for instance, I tried some syntax but I am not able to make it work, so I am trying a rerouting like this:

variable tx equal 25.0
variable ty equal 25.0
variable tz equal 0.0

variable dx atom x - v_tx
variable dy atom y - v_ty
variable dz atom z - v_tz

variable dist atom sqrt(v_dxv_dx + v_dyv_dy + v_dz*v_dz)

variable fx atom 100.0*(v_tx - x)/v_dist
variable fy atom 100.0*(v_ty - y)/v_dist
variable fz atom 100.0*(v_tz - z)/v_dist

fix 8 all addforce v_fx v_fy v_fz region skin

However, I always get the error “ERROR: Illegal variable command: expected 3 arguments but found 5 (src/variable.cpp:462) / Last command: variable dx atom x - v_tx” , but I don’t understand which ones are the wrong arguments.

Two things:
a) you are not quoting your input correctly and thus it is not typeset in a readable way. Please see the guidelines and suggestions post for details
b) you are not paying attention to the documentation. The documentation for the variable command clearly says that there either must not be any whitespace in the expressions or you need to use quotes.

As for difficulties making sense of the documentation. You cannot expect commented examples that explain everything. Do you have an idea how much effort this is and how few people are actually working on maintaining LAMMPS and how much effort it is as it is? The LAMMPS documentation is mostly a command reference. So to understand better you have to apply “the scientific method™”, i.e. create minimal inputs, make a hypothesis of what should the result be, check the result and compare, then either adjust your hypothesis after re-reading the documentation, or move on to the next step.

As a point of reference, I am providing a simple test example, to provide a proof of concept. This input is based on the “melt” example from the LAMMPS distribution and defines cubic region in the center where atoms are pulled toward the center of that region with fix addforce. It runs for a few steps and then outputs the various intermediate per-atom properties as well as the total force.


units           lj
atom_style      atomic

lattice         fcc 0.8442
region          box block 0 10 0 10 0 10
create_box      1 box
create_atoms    1 box
mass            1 1.0

# compute center
variable cx equal 0.5*lz
variable cy equal 0.5*lz
variable cz equal 0.5*lz

# compute distance from center
variable dx atom x-v_cx
variable dy atom y-v_cy
variable dz atom z-v_cz
variable dist atom sqrt(v_dx*v_dx+v_dy*v_dy+v_dz*v_dz)

# region in which to apply force
region force block 2 8 2 8 2 8

# total force and projections on x, y, and z
variable force atom -rmask(force)*2.0*v_dist
variable fx atom v_force*v_dx/v_dist
variable fy atom v_force*v_dy/v_dist
variable fz atom v_force*v_dz/v_dist

velocity        all create 3.0 87287 loop geom

pair_style      lj/cut 2.5
pair_coeff      1 1 1.0 1.0 2.5

neighbor        0.3 bin
neigh_modify    every 10 delay 0 check no

fix             1 all nve
fix             2 all addforce v_fx v_fy v_fz

run 10 post no

write_dump all custom dump.monitor id x y z v_dx v_dy v_dz v_dist v_force v_fx v_fy v_fz fx fy fz modify format int "% 4d" format float %10.5f


I apologize for the formatting, sorry. As for my outlook on LAMMPS, I shall confess that it mostly comes from the contempt toward my supervisor: after ten months spent programming in Python, and while my code was reaching the stage of “almost done” - meaning that I was studying how to get a phase diagram out of it - he told me to stop using python and switch to this program, that, honestly, I never heard of (but yeah, I did not even look for it, since I didn’t need to). Of course, now I have only a couple months to learn everything while doing stuff, and well, I am not very fond of it (and the whole situation). Sorry, I did not mean to leash out at you (or LAMMPS or the documentation), I was just sharing the origin story of my mood :stuck_out_tongue:
I will review your kind suggestion (and the mailing list documentation), thank you so much!

Cheers,
Ben

Argh, more sorry, I was writing my reply, and I didn’t notice your kind answer :frowning: thanks for the example, I will study it immediately!

If you want to have some fun, there is a way to run LAMMPS from Python: 2.3. Run LAMMPS from Python — LAMMPS documentation

I like that the introduction ends with

[…] give some examples how to use it.

:sweat_smile: