Surface tension result, correct number but negative, Godin

GROMACS version: 2021.4-Ubuntu-2021.4-2
GROMACS modification: No

Hi everybody,
I have tried to set up a simulation with which I obtain the same results regarding surface tension through thiols on a gold surface as in the publication by Godin from 2004 (https://doi.org/10.1021/la030257l).
After different variations of input parameters I got the nearly same result for the surface tension of 0.51 N/m (mine is 0.52 N/m), but my result is negative.
I noticed, that there already is a topic regarding this question, but it is not the same simulation setup and I do not know if, because I am using a solid layer, the answer is maybe a complete other.

Regarding my setup: I fixed the gold and sulfur atoms using position restraints, docked the thiols and ran a simulation at constant temperature over 50k steps.
The gold atoms are not connected to each other, since I did not manage to initialize this correctly at the beginning and since they form a plane anyway due to the position restraints and I did not know how to integrate the connections properly.

If any important information is missing, you have any idea how to solve this or further questions before you get an idea of how this could be solved, please do not hesitate to ask.

I am looking forward to possible solutions.

Kind regards,
Konsti

You do not say how you computed the surface tension. If it is the standard method of the difference of lateral and perpendicular pressure then the answer is that that method does not work in the presence of position restraints.

It might be possible to compute the surface tension difference, which I think is what is measured in the article, by taking the difference in surface tension of the system with and without thiols. But I am not 100% sure that this procedure doesn’t have issues.

I used the “gmx energy Surf*SurfTens” method from gromacs.

But I realised, that in the publication the characteristic/term ‘compressive’ is used with the surface tension. So it is a negative value, isn’t it? Would mean, that this answers my question and my results fit.

You say, that position restraints are not working. Should I then use constraints for the atoms instead?
My goal is to have a layer of gold and above are thiols. So I set the sulfur atoms as a second layer. So I had Au and S atoms with position restraints.
The layer must stay a layer. It can be irregulary, so I need to set any kind of constraints/restraints and I do not want to freeze it as a rigid layer, because this would not be natural.

Do you have any suggestions?

Yes, if it the stress is compressive the surface tension should be negative.

You have not said how you computed the stress. If you take the difference in reported surface tension between a system substrate plus thiols and substrate without tiols, I think the result is always going to be correct.

In my case, the surface tension is the stress created during the simulation that acts on the gold surface.
At the beginning, the thiols are covalently bound to the gold surface (vertically with 90° angle), during the simulation the angles change and the surface tension also changes due to attraction/repulsion between the molecules. This is the surface tension that I analyze at the end.

Or am I somehow misunderstanding your question here?

Sure, that is the physical mechanism.

But if you have position restraints, which are unphysical, the forces working on this will affect the surface tension as computed by GROMACS. If the substrate doesn’t deform much, the effects of these forces will be canceled out when you take the difference in surface tension between two simulations both with substrate but with different things on top the substrates.

So it is about the comparison between always the same basic system, only under different conditions, which will lead to a slightly different result.
For example, different position restraints, different temperatures, or thiols of different lengths.
However, the basic structure always remains the same. And in practice, it should look like a cantilever that bends either up or down due to the surface tension.
Therefore, the gold should be a solid surface. I didn’t want to freeze the atoms, because that is very unnatural, so I thought that position restraints would be a good solution.

You still didn’t answer my question.
Do you compute the stress as the surface tension from a single simulation or as a difference in surface tension between two simulations?

Ah okay, I mean I compare the calculated surface tension between many simulations.

But I compute the stress always just for one single simulation.
I run an energy minimization, an equilibration and then a production run and then in the end I let gromacs calculate the surface tension using “gmx energy #Surf*SurfTen” function.

So you compare the stress from a single simulation to experiment?

Then the position restraints can affect the computed surface tension. But you could be lucky that the effects cancel out. I would think you need to compute the surface tension for a substrate without anything on it and subtract that from all the values you got. But it could be that that value is close to zero.