Hi everybody,
I want to fill a carbon nanotube with constant density of water molecules. So, could you please what is the simplest way to do this?
Thanks,
I would just create suitably sized bulk water system (take a small water box, equilibrate, and then replicate to size) and then load a data file with the nanotube on top of that with read_data add.
Then delete overlapping water molecules and equilibrate (first by using time integration only on the water with fix nve/limit and fix langevin to quickly dissipate high-potential energies near the interfaces and then fix nvt for a bit and finally fix npt so the system can shrink to the desired density (deleting overlapping atoms tends to delete a bit too aggressively but that cannot be avoided, so the fix npt run is needed to avoid little vacuum bubbles).
Thanks for your very fruitful suggestion. Hereby, I followed your comments and prepared input and data files as follows (data files also attached).
But unfortunately I get below error after step 0. I will be appreciate any comments guide me to improve it.
ERROR: Lost atoms: original 804 current 801 (src/thermo.cpp:439)
echo both
units real
boundary p p p
neigh_modify delay 2 every 1
timestep 1
atom_style full
bond_style harmonic
angle_style harmonic
read_data config.data extra/atom/types 1 extra/bond/types 1
bond_style harmonic
bond_coeff 1 469.0 1.0 #H-O
bond_coeff 2 450.0 1.418 #C-C
angle_style harmonic
angle_coeff 1 85.0 109.47 #H-O-H
set type 1 charge 0.4238
set type 2 charge -0.8476
pair_style lj/charmm/coul/long 8.0 10.0
pair_coeff 1 1 0.0 0.0 #H-H
pair_coeff 2 2 0.1553 3.166 #O-O
pair_coeff 3 3 0.07 3.9848 #C-C
pair_coeff 1 2 0.0 0.0 #H-O
pair_coeff 1 3 0.0 1.9924 #H-C
pair_coeff 2 3 0.104264 3.5754 #O-C
kspace_style pppm 1e-4
neighbor 2 bin
neigh_modify delay 0 every 1 check yes
group solvent type 1 2
fix rigid solvent shake 1.0e-6 500 0 b 1 a 1
fix NVT solvent nvt temp 298.0 298.0 100
mass 3 12.00
velocity solvent create 298.0 12345678 dist uniform
thermo 100
thermo_style custom step temp pe ke etotal
log logequi.txt
dump 1 all custom 10000 dump.lammpstrj id type x y z
run 100 # 50000
undump 1
##################
reset_timestep 0
read_data cc.data add merge group cnt
pair_coeff 3 3 0.07 3.9848 #C-C
pair_coeff 1 3 0.0 1.9924 #H-C
pair_coeff 2 3 0.104264 3.5754 #O-C
unfix rigid
region trash cylinder z 25 25 6.77 1 60 units box side out
group trash region trash
delete_atoms group trash bond yes
region cyl cylinder z 25 25 6.77 1 60 units box
group cyl region cyl
delete_atoms overlap 0.3 solvent cyl bond yes
#fix rigid2 solvent shake 1.0e-6 500 0 b 1 a 1
group H2O type 1 2
fix 1 H2O nve/limit 0.1
fix 2 H2O langevin 1.0 1.0 1000.0 699483
fix 3 H2O nvt temp 298.0 298.0 100
fix 4 H2O npt temp 298.0 298.0 100.0 iso 0.0 0.0 1000.0
thermo 1000
thermo_style custom step temp pe ke etotal
log logequi2.txt
dump 2 all custom 10000 dump2.lammpstrj id type x y z
run 5000 # 50000
config.data (1.6 MB)
cc.data (25.7 KB)
Lost atoms can have many reasons. This has been discussed many times before. In your specific kind of system the most likely explanation are close contacts. Please search through the LAMMPS forum/mailing list archives to see previous discussions and suggestions for remedies. The first step is always careful visualization around the timestep where you lose the atoms.
An alternative solution to the one proposed by Axel is to fill the nanotube will water using grand canonical Monte Carlo (fix gcmc), as done here:
The advantage is that you impose the chemical potential, the disadvantage is that it is slow because the probability of successful insertion is very small at large density.
1 Like
Thanks @akohlmey and @simongravelle for your comments.
I followed the procedure by @akohlmey and attached minimal
solvate.in (1.9 KB)
solvate.data (815.6 KB)
input and data files. But unfortunately I get some error
Non-numeric box dimensions - simulation unstable (src/KSPACE/pppm.cpp:1864)
and using below command to delete water molecules outside the CNT leads to delete small part of carbon (create hole in CNT). So, finally the number of water molecules decrease after 100 ps.
region trash cylinder z 20 20 6.9 1 41 units box side out #7.54
group trash region trash
delete_atoms overlap 0.5 trash cnt bond yes
In this way, I attached input and data LAMMPS file. I will be appreciate, if you help me to improve it.
Thanks in advance
I looked at your script. When you delete the water molecule, it seems that your cut-off is too small, I did it with 2 Angstrom :
delete_atoms overlap 2 solvent cnt mol yes
and the output is more reasonable.
Also, if you only want to delete the water outside the CNT, then why are you using the overlap keyword? Just delete the molecules that are not located within a certain cylinder:
delete_atoms group trash
Best,
Simon
Thanks for your cmments.
but if you use
delete_atoms group trash
you will delete some O or H atoms and hereby you get :
Bond atom missing in image check
And unfortunately reducing the radius of 6.9 for CNT makes a hole in it. So, if you visualize the dump file you will see O or H atoms which move freely.
region trash cylinder z 20 20 6.9 1 41 units box side out
Thanks for any suggestion
of course you have to use ‘mol yes’
This part does not make much sense to me. The radius of what ?
6.9 Angstrom is the radius of CNT. So, in order to delete outside water molecules, I have some problems:
If I get the radius of region trash similar to CNT radius and use the keyword ‘side out’, then
- with ‘bond yes’ → since the CNT undergoes some deformation, some carbon bonds will be omited. (plz look at the attached shot)
- with ‘mol yes’ → Since all of the carbon have the same molecule ID, so this command completely deletes the entire CNT.
- If I get the radius of
region trasha bit larger than CNT radius (To count the CNT deformation, this increased value should be at least 2 angstrom), so with this increased value some O and H atoms will remain there,
But if the atoms of the CNT are not part of the group ‘trash’, then you do not delete the CNT and everything is fine right ?
It seems yes. Although I froze C atoms, CNT was deformed here ( it can be seen from above shot) Do you think this happened because of being in contact the water molecules or not correctly freezing carbon atoms?
Thanks
It is not a matter what other people think. You are trying to do computational research and you have to know yourself what you did and that it was correct. Thus every time you get an unexpected behavior you cannot rely on others telling you what to do and then push forward again. Instead you need to step back and start over and re-build you input in steps and check at each step that the commands are doing what they are supposed to be doing and that the system still behaves the way you expect it to behave.
You have now started multiple discussion threads where you keep getting stuck in details and have failed to confirm that basic steps are correct and at the same time are rushing forward without much regard of warnings or issues that should alert you and with performing steps that are contradicting common sense or the documented behavior of LAMMPS. I see a pattern forming and is not looking good for you. If you do not manage to get your process of doing research under control and be more careful, vigilant and - most of all - more systematic, you will quickly reach the point where people are getting tired of dealing with your questions about details while you still have not gotten the basics right.
If that is the diameter of your CNT, then the region capturing the water atoms outside must obviously have a larger diameter.
Besides, when you delete atoms from a group, it is very straightforward to create a group that is the intersection of the atoms that form the water molecules and atoms that are in the selected region. By definition there cannot be any CNT atoms in that group and thus will not be deleted (unless you messed up your molecule IDs).
Both of these points are common sense.
I see that many interesting solutions have been proposed.
There is one that I would like to comment on, which is quite simple, and that is to use packmol and prepare the structures first. In it, you would enter the nanotube file (a .pdb for example) and follow another file with a single water molecule. In packmol you set up a very simple script so that the water molecules (with the number of your choice) are inserted into a cylindrical shape. There a preliminary calculation is already made so that the water molecules are not too close together.
Hope this helps.
Using packmol is inferior in this case because it tends to create too low densities (a side effect of avoiding close contacts) and the water structure will need significant equilibration time. It is certainly a good method for more complex geometries that are difficult to describe by other means.