I want to write code to simulate the impingement of any molecule on a surface at different (r, theta, and phi) angles in XYZ format.

  • $\begingroup$ Could you be more specific in you question? That might make it more clear for others to provide you an answer $\endgroup$ Jan 11 at 11:14
  • $\begingroup$ to explore surface dynamics by simulating the dissociation of H2 molecules on a Pt surface. To initiate the study, I need to generate initial molecular geometries (trajectories) with varying angles (theta, phi) and distances between H2 atoms so, how to make initial trajectories (via python code or via ASE) for VASP? $\endgroup$ Jan 11 at 15:22

1 Answer 1


I don't totally understand what OP is trying to do, but to my understanding, he/she is trying to make .xyz files for of a molecule moving on the surface of some structure, where its constituent atoms are moving away from each other along with rotation in their bond.

The logic behind that would be to either read or build the required surface model, define a range of values for (theta, phi), and a range for the distances between the atoms of the molecule. We then create the molecule with some initial configuration and iterate over the range of corresponding values we have already defined so that we can vary the distances and rotations, we shall save those configurations accordingly along with the surface we have already built or read in the first steps. Last step, write those configurations as xyz files.

The following is an example of my idea, which might be a good starting point for what OP is willing to achieve. The values in the script are random so you need to tweak according to your situation.

import numpy as np
from ase import Atoms, Atom
from ase.build import fcc111, molecule
from ase.io import write
# Supposed we have Pt surface 
a = 3.92  # Check the correct value of Pt lattice constant in Å
surface = fcc111('Pt', a=a, size=(2, 2, 4), vacuum=10.0, periodic=True)# randomly chosen Pt surface, check for your case or read the surface file if you aiready have one
theta_range = np.linspace(0, np.pi, num=20)  # Range of theta angles (0 to pi), num can change the sampling inside the chosen range
phi_range = np.linspace(0, 2*np.pi, num=40)  # Range of phi angles (0 to 2*pi), same thing about num 
distance_range = np.linspace(0.5, 1, num=10)  # Range of distances (value1 to value2 Å) between the H atoms
# This list is used to store the final xyz files
configurations = []
# Initial trajectories
for theta in theta_range:
    for phi in phi_range:
        for distance in distance_range:
            # Generate H2 molecule with specific distance and orientation
            h2 = molecule('H2')
            h2.set_distance(0, 1, distance)
            h2.rotate(theta, 'y')
            h2.rotate(phi, 'z')
            # Get the positions of H2 atoms relative to Pt surface
            h2_positions = h2.positions + [0.5, 0.5, 2.0 + distance] 
            combined_system = surface.copy()
            combined_system.positions[-2:] = h2_positions
# write the final structures XYZ files
for i, configuration in enumerate(configurations):
    filename = f'configuration_{i}.xyz'
    write(filename, configuration, format='xyz') 
  • $\begingroup$ Thank you @Jaafar Mehrez $\endgroup$ Jan 17 at 13:49

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .