1"""Example, in order to run you must place a pseudopotential 'Na.psf' in
2the folder"""
3
4from ase.units import Ry, eV, Ha
5from ase.calculators.siesta import Siesta
6from ase.calculators.siesta.siesta_raman import SiestaRaman
7from ase import Atoms
8import numpy as np
9
10# Define the systems
11# example of Raman calculation for CO2 molecule,
12# comparison with QE calculation can be done from
13# https://github.com/maxhutch/quantum-espresso/blob/master/PHonon/examples/example15/README
14
15CO2 = Atoms('CO2',
16            positions=[[-0.009026, -0.020241, 0.026760],
17                       [1.167544, 0.012723, 0.071808],
18                       [-1.185592, -0.053316, -0.017945]],
19            cell=[20, 20, 20])
20
21# enter siesta input
22siesta = Siesta(
23    mesh_cutoff=150 * Ry,
24    basis_set='DZP',
25    pseudo_qualifier='',
26    energy_shift=(10 * 10**-3) * eV,
27    fdf_arguments={
28        'SCFMustConverge': False,
29        'COOP.Write': True,
30        'WriteDenchar': True,
31        'PAO.BasisType': 'split',
32        'DM.Tolerance': 1e-4,
33        'DM.MixingWeight': 0.01,
34        'MaxSCFIterations': 300,
35        'DM.NumberPulay': 4,
36        'XML.Write': True,
37        'DM.UseSaveDM': True})
38
39CO2.set_calculator(siesta)
40
41ram = SiestaRaman(CO2, siesta, label="siesta", jcutoff=7, iter_broadening=0.15/Ha,
42        xc_code='LDA,PZ', tol_loc=1e-6, tol_biloc=1e-7, freq = np.arange(0.0, 5.0, 0.05))
43ram.run()
44ram.summary(intensity_unit_ram='A^4 amu^-1')
45
46ram.write_spectra(start=200)
47