Skip to content

Moving Brushes And Janus Particles

Jump to bottom
jspngler edited this page Aug 8, 2024 · 1 revision

Overview

The source file "generate/movingBrushesAndJanus.cpp" will create an mpd file with the following requirements:

  1. Mobile polymer brush composed of variable length polymers.
  2. Many tessellated nanoparticles with variable diameters.

As with any other generate file, the input parameters are the adhesion/repulsion of the brush monomers and nanoparticles, the density of the polymers in the brush, and the number of nanoparticles placed randomly above the brush. Compilation, requiring at least c++11, can be done by:

module load gcc/8.2.0 #Don't need this on bigblue.memphis.edu
g++ -O3 generate/movingBrushesAndJanus.cpp -o ~/bin/movingBrushesAndJanus

The input arguments are:

movingBrushesAndJanus name seed nBrushes arealDensity UminBrushNano UmaxBrushNano nNanoparticles nanoRadius UminBrushBrush UmaxBrushBrush nMonomers bondLength temperature

The test input, generating test.mpd and test.xyz, that was used during testing (as a starting point) is:

movingBrushesAndJanus test 448588 4000 0.2 -2.0 200 25 6 0 100 50 0.7 3.0

The above should create a 25 nanoparticle, 4000 polymer brush system with an average of 50 monomers per polymer (using poisson distribution, roughly related to Flory-Shulz distribution) and nanoparticles with an average radius of 6 r_min (linear relation with the number of nanoparticles, none are the same size, this part needs to change). The nanoparticles have a strong adhesion (U_min=-2.0) to the polymer monomers. The polymers are repulsive to one another and the core of the tessellated nanoparticles.

I recommend opening the resulting test.mpd file and modifying the "finalTime" variable to something like 100000 or more tau (timesteps*deltaT).

