zero point energy

[SMA]

Hi
How to calculate zero point energy of H2 and D2 molecule

[mlee]

You can distinguish the H2 and D2(see the H in the periodic table) when you build the molecules. You will need to perform the vibration modes of the atoms. You can add the VibrationalMode in the analysis. https://docs.quantumatk.com/manual/Types/VibrationalMode/VibrationalMode.html

[SMA]

Thank mlee
I have made one input for D2 molecule. Is it correct
# Define elements
elements = [Hydrogen, Hydrogen]

# Define coordinates
cartesian_coordinates = [[ 0.  ,  0.  ,  0.  ],
                         [ 0.  ,  0.  ,  0.75]]*Angstrom

# Set up configuration
molecule_configuration = MoleculeConfiguration(
    elements=elements,
    cartesian_coordinates=cartesian_coordinates
    )

# Define a new class where the H atom has mass 2.00
class H2(Hydrogen):
    @staticmethod
    def atomicMass():
        return 2.00*Units.atomic_mass_unit

PeriodicTable.ALL_ELEMENTS.append(H2)

# Replace all atoms tagged "H2" with hydrogen-2
for index in molecule_configuration.indicesFromTags("H2"):
    molecule_configuration.elements()[index] = H2

# Add tags
molecule_configuration.addTags('H2', [0,1])

[mlee]

Deuterium is Hydrogen2. Instead of Hydrogen, you can directly use Hydrogen2 without any additional part. If you want to confirm the atomic mass, type configuration.atomicMasses() in the console of builder. You can see PhysicalQuantity([2.0, 2.0],amu) for D2.

# Define elements
elements = [Hydrogen2, Hydrogen2]

[SMA]

The version I am using is QuantumATK P-2019.03-SP1.
As per ypur suggestions I used the following input .
#setup H-H geometry
elements = [ Hydrogen2, Hydrogen2]
cartesian_coordinates = [[ 0.0,0.0,0.0],
                         [ 0.0,0.0,0.75]]*Angstrom
molecule_configuration = MoleculeConfiguration(
    elements=elements,
    cartesian_coordinates=cartesian_coordinates
    )

But It is displaying error as given below.

  File "h2-vib.py", line 2, in <module>
    elements = [ Hydrogen2, Hydrogen2]
NameError: name 'Hydrogen2' is not defined
application called MPI_Abort(MPI_COMM_WORLD, 1) - process 9

[mlee]

You are right. It works from R-2020.09 version.
Before R-2020.09, the isotope can define a new element type (class) with modified mass and use that in the element list.
The below tutorial includes the C14 isotope with a new element type.
https://docs.quantumatk.com/tutorials/thermoelectrics_cnt_isotope/thermoelectrics_cnt_isotope.html

[SMA]

I have followed the tutorial. But the effect of changed mass is not reflected in the vibrational mode. It is showing the mode of H2 molecule only.
Please help to resolve the issue.
# Molecule Configuration
# -------------------------------------------------------------

# Define elements
elements = [Hydrogen, Hydrogen]

# Define coordinates
cartesian_coordinates = [[ 0.  ,  0.  ,  0.  ],
                         [ 0.  ,  0.  ,  0.75]]*Angstrom

# Set up configuration
molecule_configuration = MoleculeConfiguration(
    elements=elements,
    cartesian_coordinates=cartesian_coordinates
    )
# Define a new class where the H atom has mass 2.00
class H2(Hydrogen):
    @staticmethod
    def atomicMass():
        return 2.00*Units.atomic_mass_unit

PeriodicTable.ALL_ELEMENTS.append(H2)

# Replace all atoms tagged "H2" with hydrogen-2
for index in molecule_configuration.indicesFromTags("H2"):
    molecule_configuration.elements()[index] = H2

# Add tags
molecule_configuration.addTags('H2', [0,1])

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
numerical_accuracy_parameters = NumericalAccuracyParameters(
    density_mesh_cutoff=80.0*Hartree,
    )

