Initial condition:
a) I have a solid surface A with elements a, b, c, d, e interfacing with randomly mixed liquid of B (f, g, h elements) and C (i, j, k elements), both are prepared with proper density.
b) X & Y direction are periodic while Z-direction is non-periodic, but both edges of Z-dir have reflective wall with reflect/stochastic wall.
c) The bottom 3 angstroms of A are frozen (not applying any thermostat) to anchor the A in the place.
d) Total ~11k atoms, and total z-dir length is about ~6nm, where the solid and liquid take half, so 3nm thickness is A and 3nm thickness takes liquid B and C.
MD simulation:
a) I run NVT MD at high temperature close to MP of A.
b) During MD simulation, A and B, C reacts with each other, creating reaction products and the reaction consumes the amount reactants, B and C.
Problems:
a) The first problem is that I see a lot of reaction products accumulating on the surface, prohibiting further reactions. I think this is because of reflective wall and short z-dir length, products have nowhere to go, reflected from the wall and went down and react back to the surface. Maybe this is something really happens in experiment, but I doubt that based on experiment results.
b) The second problem is that the MD looks like interdiffusion MD, not the surface reaction MD. I get the density profile of B and C, and all molecule numbers (and density too) decline as MD proceed because they are consumed by interfacial reaction. And this will impact the rate of reaction.
Here are my questions:
Easy solution is making liquid phase (B and C) a lot longer in z-dir with a lot more molecules. But due to limit on computational cost, we can’t do that. Then, what would be the better options to solve two problems?
For problem 3a):
I want to delete product molecules whenever they are created in the system after checking every iteration (or every 5 or N iterations)
I know LAMMPS has delete atoms command, but to use that, I need to know which atoms will become product at precise point, but how do I know that before running MD? So I can’t use that.
I was thinking to check and delete reaction product molecules (any molecules that contains any elements among a~e that is not a part of A). I can do that as post processing using depth first search using fortran, but I don’t know how to do this in-situ, in real-time during MD from LAMMPS.
Can I use group with dynamics flag? But still, I don’t see how to selectively find the reaction product molecules for every iteration… I can use DFS algorithm for dump file for post-processing, but I have no clue how to do that with dynamic group…
For problem 3b)
I want to supply new reactant molecules whenever the number of B or C declines after checking every iteration (or every 5 or N iterations)
I know LAMMPS has fix deposit command, and I saw people use that for gas phase. But now my system is liquid at high T, and there is no empty space to put new reactant molecules. So, I can’t use that. Also, counting B and C as postprocessing is easy, but how can I implement that in the LAMMPS?
I was thinking to fix deposit above the reflect/stochastic wall and push them into the MD system, but I don’t even know if that is possible and if that is good idea or not. Even if I do that, can I do that (create and push) for every iteration for every new molecule? So, at this point, I really have no clue…
Can please anyone share some tips, idea, experiences, and scripts for this problem? I think there should be a lot of people who ran interfacial MD between solid and liquid, so some people should’ve thought about this…
It is not at all uncommon to have chemical reactions that are diffusion limited. Classical MD by itself does not offer access to the time and length scales necessary to handle such a situation in any realistic way. This is outside my area of expertise but I would suggest that you do some more research in the published literature. All I can point to are two hints of potential pathways. a) a long time ago I was asked to be on a team to submit a proposal to implement a multi-scale modeling environment that would combine kinetic Monte Carlo simulations with classical MD. kMC is one way of modeling chemical reactions at interfaces, but you need to define the reactions and their reactions constants etc. and b) the DPD-REACT package may offer some framework that may be possible to adapt for your needs. That would require some significant degree of coarse graining as one component to overcome the size and time limits.
Thanks, and yeah, I know KMC but you already pointed out it requires proper reaction constants. And getting reaction constants for all possible interactions of that many element types seems not trivial. Also, I think KMC makes sense if the problem is just for dissolution, but I’m not sure if KMC is a good idea for the description/prediction of reactions and several reaction products (especially predict some potentially possible and previously unknown products) with those reactions. But I might be wrong.
Let me check about that package, thanks.
===
Meanwhile, could LAMMPS Python module help to resolve two issues; removing recording and removing product molecules, and add new reactant molecules to liquid with high T, during mD? I never used Python module, so I don’t know if I can control things during MD with Python module.
If everything isn’t working, then I will look into other method you mentioned, or look for computational resources to make B&C liquid much larger from this MD setting…
You have titled your topic that you want to simulate your system more realistically. What you are asking about is the opposite. There has to be a scientific justification for the removal and addition of molecules and a “protocol” that allows those processes to be connected with a bulk system behavior. Otherwise, you are just doing a computer animation with no scientific value; just to have something that looks like you want it to look and not what is realistic. How do you know that your intuition is correct? How can you justify and quantify your results if you make what are effectively arbitrary changes to the system?
As I already mentioned, at the core of your question is both a length-scale and a time-scale problem. Your statement, that you do not want to spend more computational effort is in my opinion a very unscientific one, for two reasons: 1) there are plenty of computing resources available (for example in regional or national supercomputing centers), it is just a matter of investing the effort to get access. Your choice of simulation settings should be determined by the science and not convenience. 2) if you want to accelerate time-scales or bridge length-scales, you need a solid scientific footing for that. As I already stated MD as such is not capable of that. You have not mentioned what are the time and length scales of your experimental data and how this compares to your simulations? Not to mention concentrations of components.
This is why I suggested to do your homework and dig into the published literature to learn what others have done and what quality of results they have obtained. Ultimately, you have to decide whether you are chasing after an Oscar or a Nobel Prize.
Thanks for reply. I will take your reply as “no” on my questions regarding Python module.
Regarding your concerns, sorry, I cannot disclose details, but we have pretty good reasons for why checking from all-atom scale is necessary and why we can’t use public computing resources. This is already a part of multiscale effort, and we believe science part is well justified.
About MD, the temperature is high near melting point of A, so the reaction is fast, that is good thing for us. But short z-dir thickness and limited number of reactant molecules are likely to be the problem based on our current observations. We just wanted to check and make sure about these points and that is why I’m interested if those two problems can be resolved from MD or not. Hope this answers some of your concerns, and I appreciate your answer. I will check DPD module and we will see what we can do. Thanks.
That is not correct. I simply chose to not discuss something that is in my opinion distracting for the problem at hand.
That is your choice of course, but you are effectively telling me that you don’t care about getting advice about making your model realistic like your topic says.
But those problems are self-imposed and you cannot convince me by just restating that you believe in your approach being correct that it actually is correct and and realistic.
Not really, all I am reading from your responses is that you want a) somebody that confirms your choices (I won’t) and b) some suggestions to implement your strategy (I won’t either because, as I stated before, this will turn your computer simulation into a computer animation).
Not at all. It just tells me, that I have wasted my time trying to raise scientifically motivated concerns. But if you neither willing to convince me nor willing to follow the advice given (DPD-REACT has the same limitations as kMC so I don’t think it will help you when you already dismissed kMC), there is nothing left to discuss.