Keep 2 bond lengths the same

Questions
#1

For teaching purposes, (to show what a PES is, what cuts through the PES look like and why a single cut does not necessarily get you to the minimum energy configuration, etc.), I would like to do a 2D scan of water, keeping the two OH bonds the same, scanning over the HOH angle and the OH bond length (with OH1 = OH2).

As it is for teaching, it would be nice if there were a possibility to do this via a 2D scan, as this would provide me with all results all at once. However, I did not find a possibility to either restrict the second OH bond to the value of the first, or to scan through a variable using a Z-matrix representation with this variable defining the OH bond length. Is there a possibility to do such a 2D scan while restricting the two OH bonds to be the same?

#2

I don’t think it’s possible to constrain equal bond distances, but how about scanning the H-H distance? At each point of the scan, the remaining degrees of freedom are relaxed, which naturally enforces the O-H bonds to be equal. For presentation you will also have to back-convert H-H distances to O-H distances using the H-O-H angle.

$molecule
0 1
O
H 1 0.98
H 1 0.98 2 104.5
$end

$rem
jobtype = pes_scan
method = wb97x
basis = def2-svp
sym_ignore = true
$end

$scan
bend 2 1 3 103.0 106.5 1.5
stre 3 2   1.4 1.6 0.1
$end
#3

I had already thought about this “hack” and had hoped not to need it as it might confuse the students. But as the “direct” way is not possible, this is a very nice solution. Thank you for your answer!

#4

Then the other obvious way that you’ve probably considered, but which I’ll add here for future reference, is to do a sequence of constrained geometry optimizations. This would require the student to set up a compound input file, one input for each data point on the potential curve, e.g.:

$molecule
0 1
O
H 1 1.0
H 1 1.0 2 104.5
$end

$rem
jobtype = opt
method = wb97x
basis = def2-svp
sym_ignore = true
$end

$opt
CONSTRAINT
stre 1 2 1.1
stre 1 3 1.1
ENDCONSTRAINT
$end

@@@

$molecule
read
$end

$rem
scf_guess = read
jobtype = opt
method = wb97x
basis = def2-svp
sym_ignore = true
$end

$opt
CONSTRAINT
stre 1 2 1.2
stre 1 3 1.2
ENDCONSTRAINT
$end

These will work better (less chance of failure) if you do them sequential with small incremental changes in the bond length. The geometry optimizer does not require the constraints to be satisfied in the initial structure but it can have trouble if they are way off.

#5

Thank you for your answer! I am confused though. (Now it is me being confused - not the students :rofl:) Why would I need an optimization in this case? Wouldn’t that give me the optimal angle instead of a 2D surface (angle + OH bond length)? (And if I were to specify 3 constraints (OH1, OH2 and HOH), I could just as well specify the geometry of water itself.)
Clearly, the latter is the easiest for the students: A single single-point calculation. If I have each of them computing 3 points, all together, that should give a very decent PES … Maybe it is actually fun if they can fill in their values into a grid in a shared document.

#6

I just assumed that if you were wanting to map a potential surface in reduced dimensions that you would be optimizing (relaxing) in the other degrees of freedom. Obviously you don’t need to do this, but it is what JOBTYPE = PES_SCAN does. If you want a fixed-angle scan at different OH bond lengths, then you should do this with a compound input file, using a sequence of jobs with Z-matrices to specify the internal coordinates. (You could turn it into a scripting exercise, because that’s what it is in this case.)