# Using "soft" constraints / bias in geometry optimization?

I study molecules weakly bound to nanoclusters and surfaces. During geometry optimization, I often run into problems like for example a small metal cluster behave nicely and has a given, stable geometry for 90% of the calculations, then in the remaining 10% of cases, it reorganizes to some other geometry. On one hand, these problems are expected as there are several different stability points around, often with relatively low transition barrier, on the other hand, it is rather annoying when changing a functional or basis set, or something minor chemical change induces such a transformation and makes impossible any comparison between geometries. One can use a geometry constrain to fix the geometry every time, but I feel kind of uneasy about it.

I was wondering if there is an easy way to introduce an extra force, similarly like biases can be introduced in MD simulations to do meta-dynamics, so it smoothly lowers the energy of a given geometry and ideally makes it stable. Technicalities: I am using molecular codes (Gaussian, ORCA) for such calculations at the DFT level.

• What will be the difference between fixing the geometry and a adding a (fictional) extra force? For me, both are the same.
– Camps
Mar 5, 2021 at 11:51
• I guess I just do not like to fix cartesian coordinates when e.g. I expect changes in bond length (which can happen when you change functional or basis set).
– Greg
Mar 5, 2021 at 12:01
• I think you will find it nearly impossible to do this in a useful way. You may have some luck with applying a constraint, optimizing, then releasing the constraint. This lets you start at a more optimum structure, closer to the local minimum you are trying t o hit Mar 5, 2021 at 15:20
• In Orca, if you are starting the optimization with a guess hessian, you can set the force constant of a particular coordinate (bond, angle, or dihedral or Cartesian). For example %opt Hess_Internal {A 2 1 0 D 2.0} reset 5 end will set the diagonal hessian value for the angle between 2nd, 1st and 0th atom to 2 Eh/Bohr^2 and reset it after 5 optimization steps. You can increase the value to exert more force. Mar 5, 2021 at 17:51
• Mar 7, 2021 at 14:49

If you are adamant on adding forces during geometry optimization, you could borrow formalisms from mechanochemistry. We have implemented such formalism using Python (see OpenMechanochem) and it is designed to work with Atomic Simulation Environment (ASE). This would allow you to patch quantum chemical software (molecular mechanics code too), including Gaussian and ORCA. But personally, I think fixing some important coordinates in your cluster can do the trick.

• Thanks for joining to add an answer!
– Tyberius
Nov 16, 2021 at 5:25
• Thank you! I am checking it out.
– Greg
Nov 16, 2021 at 15:23