System | Charge | Multiplicity |
---|---|---|
neutral molecule | 0 | 1 |
protonated PAHs C16H15+ | 1 | 1 |
cation H2O+ | 1 | 2 |
radical cation PAHs C16H14+ | 1 | 2 |
- Linear molecules have to align with z-axis for vibrational analysis |
$$\text{multiplicity} = \text{number of } \alpha \text{ electrons} - \text{number of } \beta \text{ electrons} + 1$$ - When #electrons is odd, then multiplicity must be even (2, 4, 6, ...). - When #electrons is even, then multiplicity must be old (1, 3, 5, ...).
By default, TRAH
is launched automatically if standard SCF (DISS/SOSCF) fails. It is possible to turn this option off by adding
! NOTRAH
or
%scf
AutoTRAH false
end
Unit conversion for electric field is 1 a.u. = 51.4220674763 V/Å (1 * ase.units.Ha / ase.units.Bohr
)
%scf efield X-field-strength-in-a.u. , Y-field-strength-in-a.u. , Z-field-strength-in-a.u. end
NOPOP
: Turns off all population analysis (Mulliken, Loewdin, Mayer)M06-2x
should use D3zero
instead of D3BJ
D3不应当理解为是用来描述色散作用的,而是对你选择的泛函基础上进行校正以提高精度,只不过是对于色散作用描述极差的诸如B3LYP,D3主要作用才是起到描述色散作用的作用。D3的参数对每种泛函是不同的。M06-2X已经能很好描述色散作用,所以加了D3只会稍微改进一点精度,完全不会double-counting。明尼苏达系列泛函不能用BJ阻尼的D3是因为BJ阻尼下必然导致double-counting。这些问题Grimme都已经充分考虑了。 - [量化理论] M06-2X加D3校正的问题
Grid6
for minnesota meta-GGAUsing MN15 functional in ORCA 4.2.1 from libxc
! RIJCOSX def2-TZVP def2/J Opt nopop %scf maxiter 600 end %method method dft exchange hyb_mgga_x_mn15 correlation mgga_c_mn15 end %maxcore 1000 %pal nprocs 16 end *xyz 0 1 C -0.047131 0.664389 0.0 O -0.047131 -0.758551 0.0 H -1.092995 0.969785 0.0 H 0.878534 -1.048458 0.0 H 0.437145 1.080376 0.891772 H 0.437145 1.080376 -0.891772 *
Complete Active Space Self-Consistent Field (CASSCF) - CASSCF(N,M): N-electrons are active in M orbitals.
%casscf nel N norb M end
ORCA 5.0.0 and 5.0.1 has a bug in interpreting the second Colvar. This can be solved by adding 3 more zeros in the end.
Metadynamics Colvar 1 Scale 1.0_A Wall Lower 0.0 50.0 Upper 10.0 50.0 Range 0.0 10.0 100 Metadynamics Colvar 2 Scale 1.0_A Wall Lower 0.0 50.0 Upper 10.0 50.0 Range 0.0 10.0 100 0 0 0
Alternatively, putting two Colvars in the same line also solves this bug
Metadynamics Colvar 1 Scale 1.0_A Wall Lower 0.0 50.0 Upper 10.0 50.0 Colvar 2 Scale 1.0_A Wall Lower 0.0 50.0 Upper 10.0 50.0 Range 0.0 10.0 100 0.0 10.0 100
For a big system (e.g. around 200 atoms), the local optimization in new redundant internal coordinates can be very slow.
For example, updating positions could take several minutes.
Diagonalizing the G-matrix is also slow.
- %geom MaxIter 200 end
sets the maximum optimization cycle to 200.
Observation | Parameter | TS OPTIMIZATION |
---|---|---|
Energy Change | TolE | 5.0000e-06 Eh |
Max. Gradient | TolMAXG | 3.0000e-04 Eh/bohr |
RMS Gradient | TolRMSG | 1.0000e-04 Eh/bohr |
Max. Displacement | TolMAXD | 4.0000e-03 bohr |
RMS Displacement | TolRMSD | 2.0000e-03 bohr |
In a TS optimization we are looking for a first order saddle point, and thus for a point on the PES where the curvature is negative in the direction of the TS mode (the TS mode is also called transition state vector, the only eigenvector of the Hessian at the TS geometry with a negative eigenvalue). As transition state finder we implemented the quasi-Newton like hessian mode following algorithm.This algorithm maximizes the energy with respect to one (usually the lowest) eigenmode and minimizes with respect to the remaining 3N - 7(6) eigenmodes of the Hessian.
OptTS
uses the eigenvector following algorithm to find the nearest stationary point on PES. This can be local minimum or transition state.
PreOpt_Ends
conflicts with NEB_Restart_XYZFile
.
You cannot use them at the same time.
NEB-TS
doesn't calculate the Hessian (Initial Hessian is stored in orca_TSOpt.appr.hess
, no idea where this if from)zoomNEB is similar to AutoNEB to some extents. Both methods insert images around the saddle point region on the fly. But zoomNEB inserts several images as many as initial images once identifying the saddle point region, while AutoNEB insert images one by one until a maximal pre-defined number.
In orca.hess
- $hessian
has the 2nd energy derivatives w.r.t cartesian coordinates (not mass-weighted yet), with the unit of Ha / Bohr**2
and a shape of (3N, 3N)
- $dipole_derivatives
has the dipole derivative in cartesian coordinates, with the unit of e*Bohr / Bohr
and a shape of (3N, 3)
import numpy as np import ase.units as units from ase import Atoms def read_orca_hess(filename='orca.hess'): ''' Read Hessian and dipole derevatives from the orca hess output file ''' n_column = 5 with open(filename, 'r') as f: for line in f: if '$hessian' in line: n_dof = int(next(f).strip().split(' ')[0]) n_block = np.ceil(n_dof/n_column).astype('int') hessian = np.zeros((n_dof, n_dof)) for idx_block in range(n_block): next(f) for i in range(n_dof): content = next(f).strip().split(" ")[1:] content = [float(item) for item in content if item != ""] hessian[i, idx_block*n_column:(idx_block+1)*n_column] = content hessian *= units.Ha / (units.Bohr)**2 if '$atoms' in line: n_atoms = int(next(f).strip().split(' ')[0]) symbols = [] masses = [] positions = [] for idx_atom in range(n_atoms): content = next(f).strip().split(" ") symbol = content[0] symbols.append(symbol) info = [float(item) for item in content[1:] if item != ""] mass = info[0] xyz = info[1:] masses.append(mass) positions.append(xyz) positions = np.array(positions) * units.Bohr atoms = Atoms(symbols, positions, masses=masses) if '$dipole_derivatives' in line: n_dof = int(next(f).strip().split(' ')[0]) dipole_derivatives = np.zeros((n_dof, 3)) for idx in range(n_dof): content = next(f).strip().split(" ") content = [float(item) for item in content if item != ""] dipole_derivatives[idx] = content n_atoms = len(atoms) masses = atoms.get_masses() mass_weights = np.repeat(masses**-0.5, 3) omega2, vectors = np.linalg.eigh(mass_weights * hessian * mass_weights[:, np.newaxis]) freq_unit_in_SI = np.sqrt(units._e/units._amu) * units.m freq_unit_in_eV = units._hbar * freq_unit_in_SI * units.J freq_in_invcm = omega2.astype(complex)**0.5 * freq_unit_in_eV / units.invcm freq_in_invcm = freq_in_invcm.real * np.sign(omega2) # ASE's version #unit_conversion = units._hbar * units.m / np.sqrt(units._e * units._amu) #energies = unit_conversion * omega2.astype(complex)**0.5 #energies /= units.invcm modes = vectors.T * mass_weights intensities = np.einsum('ab,bg->ag', modes, dipole_derivatives) intensities = np.sum(intensities**2, axis=1).real / units.Debye**2 intensities *= 42.255 for i in range(3*n_atoms): print(f'{freq_in_invcm[i].real:8.2f} cm^-1 {intensities[i]:6.2f} km/mol')
You can get results form a thermochemical analysis (based on ideal gas assumption) once freq
is turned on. The most useful quantity for me is the Gibbs free energy G
which includes contributions from electronic energies, vibrational frequencies and temperatures.
G = H - T*S
VPT2 module is available again since ORCA 5.0.0.
It is turned on by the keyword NEARIR
.
Since the evaluation of third and fourth derivatives are expensive, ORCA uses XTB2 to get high order derivatives by default.
The default settings are writen under %FREQ
:
%FREQ XTBVPT2 True DELQ 0.5 END
DELQ
is the step during the numerical differentiation in dimensionless normal mode units.
If one wants high order derivatives in the theory level as same as the single point calculation, XTBVPT2 False
can be used.
The total number of single point calculations in ORCA-VPT2 is (2N_mode+1)(N_mode+1)
where N_mode = 3*N_atoms - 6
for non-linear molecules.
orca_pltvib
is used for generating animations of vibrational modes.
It takes orca.out
or orca.hess
as arguments.
When analytical frequency is performed, orca.out
is not compatible with orca_pltvib
in ORCA 4.2.1.
grep 'The optimization did not converge' orca.out
One can use the following script to delete all .gbw
, .densities
and .tmp
files in a directory and its subdirectories.
import os rootdir = 'calculations/orca' for subdir, dirs, files in os.walk(rootdir): for file in files: if file.endswith('.gbw') or file.endswith('.densities') or file.endswith('.tmp'): full_file_path = os.path.join(subdir, file) os.remove(full_file_path) #print(full_file_path)
One has to create JOBNAME.nodes
in the running directory to run on multiple nodes.
echo $SLURM_NODELIST >> JOBNAME.nodes
will create s31n[07-09]
which cannot be recognized by ORCA.
The following script gives the correct multiple nodes setting.
#SBATCH --nodelist=s31n[07-09] #SBATCH --ntasks=36 #SBATCH --ntasks-per-node=12 #SBATCH --cpus-per-task=1 # change directory to the local scratch-directory, and run: cd /scratch/$SLURM_JOB_ID cat > orca.nodes <<! s31n07 slots=12 s31n08 slots=12 s31n09 slots=12 !
WARNING: For RIJCOSX frequencies GRIDX4 or larger grids are recommended
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!! FATAL ERROR ENCOUNTERED !!! !!! ----------------------- !!! !!! I/O OPERATION FAILED !!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
I/O OPERATION FAILED
is printed out when running out of disks.
It is better to restart a job at a node with a larger disk.
%method Z_Solver CG end
[file orca_main/mainchk.cpp, line 1630]: ERROR: RI-MP2 needs an AuxC basis but none was defined!
An AuxC basis has to be defined explicitly in the input when using RI-MP2, such as RI-MP2 ANO-pVDZ cc-pVTZ/C TightSCF Engrad
.
[file orca_main/mainfcts.cpp, line 625]: Error : multiplicity (2) is even and number of electrons (204) is even -> impossible
---------------------------------------------------------------------------- ERROR !!! The optimization did not converge but reached the maximum number of optimization cycles. As a subsequent Frequencies calculation has been requested ORCA will abort at this point of the run. Please restart the calculation with the lowest energy geometry and/or a larger maxiter for the geometry optimization. ---------------------------------------------------------------------------- [file orca_main/run.cpp, line 10680]: ORCA finished with error return - aborting the run