/
index.html
26 lines (25 loc) · 14.1 KB
/
index.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Input and output formats · DFTK.jl</title><link rel="canonical" href="https://docs.dftk.org/stable/guide/input_output/"/><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="../.."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="../../assets/documenter.js"></script><script src="../../siteinfo.js"></script><script src="../../../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="../../assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="../../assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="../../assets/themeswap.js"></script><link href="../../assets/favicon.ico" rel="icon" type="image/x-icon"/></head><body><div id="documenter"><nav class="docs-sidebar"><a class="docs-logo" href="../../"><img src="../../assets/logo.png" alt="DFTK.jl logo"/></a><div class="docs-package-name"><span class="docs-autofit">DFTK.jl</span></div><form class="docs-search" action="../../search/"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="../../">Home</a></li><li><a class="tocitem" href="../../school2022/">DFTK School 2022</a></li><li><span class="tocitem">Getting started</span><ul><li><a class="tocitem" href="../installation/">Installation</a></li><li><a class="tocitem" href="../tutorial/">Tutorial</a></li><li class="is-active"><a class="tocitem" href>Input and output formats</a><ul class="internal"><li><a class="tocitem" href="#Reading-input-formats-supported-by-ASE"><span>Reading input formats supported by ASE</span></a></li><li><a class="tocitem" href="#Writing-VTK-files-for-visualization"><span>Writing VTK files for visualization</span></a></li><li><a class="tocitem" href="#Writing-and-reading-JLD2-files"><span>Writing and reading JLD2 files</span></a></li></ul></li><li><a class="tocitem" href="../parallelization/">Timings and parallelization</a></li><li><a class="tocitem" href="../periodic_problems/">Introduction to periodic problems</a></li><li><a class="tocitem" href="../density_functional_theory/">Density-functional theory</a></li></ul></li><li><span class="tocitem">Examples</span><ul><li><a class="tocitem" href="../../examples/metallic_systems/">Temperature and metallic systems</a></li><li><a class="tocitem" href="../../examples/pymatgen/">Creating supercells with pymatgen</a></li><li><a class="tocitem" href="../../examples/ase/">Creating slabs with ASE</a></li><li><a class="tocitem" href="../../examples/collinear_magnetism/">Collinear spin and magnetic systems</a></li><li><a class="tocitem" href="../../examples/geometry_optimization/">Geometry optimization</a></li><li><a class="tocitem" href="../../examples/scf_callbacks/">Monitoring self-consistent field calculations</a></li><li><a class="tocitem" href="../../examples/scf_checkpoints/">Saving SCF results on disk and SCF checkpoints</a></li><li><a class="tocitem" href="../../examples/polarizability/">Polarizability by linear response</a></li><li><a class="tocitem" href="../../examples/gross_pitaevskii/">Gross-Pitaevskii equation in one dimension</a></li><li><a class="tocitem" href="../../examples/gross_pitaevskii_2D/">Gross-Pitaevskii equation with external magnetic field</a></li><li><a class="tocitem" href="../../examples/cohen_bergstresser/">Cohen-Bergstresser model</a></li><li><a class="tocitem" href="../../examples/arbitrary_floattype/">Arbitrary floating-point types</a></li><li><a class="tocitem" href="../../examples/forwarddiff/">Polarizability using automatic differentiation</a></li><li><a class="tocitem" href="../../examples/custom_solvers/">Custom solvers</a></li><li><a class="tocitem" href="../../examples/custom_potential/">Custom potential</a></li><li><a class="tocitem" href="../../examples/wannier90/">Wannierization using Wannier90</a></li><li><a class="tocitem" href="../../examples/error_estimates_forces/">Practical error bounds for the forces</a></li></ul></li><li><span class="tocitem">Advanced topics</span><ul><li><a class="tocitem" href="../../advanced/conventions/">Notation and conventions</a></li><li><a class="tocitem" href="../../advanced/data_structures/">Data structures</a></li><li><a class="tocitem" href="../../advanced/useful_formulas/">Useful formulas</a></li><li><a class="tocitem" href="../../advanced/symmetries/">Crystal symmetries</a></li></ul></li><li><a class="tocitem" href="../../api/">API reference</a></li><li><a class="tocitem" href="../../publications/">Publications</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li><a class="is-disabled">Getting started</a></li><li class="is-active"><a href>Input and output formats</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href>Input and output formats</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://github.com/JuliaMolSim/DFTK.jl/blob/master/docs/src/guide/input_output.jl" title="Edit on GitHub"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitHub</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="Input-and-output-formats"><a class="docs-heading-anchor" href="#Input-and-output-formats">Input and output formats</a><a id="Input-and-output-formats-1"></a><a class="docs-heading-anchor-permalink" href="#Input-and-output-formats" title="Permalink"></a></h1><p><a href="https://mybinder.org/v2/gh/JuliaMolSim/DFTK.jl/gh-pages?filepath=v0.5.0/guide/input_output.ipynb"><img src="https://mybinder.org/badge_logo.svg" alt/></a> <a href="https://nbviewer.jupyter.org/github/JuliaMolSim/DFTK.jl/blob/gh-pages/v0.5.0/guide/input_output.ipynb"><img src="https://img.shields.io/badge/show-nbviewer-579ACA.svg" alt/></a></p><p>This section provides an overview of the input and output formats supported by DFTK, usually via integration with a third-party library.</p><h2 id="Reading-input-formats-supported-by-ASE"><a class="docs-heading-anchor" href="#Reading-input-formats-supported-by-ASE">Reading input formats supported by ASE</a><a id="Reading-input-formats-supported-by-ASE-1"></a><a class="docs-heading-anchor-permalink" href="#Reading-input-formats-supported-by-ASE" title="Permalink"></a></h2><p>ASE is short for the <a href="https://wiki.fysik.dtu.dk/ase/index.html">atomistic simulation environment</a>, a Python package to simplify the process of setting up, running and analysing results from atomistic simulations across different programs. If it is installed it is automatically used by DFTK in order to read a wide range of input files, such as xyz, CIF, input files of various other codes (e.g. Quantum Espresso, VASP, ABINIT, CASTEP, ...). The full list of formats can be found in the <a href="https://wiki.fysik.dtu.dk/ase/ase/io/io.html">ASE IO documentation</a>.</p><p>As an example we start the calculation of a simple antiferromagnetic iron crystal using a Quantum-Espresso input file, <a href="../Fe_afm.pwi">Fe_afm.pwi</a>. From this file the lattice, atomic positions and the initial magnetisation are read. For more details about calculations on magnetic systems using collinear spin, see <a href="../../examples/collinear_magnetism/#Collinear-spin-and-magnetic-systems">Collinear spin and magnetic systems</a>.</p><p>First we parse the Quantum Espresso input file to DFTK datastructures:</p><pre><code class="language-julia">using DFTK
qe_input = "Fe_afm.pwi"
atoms = load_atoms(qe_input)
positions = load_positions(qe_input)
lattice = load_lattice(qe_input)
magnetic_moments = load_magnetic_moments(qe_input)</code></pre><pre class="documenter-example-output">2-element Vector{StaticArrays.SVector{3, Float64}}:
[0.0, 0.0, -9.6]
[0.0, 0.0, -9.6]</pre><p>At this point a file of any format supported by ASE could be passed instead, e.g. an <code>xyz</code> file or an ABINIT input file. Behind the scenes ASE takes care to select the right parser and extract the required structural information.</p><p>Next we attach the pseudopotential information, since this information is currently not exposed inside the ASE datastructures. See <a href="../../examples/ase/#Creating-slabs-with-ASE">Creating slabs with ASE</a> for more details.</p><pre><code class="language-julia">atoms = map(atoms) do el
@assert el.symbol == :Fe
ElementPsp(:Fe, psp=load_psp("hgh/pbe/fe-q16.hgh"))
end;</code></pre><p>Finally we run the calculation.</p><pre><code class="language-julia">model = model_LDA(lattice, atoms, positions; magnetic_moments, temperature=0.01)
basis = PlaneWaveBasis(model; Ecut=10, kgrid=(2, 2, 2))
ρ0 = guess_density(basis, magnetic_moments)
scfres = self_consistent_field(basis, tol=1e-4, ρ=ρ0, mixing=KerkerMixing());</code></pre><pre class="documenter-example-output">n Energy log10(ΔE) log10(Δρ) Magnet Diag
--- --------------- --------- --------- ------ ----
1 -223.7481691487 0.22 -6.321 4.7
2 -224.1772742481 -0.37 -0.23 -3.269 1.6
3 -224.2162796534 -1.41 -1.07 -1.801 2.7
4 -224.2192878595 -2.52 -1.36 -1.353 2.4
5 -224.2206002824 -2.88 -1.68 -0.898 1.2
6 -224.2205609626 + -4.41 -1.78 -0.634 1.0</pre><div class="admonition is-info"><header class="admonition-header">DFTK and ASE</header><div class="admonition-body"><p>DFTK also supports using ASE to setup a calculation and provides two-way conversion routines between key DFTK and ASE datastructures. See <a href="../../examples/ase/#Creating-slabs-with-ASE">Creating slabs with ASE</a> for details.</p></div></div><h2 id="Writing-VTK-files-for-visualization"><a class="docs-heading-anchor" href="#Writing-VTK-files-for-visualization">Writing VTK files for visualization</a><a id="Writing-VTK-files-for-visualization-1"></a><a class="docs-heading-anchor-permalink" href="#Writing-VTK-files-for-visualization" title="Permalink"></a></h2><p>For visualizing the density or the Kohn-Sham orbitals DFTK supports storing the result of an SCF calculations in the form of VTK files. These can afterwards be visualized using tools such as <a href="https://www.paraview.org/">paraview</a>. Using this feature requires the <a href="https://github.com/jipolanco/WriteVTK.jl/">WriteVTK.jl</a> Julia package.</p><pre><code class="language-julia">using WriteVTK
save_scfres("iron_afm.vts", scfres; save_ψ=true);</code></pre><p>This will save the iron calculation above into the file <code>iron_afm.vts</code>, using <code>save_ψ=true</code> to also include the KS orbitals.</p><h2 id="Writing-and-reading-JLD2-files"><a class="docs-heading-anchor" href="#Writing-and-reading-JLD2-files">Writing and reading JLD2 files</a><a id="Writing-and-reading-JLD2-files-1"></a><a class="docs-heading-anchor-permalink" href="#Writing-and-reading-JLD2-files" title="Permalink"></a></h2><p>The full state of a DFTK self-consistent field calculation can be stored on disk in form of an <a href="https://github.com/JuliaIO/JLD2.jl">JLD2.jl</a> file. This file can be read from other Julia scripts as well as other external codes supporting the HDF5 file format (since the JLD2 format is based on HDF5).</p><pre><code class="language-julia">using JLD2
save_scfres("iron_afm.jld2", scfres);</code></pre><p>Since such JLD2 can also be read by DFTK to start or continue a calculation, these can also be used for checkpointing or for transferring results to a different computer. See <a href="../../examples/scf_checkpoints/#Saving-SCF-results-on-disk-and-SCF-checkpoints">Saving SCF results on disk and SCF checkpoints</a> for details.</p><p>(Cleanup files generated by this notebook)</p><pre><code class="language-julia">rm("iron_afm.vts")
rm("iron_afm.jld2")</code></pre></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="../tutorial/">« Tutorial</a><a class="docs-footer-nextpage" href="../parallelization/">Timings and parallelization »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Sunday 20 March 2022 18:55">Sunday 20 March 2022</span>. Using Julia version 1.7.2.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>