Hello all,
I’ve recently begun using Q-Chem EOM-EE-CCSD to to perform excited state geometry optimizations of open-shell diatomic clusters (eg CN-CO). Ideally I would like to use one z-matrix to sweep through Jacobi coordinates for the diatomics (r is the distance between the centers of mass between each pair, and theta being the angle created by CO and the C in CN). I thought I could put a dummy atom at the COM of the CO molecule, and another above the plane created by the CCO atoms (to keep all angles from approaching 180 degrees), but I have not had much luck with my DUMMY section. The below z-matrix (and other attempts) resulted in error messages stating I needed to list referenced centers in a dummy section. Do you have examples of z-matrices that would work better or examples of using the types of dummy atoms?
Thanks,
Benj
O
xx 1 xxo2
C 1 co3 2 coxx3
xx 2 xxxx4 3 xxxxc4 1 dih4
C 3 cc5 1 occ5 4 dih5
N 5 cn6 2 cnxx6 4 dih6
xxo2 0.550000
co3 1.100000
coxx3 1.000
xxxx4 1.000000
xxxxc4 90.000
dih4 179.000
cc5 2.000000
occ5 160.000
dih5 179.000
cn6 1.340000
cnxx6 60.000
dih6 0.000
Can you provide a complete input file? I will try and help with your issues.
Hello,
After reading the manual yet one more time and looking through the output of some failed runs, I think I have a better understanding of how dummy atoms are processed by Q-Chem. I had thought that the atoms were numbered based on their appearance in the z matrix. My current understanding is the real atoms are numbered by their positions in the z matrix, and the dummy atoms then start at one more than the number of real atoms (5 in this case). Below I understood the dummy section to create atom 5 using two atoms, C and O, and the geometric mean distance between them.
Then when I specify constraints, I have to reference the atoms based on this new ordering and not on their orders in the original z matrix. Does this sound correct?
My below z matrix worked properly in the constrained optimization, which is definitely an improvement over previous attempts. The downside of the z matrix I used is it fails near the T shaped structure where the C, N, dummy atom form an angle of 180 degrees. Is there a better z matrix that bypasses having an angle approaching 180 degrees while mapping out the approximate Jacobi angle theta?
$comment
Uses constrained optimizations and zmat f
$end
$molecule
0 2
x
x 1 xxxx2
c 1 cxx3 2 cxxxx3
n 3 nc4 1 ncxx4 2 dih4
c 1 cxx5 3 cxxo5 2 dih5
o 5 co6 1 coxx6 2 dih6
xxxx2 1.000000
cxx3 0.5750000
cxxxx3 90.000
nc4 1.150000
ncxx4 0.000
dih4 0.000
cxx5 3.00000
cxxo5 50.000
dih5 180.000
co6 1.220000
coxx6 180.000
dih6 180.000
$end
$rem
mem_static 1000
mem_total 4000
JOBTYPE Opt
GEOM_OPT_MAX_CYCLES 20
SYM_IGNORE True
SYMMETRY FALSE
CC_SYMMETRY False
METHOD EOM-CCSD
BASIS aug-cc-pvdz
SCF_convergence 8
SCF_algorithm gdm
MaxSCF 1000
SCF_GUESS SAD
gen_scfman true
gen_scfman_final true
CC_T_CONV 9 CCSD amplitudes convergence
EOM_DAVIDSON_CONVERGENCE 9 EOM amplitude convergence
EE_STATES [3]
CC_STATE_TO_OPT [1,1]
EOM_NGUESS_SINGLES 3
$end
$opt
DUMMY
5 1 2 3 4
ENDDUMMY
CONSTRAINT
stre 1 2 1.15
stre 2 3 3.000
stre 3 4 1.22
bend 3 5 1 40.000
ENDCONSTRAINT
$end
Thank you for the input this helps me. I already see several issues that may help your calculation.
First when running the job this portion of the output is useful for you.
**WARNING** Cannot use Internal Coordinates with Fixed Atoms
Cartesian Optimization Enforced
**WARNING** Constraint 4 Involves a Non-Fixed Dummy Atom
Such Constraints are NOT Guaranteed to Converge!
The Best Way to Define a Constraint with Respect to a Dummy Atom
is Via a Z-Matrix
** CONSTRAINED OPTIMIZATION IN CARTESIAN COORDINATES **
Searching for a Minimum
Optimization Cycle: 1
Coordinates (Angstroms)
ATOM X Y Z
1 C 0.5750000000 0.0000000000 0.0000000000
2 N -0.5750000000 0.0000000000 0.0000000000
3 C 1.9283628291 0.0000000000 -2.2981333294
4 O 2.7125637129 0.0000000000 -3.2327075500
5 Du 2.3204632710 0.0000000000 -2.7654204397
Point Group: c1 Number of degrees of freedom: 6
Energy is -205.469479943
Constraints and their Current Values
Value Constraint
Distance(Angs): 1 2 1.150000 1.150000
Distance(Angs): 2 3 3.398270 3.000000
Distance(Angs): 3 4 1.220000 1.220000
Angle: 3 5 1 7.741 40.000
This is printed at the start of the geometry optimization portion. So first for clarity in your input you should add geom_opt_coords = 0 this switches to Cartesian optimization, since you can not do internal coordinates optimization with fixed atoms. Second, the warning about constraint 4 deals with the constraint with a dummy atom, states you are not guaranteed to converge.
But I think you might be over complicating the input. So for clarity I want to try and understand the type of results you are trying to obtain. So you have these two diatomics and you want to 1) optimize their distance from each other (COM-COM or C-C distance) and 2) optimize the angle between the two diatomics ( O-C-C)? Next, what are you expecting the optimized position to be? A T-Shape, a linear molecule (whether it is O-C - C-N or O-C - N-C), or something in-between (this can be were there is an angle between O-C - [COM of C-N] or if the O-C - COM is co-linear and then there exists some dihedral angle between O-C - [COM of C-N] - N)? Also you utilize a constraint for the C-C distance of 3.0, it is not obvious to me why you want this constraint when you are trying to optimize
Adding this bond constraint to me seems like it will cause problems since to me it seems like you are using it to force the optimizer to find an angle for you about this bond.
I think the first simplification you can do is to keep one of the diatomics fixed, for example the C-O. It is not obvious to me why this molecule shouldn’t be fixed using
FIXED
3 XYZ
4 XYZ
ENDFIXED
I do not see an example that you can use. But after you answer my questions above I should be able to help with the input section for this job.
@pfmclaug thank you very much for the detailed response. I will definitely add your suggestions into my input file.
My end goal is to create a map of the COM distances from about 2 to 5 A and the angles (say C-N-C) from 0 to 180 degrees for both the X and the A state. My thought is to use constrained optimizations with different input configurations (T, H, X, L, I), a fixed distance between the diatomics, a fixed C-N-C angle, to find which parts of the PES have minima corresponding to those configurations. Starting with the different configurations seems necessary to ensure the lowest energy geometry is found for a given pair of constraints. The second step is to use this constrained geometry to perform a scan (say from 2 to 3 A between the diatomics). This last step would then be repeated for all of the C-N-C angles to build the surface.
I 100 % agree one of the diatomics can be fixed given I am holding those bond distances constant, thank you very much for the suggestion.
Ok hopefully my suggestions will help. I do want to point out the fact that comparing PES from different constrained optimization configurations with each other is not trivial. If you are only changing the electronic state then you shouldn’t have a problem.
I also think you could simplify your starting position by performing a ground state optimization with minimal constraints. Then with that structure perform the excited state optimizations from that initial guess with your desired constraints. That will speed up your job workflow since you can first obtain a good quality initial geometry guess. Then add in the costly excited state calculations.
1 Like
Thank you very much for all of your suggestions, they were very helpful. I used pes_scan to generate the ground-state surface. Most of the diatomic bond lengths are very similar, and I would like to use pes_scan to map the A surface with the diatomic bond lengths set to those from the ground state. This way I can compare the vertical excitation energies to a more adiabatic representation. One follow up question, is it possible to use constraints with the pes_scan job type?
If you follow Example 9.10 in the Q-Chem manual it does a one-dimensional torsional scan of butane. The job information is provided in section 9.5.
Looking through section 9.5 and example 9.10 definitely fits for what I did with the ground state, which worked great. My apologies for not being clear, in addition to stepping through two variables (distance between N in diatomic 1 and C in diatomic 2, and the angle between the first 3 atoms) using the $scan section I need to constrain two other variables (the bond distances in each diatomic). Hopefully the below example helps clarify what I am after (that does not work because CONSTRAINT does not appear to be available in $scan). I know I could do it if I only needed two variables (instead of four) or if pes_scan could handle four.
$scan
stre 2 3 3.000 4.000 0.1
bend 3 5 1 40.000 50.000 10.0
CONSTRAINT
stre 1 2 1.15
stre 3 4 1.22
ENDCONSTRAINT
$end
So you can only specify a maximum of two stretches/bends/torsions for a PES scan.
The CONSTRAINT option is for $opt section which is not read during the PES scan job.
I think you want to use the FROZEN_SCAN TRUE rem variable which is described in section 9.5, the only stipulation is that the $molecule needs to be in a Z-Mat format.
@pfmclaug thank you so much for all of your help! I was hoping Q-Chem had a secret way to add in constraints to $scan, akin to what can be done in $opt. Oh well. I will continue on with all of your advice, constrained optimizations, and frozen scans. Thanks again!