3. Vibrational spectra of Si#
Here we compute the infrared and Raman spectra of silicon using aiida-vibroscopy.
3.1. Finite electric fields#
In the Placzek approximations (good for insulators), all the needed quantities for Raman can be computed as finite difference using the electric fields. Now we a more challenging task, as a third order derivative is needed: the Raman tensor. This is defined as:
(1)#\[\begin{equation}
\frac{\partial \chi_{ij}}{\partial \tau_{K,k}}
=
\frac{1}{\Omega}
\frac{\partial^2 F_{K,k}}{\partial \mathcal{E}_i \partial \mathcal{E}_j}
\end{equation}\]
The same theory we saw in previous tutorials applies here. Note that for computing second order derivatives of forces we need a slightly different formula:
(2)#\[\begin{equation}
\left . \frac{\partial^2 f(x)}{\partial x^2} \right|_{x=0}
=
\frac{1}{h^2}
\left [
f(h) -2f(0) +f(-h)
\right ]
+ \mathcal{O}(h^2)
\end{equation}\]
Have a look at the finite difference coefficients for coefficients of any accuracy order.
Show cell content
from local_module import load_temp_profile
# If you download this file, you can run it with your own profile.
# Put these lines instead:
# from aiida import load_profile
# load_profile()
data = load_temp_profile(
name="raman-tutorial",
add_computer=True,
add_pw_code=True,
add_sssp=True,
)
/opt/conda/lib/python3.8/site-packages/paramiko/transport.py:219: CryptographyDeprecationWarning: Blowfish has been deprecated
"class": algorithms.Blowfish,
3.2. The IRamanSpectraWorkChain#
For computing the spectra we need both phonons and Raman tensors.
Let’s import the WorkChain and run it! We use the get_builder_from_protocol to get a prefilled builder with all inputs.
Note
These inputs should be considered not as converged parameters, but as a good starting point. You may also need to tweak some parameters, e.g. add magnetization etc., depending on your case.
Show cell content
from aiida.plugins import DbImporterFactory
CodDbImporter = DbImporterFactory('cod')
cod = CodDbImporter()
results = cod.query(id='1526655') # Si 1526655
structure = results[0].get_aiida_structure() # it has 8 atoms
from aiida.plugins import WorkflowFactory
from aiida.engine import run_get_node
IRamanSpectraWorkChain = WorkflowFactory("vibroscopy.spectra.iraman")
builder = IRamanSpectraWorkChain.get_builder_from_protocol(
code=data.pw_code,
structure=structure,
protocol="fast",
)
results, calc = run_get_node(builder)
Show cell output
08/25/2023 09:57:41 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [106|IRamanSpectraWorkChain|run_spectra]: submitting `HarmonicWorkChain` <PK=108>
08/25/2023 09:57:41 PM <1014064> aiida.engine.processes.functions: [WARNING] function `generate_preprocess_data` has invalid type hints: unsupported operand type(s) for |: 'AbstractNodeMeta' and 'NoneType'
08/25/2023 09:57:41 PM <1014064> aiida.engine.processes.functions: [WARNING] function `generate_phonopy_data` has invalid type hints: unsupported operand type(s) for |: 'AbstractNodeMeta' and 'NoneType'
08/25/2023 09:57:42 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [108|HarmonicWorkChain|run_phonon]: submitting `PhononWorkChain` <PK=114>
08/25/2023 09:57:42 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [108|HarmonicWorkChain|run_dielectric]: submitting `DielectricWorkChain` <PK=118>
08/25/2023 09:57:43 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_base_scf]: launching base scf PwBaseWorkChain<131>
08/25/2023 09:57:44 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|run_process]: launching PwCalculation<135> iteration #1
08/25/2023 09:57:44 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [114|PhononWorkChain|run_base_supercell]: launching base supercell scf PwBaseWorkChain<140>
08/25/2023 09:57:44 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|run_process]: launching PwCalculation<143> iteration #1
08/25/2023 09:57:54 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|sanity_check_insufficient_bands]: PwCalculation<135> run with smearing and highest band is occupied
08/25/2023 09:57:54 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|sanity_check_insufficient_bands]: BandsData<147> has invalid occupations: Occupation of 0.009307178053852729 at last band lkn<0,0,20>
08/25/2023 09:57:54 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|sanity_check_insufficient_bands]: PwCalculation<135> had insufficient bands
08/25/2023 09:57:54 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|sanity_check_insufficient_bands]: Action taken: increased number of bands to 24 and restarting from the previous charge density.
08/25/2023 09:57:54 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|inspect_process]: PwCalculation<135> finished successfully but a handler was triggered, restarting
08/25/2023 09:57:54 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|run_process]: launching PwCalculation<156> iteration #2
08/25/2023 09:57:56 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|sanity_check_insufficient_bands]: PwCalculation<143> run with smearing and highest band is occupied
08/25/2023 09:57:56 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|sanity_check_insufficient_bands]: BandsData<151> has invalid occupations: Occupation of 0.009307178051672582 at last band lkn<0,0,20>
08/25/2023 09:57:56 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|sanity_check_insufficient_bands]: PwCalculation<143> had insufficient bands
08/25/2023 09:57:56 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|sanity_check_insufficient_bands]: Action taken: increased number of bands to 24 and restarting from the previous charge density.
08/25/2023 09:57:56 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|inspect_process]: PwCalculation<143> finished successfully but a handler was triggered, restarting
08/25/2023 09:57:56 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|run_process]: launching PwCalculation<160> iteration #2
08/25/2023 09:58:03 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|results]: work chain completed after 2 iterations
08/25/2023 09:58:03 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [131|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:58:04 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|results]: work chain completed after 2 iterations
08/25/2023 09:58:04 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [140|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:58:04 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_nscf]: launching base scf PwBaseWorkChain<173>
08/25/2023 09:58:04 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [114|PhononWorkChain|run_forces]: submitting `PwBaseWorkChain` <PK=176> with supercell n.o 1
08/25/2023 09:58:05 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [173|PwBaseWorkChain|run_process]: launching PwCalculation<179> iteration #1
08/25/2023 09:58:06 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|run_process]: launching PwCalculation<182> iteration #1
08/25/2023 09:58:14 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [173|PwBaseWorkChain|results]: work chain completed after 1 iterations
08/25/2023 09:58:14 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [173|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:58:16 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_null_field_scfs]: launching PwBaseWorkChain<201> with null electric field 0
08/25/2023 09:58:16 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_null_field_scfs]: launching PwBaseWorkChain<202> with null electric field 1
08/25/2023 09:58:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [201|PwBaseWorkChain|run_process]: launching PwCalculation<205> iteration #1
08/25/2023 09:58:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [202|PwBaseWorkChain|run_process]: launching PwCalculation<208> iteration #1
08/25/2023 09:58:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|sanity_check_insufficient_bands]: PwCalculation<182> run with smearing and highest band is occupied
08/25/2023 09:58:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|sanity_check_insufficient_bands]: BandsData<193> has invalid occupations: Occupation of 0.008059168786845428 at last band lkn<0,0,20>
08/25/2023 09:58:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|sanity_check_insufficient_bands]: PwCalculation<182> had insufficient bands
08/25/2023 09:58:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|sanity_check_insufficient_bands]: Action taken: increased number of bands to 24 and restarting from the previous charge density.
08/25/2023 09:58:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|inspect_process]: PwCalculation<182> finished successfully but a handler was triggered, restarting
08/25/2023 09:58:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|run_process]: launching PwCalculation<211> iteration #2
08/25/2023 09:58:35 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|results]: work chain completed after 2 iterations
08/25/2023 09:58:36 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [176|PwBaseWorkChain|on_terminated]: cleaned remote folders of calculations: 182 211
08/25/2023 09:58:37 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [114|PhononWorkChain|on_terminated]: cleaned remote folders of calculations: 143 160 182 211
08/25/2023 09:58:40 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [201|PwBaseWorkChain|results]: work chain completed after 1 iterations
08/25/2023 09:58:40 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [201|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:58:40 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [202|PwBaseWorkChain|results]: work chain completed after 1 iterations
08/25/2023 09:58:41 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [202|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:58:42 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_electric_field_scfs]: launching PwBaseWorkChain<232> with electric field index 2 and sign 1.0 iteration #0
08/25/2023 09:58:43 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_electric_field_scfs]: launching PwBaseWorkChain<235> with electric field index 3 and sign 1.0 iteration #0
08/25/2023 09:58:44 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [232|PwBaseWorkChain|run_process]: launching PwCalculation<238> iteration #1
08/25/2023 09:58:44 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [235|PwBaseWorkChain|run_process]: launching PwCalculation<241> iteration #1
08/25/2023 09:59:15 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [232|PwBaseWorkChain|results]: work chain completed after 1 iterations
08/25/2023 09:59:15 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [232|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:59:17 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [235|PwBaseWorkChain|results]: work chain completed after 1 iterations
08/25/2023 09:59:18 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [235|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:59:19 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_electric_field_scfs]: launching PwBaseWorkChain<254> with electric field index 2 and sign 1.0 iteration #1
08/25/2023 09:59:20 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_electric_field_scfs]: launching PwBaseWorkChain<257> with electric field index 3 and sign 1.0 iteration #1
08/25/2023 09:59:21 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [254|PwBaseWorkChain|run_process]: launching PwCalculation<260> iteration #1
08/25/2023 09:59:21 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [257|PwBaseWorkChain|run_process]: launching PwCalculation<263> iteration #1
08/25/2023 09:59:42 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [254|PwBaseWorkChain|results]: work chain completed after 1 iterations
08/25/2023 09:59:42 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [254|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:59:56 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [257|PwBaseWorkChain|results]: work chain completed after 1 iterations
08/25/2023 09:59:56 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [257|PwBaseWorkChain|on_terminated]: remote folders will not be cleaned
08/25/2023 09:59:58 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|run_numerical_derivatives]: launching NumericalDerivativesWorkChain<281> for computing numerical derivatives.
/aiida-plugins/aiida-vibroscopy/src/aiida_vibroscopy/calculations/symmetry.py:163: UserWarning: Symmetry of dChi/dr tensors is largely broken.
The max difference is:
22.241965364020736
warnings.warn('\n'.join(lines))
/aiida-plugins/aiida-vibroscopy/src/aiida_vibroscopy/calculations/symmetry.py:163: UserWarning: Symmetry of dChi/dr tensors is largely broken.
The max difference is:
95.91847563234091
warnings.warn('\n'.join(lines))
/aiida-plugins/aiida-vibroscopy/src/aiida_vibroscopy/calculations/symmetry.py:163: UserWarning: Symmetry of dChi/dr tensors is largely broken.
The max difference is:
120.86345762045994
warnings.warn('\n'.join(lines))
08/25/2023 10:00:02 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [118|DielectricWorkChain|on_terminated]: cleaned remote folders of calculations: 135 156 179 205 208 238 241 260 263
08/25/2023 10:00:06 PM <1014064> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [108|HarmonicWorkChain|on_terminated]: cleaned remote folders of calculations: 143 160 182 211 135 156 179 205 208 238 241 260 263
These are the results:
results
{'output_dielectric': AttributeDict({'fields_data': AttributeDict({'field_index_2': AttributeDict({'0': <TrajectoryData: uuid: 6cdaca66-3cfb-4012-b557-65e61d348523 (pk: 246)>, '1': <TrajectoryData: uuid: 132dbb1d-ddf5-4a47-b158-9c1f07209727 (pk: 268)>}), 'field_index_3': AttributeDict({'0': <TrajectoryData: uuid: c5d3dafe-008e-4400-ba6f-762e3b15121d (pk: 250)>, '1': <TrajectoryData: uuid: 534d19a6-c13e-4d69-9086-f5a8a9988735 (pk: 272)>})}), 'tensors': AttributeDict({'numerical_accuracy_2_step_2': <ArrayData: uuid: 14d990ca-e2d7-4e03-a1cd-430ebec2bb70 (pk: 294)>, 'numerical_accuracy_2_step_1': <ArrayData: uuid: 3b4c38ee-175a-495d-bb72-d025c3e84aa9 (pk: 296)>, 'numerical_accuracy_4': <ArrayData: uuid: c4290950-7f09-4707-9dd5-8cbeffe44e72 (pk: 298)>}), 'critical_electric_field': <Float: uuid: 030bb7bd-3a2f-440f-95cd-58d6fdd023e1 (pk: 191) value: 0.00091230365735815>, 'electric_field_step': <Float: uuid: bf88df31-621a-42a6-8ba6-49fe1bb0d443 (pk: 198) value: 0.00045615182867907>, 'units': <Dict: uuid: 1b723df7-c45f-42a0-b439-9a397922a204 (pk: 292)>}),
'output_phonon': AttributeDict({'supercells': AttributeDict({'supercell_1': <StructureData: uuid: 12e71b8f-f0c1-4936-b200-4ca4a70947ba (pk: 174)>}), 'supercells_forces': AttributeDict({'forces_0': <TrajectoryData: uuid: b69e3687-9176-4e56-9963-4f459bcd99d6 (pk: 168)>, 'forces_1': <TrajectoryData: uuid: 69554bf6-766a-4124-a806-868c99b55db8 (pk: 217)>}), 'phonopy_data': <PhonopyData: uuid: 56075619-37dc-4a33-9c36-30cdd32c6f8a (pk: 221)>}),
'vibrational_data': AttributeDict({'numerical_accuracy_2_step_2': <VibrationalData: uuid: dd5581f7-843f-4dee-ac70-ac24a7e0da40 (pk: 302)>, 'numerical_accuracy_2_step_1': <VibrationalData: uuid: d67d33d4-2d61-434b-a190-306d5d9a7d63 (pk: 306)>, 'numerical_accuracy_4': <VibrationalData: uuid: 07be12d7-d563-4e75-a286-1f1a64b6f02d (pk: 310)>})}
3.3. The VibrationalData#
In aiida-vibroscopy we design a data type able to store all the phonon and dielectric properties, to ease the share of the data and make it as costistent as possible.
In this data, you can access to multiple information:
Structure (cell, sites, …)
Symmetry
Dielectric tensor (\(\epsilon^{\infty}\))
Born effective charges (\(Z^*\))
Raman tensors (\(\partial \chi / \partial \tau\))
Non-linear optical susceptibility (\(\chi^{(2)}\))
Again, Born effective charges and non-linear optical susceptibility are null in silicon, due to symmetry. In the next tutorial we will compute these quantities for AlAs, which has non vanishing tensors.
vibro = calc.outputs.vibrational_data.numerical_accuracy_4
print("The VibrationalData: ", vibro, "\n")
print("The dielectric tensor: ", "\n", vibro.dielectric.round(5), "\n")
print("The Raman tensors (in 1/Angstrom): ", "\n", vibro.raman_tensors.round(5), "\n")
The VibrationalData: uuid: 07be12d7-d563-4e75-a286-1f1a64b6f02d (pk: 310)
The dielectric tensor:
[[11.08481 -0. -0. ]
[ 0. 11.08481 0. ]
[ 0. -0. 11.08481]]
The Raman tensors (in 1/Angstrom):
[[[[-0. -0. -0. ]
[-0. -0. -0.08896]
[-0. -0.08896 -0. ]]
[[ 0. -0. -0.08896]
[ 0. -0. -0. ]
[-0.08896 -0. -0. ]]
[[-0. -0.08896 -0. ]
[-0.08896 -0. -0. ]
[-0. -0. -0. ]]]
[[[ 0. -0. 0. ]
[ 0. -0. 0.08896]
[ 0. 0.08896 0. ]]
[[ 0. -0. 0.08896]
[-0. -0. -0. ]
[ 0.08896 -0. 0. ]]
[[ 0. 0.08896 0. ]
[ 0.08896 0. 0. ]
[ 0. 0. 0. ]]]]
And we can also check how much the Raman tensors, which are third order derivatives, are sensible to the numerical accuracy chosen. Usually, the Raman tensor is expressed in Å\(^2\). To get it, one can use the convention in the code:
(3)#\[\begin{equation}
\left [
\frac{\partial \chi}{\partial \tau}
\right] (\text{Å}^2)
=
\Omega_{unitcell}
\left [
\frac{\partial \chi}{\partial \tau}
\right] (\text{Å}^{-1})
\end{equation}\]
vol = calc.outputs.vibrational_data.numerical_accuracy_4.get_unitcell().get_cell_volume()
print(
calc.outputs.vibrational_data.numerical_accuracy_4.raman_tensors[1,0,1,2]*vol,
calc.outputs.vibrational_data.numerical_accuracy_2_step_1.raman_tensors[1,0,1,2]*vol,
calc.outputs.vibrational_data.numerical_accuracy_2_step_2.raman_tensors[1,0,1,2]*vol,
)
13.860430602407677 13.89715997246768 14.007348082647715
We can see that the values are converged within 0.1 Å\(^2\) with all three steps.
3.4. Plotting the powder Raman spectra#
Now, we can get directly from the data the Raman intensities.
hh_intensities, hv_intensities, frequencies, labels = vibro.run_powder_raman_intensities(
frequency_laser=532, temperature=300)
print(frequencies, "\n", labels)
[531.91425155 531.91425155 531.91425155]
['T2g', 'T2g', 'T2g']
We have 3 degenerate modes. Quite expected for silicon, which is very symmetric. To plot the actual spectra:
from aiida_vibroscopy.utils.plotting import get_spectra_plot
total_intensities = hh_intensities + hv_intensities
plt = get_spectra_plot(frequencies, total_intensities)
plt.show()
3.5. Analysing the workflow#
Look how the complexity of computing a Raman spectra has been handled fully automatically.
%verdi process status {calc.pk}
IRamanSpectraWorkChain<106> Finished [0] [2:if_(should_run_average)]
└── HarmonicWorkChain<108> Finished [0] [4:if_(should_run_phonopy)]
├── generate_preprocess_data<109> Finished [0]
├── PhononWorkChain<114> Finished [0] [7:if_(should_run_phonopy)]
│ ├── generate_preprocess_data<119> Finished [0]
│ ├── get_supercell<132> Finished [0]
│ ├── create_kpoints_from_distance<137> Finished [0]
│ ├── PwBaseWorkChain<140> Finished [0] [3:results]
│ │ ├── PwCalculation<143> Finished [0]
│ │ └── PwCalculation<160> Finished [0]
│ ├── get_supercells_with_displacements<170> Finished [0]
│ ├── PwBaseWorkChain<176> Finished [0] [3:results]
│ │ ├── PwCalculation<182> Finished [0]
│ │ └── PwCalculation<211> Finished [0]
│ └── generate_phonopy_data<220> Finished [0]
├── DielectricWorkChain<118> Finished [0] [11:results]
│ ├── create_kpoints_from_distance<120> Finished [0]
│ ├── create_directional_kpoints<124> Finished [0]
│ ├── create_directional_kpoints<127> Finished [0]
│ ├── PwBaseWorkChain<131> Finished [0] [3:results]
│ │ ├── PwCalculation<135> Finished [0]
│ │ └── PwCalculation<156> Finished [0]
│ ├── PwBaseWorkChain<173> Finished [0] [3:results]
│ │ └── PwCalculation<179> Finished [0]
│ ├── compute_critical_electric_field<189> Finished [0]
│ ├── get_accuracy_from_critical_field<192> Finished [0]
│ ├── get_electric_field_step<197> Finished [0]
│ ├── PwBaseWorkChain<201> Finished [0] [3:results]
│ │ └── PwCalculation<205> Finished [0]
│ ├── PwBaseWorkChain<202> Finished [0] [3:results]
│ │ └── PwCalculation<208> Finished [0]
│ ├── PwBaseWorkChain<232> Finished [0] [3:results]
│ │ └── PwCalculation<238> Finished [0]
│ ├── PwBaseWorkChain<235> Finished [0] [3:results]
│ │ └── PwCalculation<241> Finished [0]
│ ├── PwBaseWorkChain<254> Finished [0] [3:results]
│ │ └── PwCalculation<260> Finished [0]
│ ├── PwBaseWorkChain<257> Finished [0] [3:results]
│ │ └── PwCalculation<263> Finished [0]
│ ├── subtract_residual_forces<276> Finished [0]
│ └── NumericalDerivativesWorkChain<281> Finished [0] [None]
│ ├── generate_preprocess_data<282> Finished [0]
│ ├── compute_nac_parameters<284> Finished [0]
│ ├── compute_susceptibility_derivatives<288> Finished [0]
│ ├── join_tensors<293> Finished [0]
│ ├── join_tensors<295> Finished [0]
│ └── join_tensors<297> Finished [0]
├── elaborate_tensors<299> Finished [0]
├── generate_vibrational_data_from_phonopy<301> Finished [0]
├── elaborate_tensors<303> Finished [0]
├── generate_vibrational_data_from_phonopy<305> Finished [0]
├── elaborate_tensors<307> Finished [0]
└── generate_vibrational_data_from_phonopy<309> Finished [0]
3.6. Convergence with k points#
As a final remark, we show here how the convergence with k points is speeded up by sampling in this directional manner, compared to the most ‘traditional’ uniform mesh.
import numpy as np
import matplotlib.pyplot as plt
# Some options to make the plot nice
plt.rcParams['font.size'] = 18
plt.rcParams['font.family'] = "serif"
plt.rcParams['text.usetex'] = False
plt.rcParams['xtick.major.size'] = 7.0
plt.rcParams['xtick.minor.size'] = 4.
plt.rcParams['ytick.major.size'] = 7.0
plt.rcParams['ytick.minor.size'] = 4.
plt.rcParams['axes.linewidth'] = 1.2
plt.rcParams['legend.fontsize'] = 15
# Directional sampling
raman_new = np.array([18.445, 19.943, 20.442, 20.442])
kpoints_new = np.array([3*6*6, 3*12*12, 3*24*24, 3*58*58])
# Uniform sampling
raman = np.array([15.610, 18.762, 20.436])
kpoints = np.array([3*3*3, 6*6*6, 12*12*12])
# Plotting
fig, ax = plt.subplots(figsize=(6.5,5)) # dichiaro il canvas e i subplots
ax.plot(kpoints_new, raman_new, color='red', marker='o', fillstyle='none', ms=10, linewidth=1.5, linestyle='-', label='Directional sampling')
ax.plot(kpoints, raman, color='black', marker='o', fillstyle='none', ms=10, linewidth=1.5, linestyle='-', label='Uniform sampling')
ax.set_xlabel('# k points')
ax.set_ylabel(r'Raman tensor ($\AA^2$)')
ax.legend(bbox_to_anchor=(0., 1.02, 1., .102), loc='lower left',
ncol=1, mode="expand", borderaxespad=0., frameon=False)
<matplotlib.legend.Legend at 0x7fa49aad0b20>