calculator = LCAOCalculator(
    numerical_accuracy_parameters=numerical_accuracy_parameters,
    )

molecule_configuration.setCalculator(calculator)
nlprint(molecule_configuration)
molecule_configuration.update()
nlsave('D2-vib.hdf5', molecule_configuration)

# -------------------------------------------------------------
# Optimize Geometry
# -------------------------------------------------------------
molecule_configuration = OptimizeGeometry(
    molecule_configuration,
    max_forces=0.05*eV/Ang,
    max_steps=200,
    max_step_length=0.2*Ang,
    trajectory_filename='D2-dos-vib.hdf',
    disable_stress=True,
    optimizer_method=LBFGS(),
)
nlsave('D2.nc.dos-vib.hdf5', molecule_configuration)
nlprint(molecule_configuration)

# -------------------------------------------------------------
# Dynamical Matrix
# -------------------------------------------------------------
dynamical_matrix = DynamicalMatrix(
    molecule_configuration,
    filename='D2-vib.hdf5',
    object_id='dynamical_matrix',
    repetitions=Automatic,
    atomic_displacement=0.01*Angstrom,
    acoustic_sum_rule=True,
    finite_difference_method=Central,
    force_tolerance=1e-08*Hartree/Bohr**2,
    processes_per_displacement=1,
    log_filename_prefix='forces_displacement_',
    use_wigner_seitz_scheme=False,
    )
dynamical_matrix.update()

# -------------------------------------------------------------
# Phonon Density Of States
# -------------------------------------------------------------
qpoint_grid = MonkhorstPackGrid()

phonon_density_of_states = PhononDensityOfStates(
    configuration=molecule_configuration,
    dynamical_matrix=dynamical_matrix,
    qpoints=qpoint_grid,
    number_of_bands=None,
    )
nlsave('D2.nc.dos-vib.hdf5', phonon_density_of_states)
nlprint(phonon_density_of_states)

# -------------------------------------------------------------
# Vibrational Mode
# -------------------------------------------------------------
vibrational_mode = VibrationalMode(
    configuration=molecule_configuration,
    dynamical_matrix=dynamical_matrix,
    kpoint_fractional=[0, 0, 0],
    mode_indices=None,
    )
nlsave('D2-vib.hdf5', vibrational_mode)

[mlee]

Unfortunately, I have a same problem in the previous version. Can you update the version of Q-2020.09?

[SMA]

How to update the version of 2019.03 .

[mlee]

If you have a valid license, you can download it in SolvNetPlus.
You need the Q-2020.09 to use Hydrogen2 for D2.

See more information:
https://www.synopsys.com/support/licensing-installation-computeplatforms.html

[SMA]

I have a valid license  for 2019.03. Can I get the Q-2020.09 without any additional charge please

[mlee]

Please contact to your regional sales.

[Tue Gunst]

Hi,
Some comments:
1) You need to define the tags before you replace the elements.
"
# Add tags
molecule_configuration.addTags('H2', [0,1])

# Replace all atoms tagged "H2" with hydrogen-2
for index in molecule_configuration.indicesFromTags("H2"):
    molecule_configuration.elements()[index] = H2
"
So changing the order of "Add tags" and "Replace all atoms tagged".

2) Version dependency:
In the 2019 version the Isotope method from the tutorial only works for a ForceFieldCalculator.
For a LCAOCalculator the new element class is not sufficient as it also requires basis set functions for the new element.
Once could instead use LowLevelEntities to get the dynamical matrix and rescale the result if what you want is the eigenmodes.
For more advanced analysis there is not existing solution for the 2019 version.

In the 2020 version we have implemented actual Isotope element classes. So one can for instance use Hydrogen2 instead of Hydrogen etc.
Therefore this version is required if you need LCAOCalculator based studies with advanced analysis options.

[SMA]

Thanks.