I would like to create a molecular dynamic simulation of interacting magnetic nanoparticles moving through a barrier.
(The possibly unusual part is that there is a coupling between the particles, the magnetization of one affects magnetization of another, item 3## below.)
The features of the simulation I need are:
- I would track N magnetic nanoparticles, each particle a rigid spherical body (e.g. 100 nm in size). So the state vector would be R = [R1;R2;…;RN] where Rj=(xj,yj,zj,wx,wy,wz), the 3D location and 3D orientation for the jth magnetic paricle.
(see the comment below regarding “fix rigid” and DATA files.)
- Each particle contains Kj iron oxide cores inside it. An iron oxide core is usually about 5 nm in diameter. So each rigid sphere contains a random placement of iron cores inside it, randomly placed. As the sphere translates and turns, the cores are displaced just by rigid body translation and rotation of that nanoparticle sphere.
There are many different versions of “fix rigid … molecule” which should do what you want.
Create a data file and assign each rigid sphere a different molecule-ID number
https://lammps.sandia.gov/doc/fix_rigid.html
As for how to create DATA files, you can either write a program to build them yourself,
https://lammps.sandia.gov/doc/2001/data_format.html
https://lammps.sandia.gov/doc/read_data.html
https://lammps.sandia.gov/doc/atom_style.html
…There is also a list of DATA-file builder tools at
https://lammps.sandia.gov/prepost.html
Either moltemplate or topotools should be general enough to build these kinds of particles. (Other molecule builders tend to be more specialized.) Topotools uses TCL. Moltemplate uses its own “object-oriented” style language.
- There are magnetic forces on each on iron oxide core. Two types of magnetic forces:
a) an force on each core due to the applied external field;
First, I found the following PDF file which looks vaguely relevant, but what they did is not as sophisticated as what you are trying to do:
http://www.rsc.org/suppdata/sm/c0/c0sm01424a/c0s
The closest thing that LAMMPS currently has to implementing an external B field is “fix addforce”, but I don’t think you can use it to apply torque to a point dipole. However if you represent each dipole by constructing it from two discrete particles, then it should be possible. Take a look at these (related) web pages:
https://lammps.sandia.gov/doc/fix_addforce.html
https://lammps.sandia.gov/doc/variable.html
If that is not flexible enough for your needs, it sounds like you will have to write your own “fix”. In LAMMPS, a “fix” refers to some code which modifies the equations of motions of the atoms.
Incidentally, you might not have to do this. Tim Sirk posted that he had written a “fix bfield” a couple years ago. He still checks the mailing list, but your more likely to get a reply from him if you mail him directly.
https://lammps.sandia.gov/threads/msg66010.html
A (slopily written) introduction to writing fix was posted a few days ago
https://sourceforge.net/p/lammps/mailman/message/36443505/
…and you should also read the official lammps documentation for customizing the code here:
https://lammps.sandia.gov/doc/Developer.pdf
b) dipole-dipole magnetic interactions between the cores.
Fixed point dipole-dipole interactions are already implemented in lammps:
https://lammps.sandia.gov/doc/pair_dipole.html
If you want to use point dipoles, then this will require that you use “atom_style hybrid full dipole” or something similar:
https://lammps.sandia.gov/doc/atom_style.html
…however, if you need polarizable dipoles, then you cannot (currently) use point dipoles, so please ignore what I just said (and the links above).
(##) Part b creates a nonlinear time-varying magnetic coupling between the cores.
If core p is magnetized, that creates a surrounding magnetic field that can change
the magnetization of cores q,r,s. The end result is that at each moment in time
there is a linear matrix equation ( [ A® ] M = B) that needs to be solved to find
the magnetizations of all the core (M = [m1x,m1y,m1z, m2x,m2y,m2z,…] ).
Once the magnetizations are known, then the force on each core can be
calculated. Once force are known, rotation is just sum of the torques on
each sphere.
An approach similar to the one used for polarizable electric dipoles could be used. (See below.) But I suspect this will need to be modified in a non-trivial way if the matrix is not diagonal. So it sounds like you are going to have to write your own fix to implement this feature.
I suggest you copy and modify the code for “fix drude”. You can obtain this code by downloading the lammps code, cd to the “src/” subdirectory, and typing “make yes-user-drude”. This will copy the “fix_drude.cpp” and “fix_drude.h” files into the “src/” directory where you can copy and edit them. Your fix will probably be more complicated. Anyway, here is the documentation for “fix drude”
https://lammps.sandia.gov/doc/Howto_polarizable.html
https://lammps.sandia.gov/doc/Howto_drude2.html
- The applied magnetic field and magnetic gradient is prescribed and is simple, is just set by an analytical formula. I don’t need to solve Maxwell’s equations.
Good.
- The magnetic particles will have diffusion and steric hindrance (so they don’t overlap each other). It’s ok to make the steric hindrance smoothed out over a long length scale, e.g. over 20% of the particle diameter, so that the governing equations are less stiff. The real important part is in the core to core interactions and how the forces that creates causes the magnetic particle to aggregate and chain up. If they aggregate with 20 nm gaps instead of 1 nm gaps between them, that’s ok (and I can correct for this by choosing particle size and core loading density appropriately).
You can use Lennard-Jones particles with large “sigma” parameters to implement as much steric hindrance as you need.
If you need to have the Lennard-Jones barrier centered exactly half-way between the two particles in the dipole (assuming you are approximating dipoles as two particles), then you can add a third particle between the other two and connect all three together with (stiff?) springs.
If you decide you don’t need polarizability and are using point-dipoles, then you can use “pair_style hybrid/overlay dipole lj/cut …” (or something similar):
https://lammps.sandia.gov/doc/pair_hybrid.html
- There will be viscous drag on the particles, both translational and rotational drag. It’s ok to keep this part simple (there are modifications to Stokes drag when particles are nearby, but it is not the central feature I am trying to simulate so is ok to keep this part as simple as possible.)
There are many ways to add viscous damping and thermostatting to LAMMPS. You can even add hydrodynamic interactions.
- If helpful, it is ok to ignore inertia (both translational and rotational inertia). These are nanoparticles with very tiny mass. So m x’’ + mu x’ = Fmagnetic + Fsteric + dW could be reduced to x’ = (Fmagnetic + Fsteric + dW)/mu where m is the mass, mu is the viscosity coefficient, ’ denote d/dt, and dW is Brownian noise. Whether LAMMPS allows a drop from second to first order dynamics, I don’t know, but it’s ok from a physics point of view.
For what it’s worth, in my own experience, I was able to get simulations to move much more quickly and efficiently using second order dynamics (with fix langevin) than writing my own code that implemented Brownian dynamics. I ranted about this topic here:
https://lammps.sandia.gov/threads/msg69500.html
For the timescales you probably care about, it does not matter if you use fix langevin.
Docs below:
https://lammps.sandia.gov/doc/fix_langevin.html
https://lammps.sandia.gov/doc/fix_nve.html
The following web sites are probably not relevant, but I mention them anyway:
https://lammps.sandia.gov/doc/pair_dpd.html
https://lammps.sandia.gov/workshops/Aug13/Denniston/Denniston_LAMMPS2013.pdf
https://lammps.sandia.gov/doc/fix_lb_fluid.html
https://lammps.sandia.gov/doc/PDF/SPH_LAMMPS_userguide.pdf
http://web.math.ucsb.edu/~atzberg/pmwiki_intranet/index.php?n=SoftwareMango.Homepage?setskin=mango_selm_Layout1
https://www.research-collection.ethz.ch/bitstream/handle/20.500.11850/153451/eth-5603-01.pdf
- There will be a barrier that the particles have to cross. The key issue is to better understand the physics of that crossing. The barrier could be a) a region of increased viscosity, or b) pores (meaning steric hindrance at a wall except for gaps in that wall). I need to see how the particles cross that barrier as the magnetic field moves and aggregates them.
You could either use “fix addforce” (as I mentioned above)
https://lammps.sandia.gov/doc/fix_addforce.html
https://lammps.sandia.gov/doc/variable.html
or create a wall by placing atoms at locations at specific locations in the barrier, and preventing them from moving. For an example of this, see:
http://moltemplate.org/visual_examples.html#translocation
http://moltemplate.org/visual_examples.html#nanotube+water
(If you don’t use “moltemplate” then you can ignore the “.lt” files and just look at the “run” files in these examples. The “run.in.nve” and “run.in.npt” files for both examples have the instructions for specifying groups of atoms which you wish to remain immobile.)
If at all possible, please let me know if this is something that is sensible/feasible to set up in LAMMPS. I don’t know how easy/hard it is to customize forces in LAMMPS. I would have to set up the kind of rigid-body multi-core particles I am interested in. I am especially concerned about how to put in the core-to-core magnetic coupling mentioned in item 3## above because that involves solving an additional (linear at each moment) system of equations at each time step within LAMMPS.
It is my first time posting to this ListServe. Any thoughts would be much appreciated.
My email is BenjaminShapiro.Shapiro@…24…
Thank you
Ben
I don’t know if this helps.
What you are doing is probably possible in LAMMPS if you are willing to write some C++ code (and/or email Tim Sirk, see link above). In that case, it’s probably easier to start with LAMMPS than writing your own MD engine from scratch. (Gone are those days, hopefully forever.)
Good luck.
Andrew