MD.cpp can be compiled like normal module load intel/2022.2; icc -fast -fopenmp MD.cpp -o ~/bin/MD,(on hpclogin.memphis.edu) or module load aocc/4.1.0; clang++ -Ofast -fopenmp MD.cpp -o ~/bin/MD (on bigblue.memphis.edu), and a submission script on the University of Memphis HPC might look like the following (we'll name it test.sh):

#!/bin/bash
#SBATCH --cpus-per-task=40
#SBATCH --partition=computeq
#SBATCH --mem=2G
#SBATCH --time=5-00:00:00

#Note: --partition should be "acomputeq" on bigblue.memphis.edu 

MD test

Then, submit with sbatch test.sh. The simulation above would run on 40 CPU-cores in the computeq partition (alternatively: shortq and wholeq) with 2 GB of memory allocated for up to 5 days. If the simulation terminates before then, just resubmit the above. The simulation doesn't need all 40 CPU-cores (between 8 and 40 is good, more==faster. Also, although 192 cores are available on bigblue.memphis.edu, 64 is maximum recommended) or all 2 GB of memory (memory allocation is pretty efficient).

Check out /home/jspngler/users/chdavis for the testing done with this particular code. Running MD on this code produces:

  1. bend_test.dat: average bending angles in cos(theta) values.
  2. frames_test.xyz: frame file for positions of all particles in the system (this can be rather large).
  3. kinetic_test.dat: total kinetic energy for all particles in the simulation over time (aka sum((3/2)_k_b_T)).
  4. meanSquareDisplacement.dat: displacement of every "molecule" type in the simulation over time.
  5. size_test.dat: simulation box size over time.
  6. flicker_test.dat: ignore this one.
  7. kEnergyDensity_test.dat: total kinetic energy probability density of simulation.
  8. lBond_test.dat: average bond length of every bond in system over time.
  9. potential_test.dat: total potential energy for all particles in simulation over time.
  10. temp_test.dat: thermostat set point over time.

Modifiable parameters in the code

Two things can be modified from this file, the inputs for various functions that modify state (something like function(system,...);) and the sequence that the system state is modified.

Inputs

variableBrush

template <typename T>
void variableBrush(Blob<T> &System, std::vector<int> polyLengths, threeVector<T> pos, T bondLength, T arealDensity, std::vector<T> constants, T aspectRatio, int bottomType, int chainType, T k, T temperature)
  • T - template parameter for floating point precision (double or float).
  • Blob &System - Reference to the system state.
  • std::vector polyLengths - Polymer lengths by monomer count. Also, number of polymers from size() member. In other functions under models/, this might be split up into a pointer, requiring input.data(), and the count, requiring input.size().
  • threeVector pos - Position offset of bottom corner, mostly for the z component.
  • T bondLength - Harmonic bond length of polymer monomers. Also describes height with polyLengths.
  • T arealDensity - X and Y plane density of the polymer. With number of polymers, it describes the area of the polymer brush.
  • std::vector constants - Harmonic bond (and optional bending) components and must contain 2 or 4 elements.
  • T aspectRatio - Anisotropy of the system, the x/y ratio for system.
  • int bottomType - Type of monomer at substrate (at bottom, AKA pos).
  • int chainType - Type of monomer that makes up the bulk of the polymer.
  • T k - Constant for metropolis algorithm that relaxes the substrate (spacing between polymers).
  • T temperature - Temperature of simulation. Used to initialize the velocities.

distributeNanoparticles

template <typename T>
void distributeNanoparticles(Blob<T> &System, std::vector<T> radii, std::vector<int> nTessellations, threeVector<T> lowerBound, threeVector<T> upperBound)
  • T - template parameter for floating point precision (double or float).
  • Blob &System - Reference to the system state.
  • std::vector radii - List of nanoparticle radii, and the number of nanoparticles. Should be 2 or more r_min units.
  • std::vector nTessellations - List of nanoparticle tessellation counts. Usually 0 to 5. Each increase will quadruple the number of nanoparticle constituent particles.
  • threeVector lowerBound - Lower boundary where nanoparticles are placed (up to the shell surface, ~radius).
  • threeVector upperBound - Upper boundary where nanoparticles are placed (up to the shell surface, ~radius).

Sequences from here: https://gitlab.com/ericspan/mpd-md/-/blob/32ba04955927b3a96f29e8a0f44aa22e668f7861/generate/movingBrushesAndJanus.cpp

main()

The main function takes a lot of arguments (see above). The molecules, particles, and boundaries can be adjusted around to any position in the main function after the system state has been initialized on lines 301-334.

From main, line 281, the general sequence of events is:

  1. Read in the parameters from the command line arguments (281-300).
  2. Initialize the system state (301-334).
  3. Initialize the polymer inputs (polymer lengths, constants, positions, etc...) (335-359).
  4. Run variableBrush function to insert the polymer brush into the system state (360).
  5. Initialize the nanoparticle inputs (boundaries and radii) (364-406)
  6. Run distributeNanoparticles function to insert nanoparticles into the system state (408).
  7. Create a FLOATING_BASE molecule (implements solid substrate as described by include/potentials/laradjiPoursoroush.h) and adding all the particles that will interact with that substrate (411-444).
  8. Add a harder boundary that prevents molecules from entering on either side near the z=0 plane (446-461).
  9. Initialize and add non-bonded constants (implementation described by include/potentials/laradjiRevalee.h) (469-509).
  10. Store the simulation (513-515).
  11. Store initial frame for vmd (517-530).

variableBrush()

This particular function will generate variable length polymers on a grid, then relax the configuration using a Metopolis algorithm. The main things to modify are the polymer densities and polymer lengths. The only downside to this particular configuration is that the substrate, from main() on lines 411-444, isn't particularly strong (from -1 to -6 Umin) when the polymer density (from ~0 to 0.5 polymers per rmin) is high. Increasing the strength, by lowering Umin, can cause the brush to shed more.

From the definition, on line 13, the sequence to place brush polymers is:

  1. Error checking the constants (17-21).
  2. Count the number of particles, bonds, and bends in polymers (24-34).
  3. Initialize the bonding and bending molecule components (37-53).
  4. Adjust the box size to place all the polymers (56-77).
  5. Set the rms velocity for the monomers using the system temperature (80-83).
  6. Pre-allocate the particles in the system (86).
  7. Set up the base lattice for polymers (89-117).
  8. Relax the lattice using Metropolis algorithm (119-121).
  9. Place the polymer monomers, bond indices, and bend indices into the system (124-171).

distributeNanoparticles()

The nanoparticles added to the simulation are disbursed into a box with defined boundaries. If the boundary is too small, the function will throw an error saying "10,000 nanoparticle placement attempts exceeded". You should probably make a best guess on the box dimensions based on how many nanoparticles of an average size should fit inside this box with at least 1/2 or more of the volume being empty. Another note, you can take the NANOCORE particle out of the nanoparticles, but the longer polymers will ingress into the interior of the janus shell.

The janus shell itself is based on a icosahedron with tessellation increasing the mesh density. You should probably aim for 1 to 10 particles per unit area. Also, note that each particle added to the shell increases the repulsive/attractiveness of the nanoparticle to nearby nanoparticles/monomers. A better modification would be to have a different type for each nanoparticle (assuming a homogeneous shell) and scale the strength based on the diameter. The interior nanoparticle is already density scaled (to ~5.88 particles per rmin^2 on line 253).

From the definition, on line 175, the sequence to place janus nanoparticles is:

  1. Set up the dimensions of the box in which nanoparticles will be placed (178).
  2. Set up the vector of nanoparticle coordinates (182).
  3. Set up the random number generator (185).
  4. Loop through all the particle radii, placing the nanoparticles randomly while looking at previous nanoparticle placements to determine if there is any overlap. There is a hard coded 10,000 attempts at line 221. (187-230)
  5. Initialize NANOCORE type to prevent polymers from entering janus shell (232-236).
  6. Loop through all the radii again to build the janus shells, using a function that adds them to the system in a similar manner as this function, and incorporate the NANOCORE to the anchoring particle (238-274).
  7. Add the NANOCORE molecule to the system (277).

"Hardcoded" values in the code

Feel free to change these to anything you wish. I've only chosen them as a solid starting point. I've run about ~45 simulations to narrow the ranges a little bit (in addition to the work Yu Zhu and Asama Poursoroush did over hundreds of simulations).

  1. BB, BM and NANOTYPE are just numerical representations of the Brush Base (substrate), Brush Monomer, and nanoparticle shell particle types. They can be removed and replaced with anything.
  2. Max Z height offset (1.0) for the brush base beads on line 58.
  3. Relax substrate has 1000 iterations on line 121.
  4. 10000 maximum attempts at placing a single nanoparticle on line 221.
  5. Yu Zhu's Janus function, called on line 246, takes 1200 for kBond between shell beads (aka edges), 45 kBond for bonds between the center anchor and vertex beads, 250 kBend and 50 particleBond isn't used (oops), and 0 is the janus particle's ratio from one type to the other (so unused, but can be modified on a per particle basis).
  6. Eric and Laradji's constants for an impermeable shell particle, called on line 253, takes a CUTOFF of 2.0, RMIN of 1.0, radius 1.3 rmin smaller than the janus shell, 5.88 particles per unit area scaling, 30.0 Umax, and 0.0 Umin (making a repulsive interaction to all particles).
  7. 1.0 is the damping for the langevin thermostat on line 304.
  8. 4 types are defined on line 305 for 0=immobile, 1=BM, 2=NANOTYPE, and 3=BB particles. Make this as big as you want, but when you exceed 100 types, viewing the frames file in vmd is a little difficult.
  9. 2.0 for cutoff on line 314. This is from a simulation that Laradji did with Joel Revalee in ~2006, I believe.
  10. Initial time of 0 on line 322.
  11. Final time of 10000 on line 323.
  12. DeltaT of 0.02 on line 324 (same paper as cutoff above). This can be lowered to increase accuracy and stability, but don't expect the simulation to run very fast.
  13. Store interval of 1000 (1000/0.02=50000 time steps) per frame/restart interval on line 325. We've used this between 100 and 5000. It doesn't affect the simulation itself, but may affect parameters you measure after the simulation completes. Also, low values can take a significant amount of storage.
  14. Measure interval of 100 (100/0.02=5000 time steps) per measurement interval for potential and kinetic energy on line 326. I recommend modifying the MD.cpp code to take new measurements under line ~514, in the main loop of the simulation under comment "Measurements are output here".
  15. Offset of 3.0 on line 344 is for the poisson distribution. Basically the average polymer monomers minus 3 reserved for minimum lengths.
  16. 100 kBond for polymers is pretty ridged for the harmonic bond between monomers on line 351.
  17. On line 361-362 variableBrush has 1.0 for x/y plane aspect ratio, 100.0 for k in metropolis algorithm, and 3.0 for temperature (I probably should have used System.readTemperature() function!).
  18. 0.2 scaling and 0.5 mean offset for the nanoparticles on line 391. Read the comment in this section of code as it is intended to be modified.
  19. Umin and Umax for repulsive and attractive components for particles near the floating boundary on lines 430 to 436. This might need tweaking to get the polymers to remain adhered.
  20. 5.88 for floating boundary surface density on line 438. Scaling constant for boundary potential.
  21. Position of hard boundary along the z component on line 456 could be adjusted within System.readSize().z and 0.
  22. 200 for repulsiveness of hard boundary on line 458. Higher makes the boundary harder. And don't reduce it to zero or lower (divide by zero error when particle lands on boundary).
  23. Umin=0 and Umax=100 for non-bonded forces are repulsive by default on line 484-485.
Clone this wiki locally