【CP2K教程(二)】WO3的投影态密度和能带结构

1. WO3的投影态密度和能带结构

2. 获得WO3晶格的能带结构

1. WO3的投影态密度和能带结构

在本练习中,您将使用立方晶格WO3的K点采样进行态密度(DOS)和能带结构计算。在本文中可以找到参考DOS和能带结构。

获取PDOS:

在下面的练习中,我们将研究WO3的态密度:

与上一个练习类似,我们根据unit cell编写坐标:

WO3-PDOS.inp

&GLOBAL
   PROJECT WO3-pdos
   RUN_TYPE ENERGY
   PRINT_LEVEL MEDIUM
&END GLOBAL

&FORCE_EVAL
   METHOD Quickstep
   &DFT
      BASIS_SET_FILE_NAME  BASIS_MOLOPT
      POTENTIAL_FILE_NAME  POTENTIAL

      &POISSON
         PERIODIC XYZ
      &END POISSON

      &SCF
         SCF_GUESS ATOMIC
         EPS_SCF 1.0E-6
         MAX_SCF 300
         ADDED_MOS 100
         &DIAGONALIZATION
            ALGORITHM STANDARD
            EPS_ADAPT 0.01
         &END DIAGONALIZATION
         &SMEAR  ON
            METHOD FERMI_DIRAC
            ELECTRONIC_TEMPERATURE [K] 300
         &END SMEAR

         &MIXING
            METHOD BROYDEN_MIXING
            ALPHA 0.2
            BETA 1.5
            NBROYDEN 8
         &END MIXING

      &END SCF
      &XC
         &XC_FUNCTIONAL PBE
         &END XC_FUNCTIONAL
      &END XC
      &PRINT
         &PDOS
            # print all projected DOS available:
            NLUMO -1
            # split the density by quantum number:
            COMPONENTS
         &END
      &END PRINT
   &END DFT

   &SUBSYS
      &CELL
         ABC [angstrom] 3.810000 3.810000 3.810000
         PERIODIC XYZ
         MULTIPLE_UNIT_CELL 2 2 2
      &END CELL
      &TOPOLOGY
         MULTIPLE_UNIT_CELL 2 2 2 
      &END TOPOLOGY
      &COORD
         SCALED
         W 0.0 0.0 0.0
         O 0.5 0.0 0.0
         O 0.0 0.5 0.0
         O 0.0 0.0 0.5
      &END
      &KIND W
         ELEMENT W
         BASIS_SET DZVP-MOLOPT-SR-GTH
         POTENTIAL GTH-PBE
      &END KIND
      &KIND O
         ELEMENT O
         BASIS_SET DZVP-MOLOPT-SR-GTH
         POTENTIAL GTH-PBE
      &END KIND

   &END SUBSYS

&END FORCE_EVAL

unit cell的扩增是必要的,因为除非另有指示,否则程序仅在Γ点采样,我们难以获得态密度的有意义采样(例如,布里渊区上的网格将过于粗糙)。另一个选项(我们将在下一个练习中研究)是在k点上采样。

除了输出文件之外,您还将得到一个名为WO3_pdos-k1-1.pdos的文件(确切地说,每个原子类型您将获得一个这样的文件,但这里我们只有一个,碳),其内容类似于:

# Projected DOS for atomic kind W at iteration step i = 0, E(Fermi) =     0.144475 a.u.
#     MO Eigenvalue [a.u.]      Occupation                 s                py                pz                px               d-2               d-1                d0               d+1               d+2               f-3               f-2               f-1                f0               f+1               f+2               f+3
       1         -2.621088        2.000000        0.87115225        0.00000000        0.00000000        0.00000000        0.00000000        0.00000000        0.00000256        0.00000000        0.00000256        0.00000000        0.00000000        0.00000000        0.00000000        0.00000000        0.00000000        0.00000000
       2         -2.621080        2.000000        0.85006340        0.00131878        0.00166813        0.00134861        0.00000000        0.00000000        0.00515023        0.00000000        0.00441289        0.00006412        0.00000000        0.00003847        0.00012977        0.00003934        0.00000000        0.00006557

[...]

这些列对应于基组中存在的轨道(因此是投影DOS)。现在,您可以使用高斯曲线绘制卷积图,以获得平滑DOS

对于复杂的DOS,您可能需要查看此网站。这里提供了两个Python脚本来进行卷积。下载 pdos.py 和 get-smearing-pdos.py两个文件到您的文件夹。并使用以下命令执行python脚本

python get-smearing-pdos.py file.pdos

 pdos.py

#! /usr/bin/env python

from math import pi, sqrt	
import numpy as np

class pdos:
    """ Projected electronic density of states from CP2K output files

        Attributes
        ----------
        atom: str 
            the name of the atom where the DoS is projected
        iterstep: int
            the iteration step from the CP2K job
        efermi: float
            the energy of the Fermi level [a.u]
        e: float
            (eigenvalue - efermi) in eV
        occupation: int
            1 for occupied state or 0 for unoccupied
        pdos: nested list of float
            projected density of states on each orbital for each eigenvalue
            [[s1, p1, d1,....], [s2, p2, d2,...],...]
            s: pdos in s orbitals
            p: pdos in p orbitals
            d: pdos in d orbitals
            .
            .
            .
        tpdos: list of float
            sum of all the orbitals PDOS
            
        Methods
        -------
        smearing(self,npts, width)
            return the smeared tpdos 
    """
    
    def __init__(self, infilename):
        """Read a CP2K .pdos file and build a pdos instance

        Parameters
        ----------
        infilename: str
            pdos output from CP2K. 

        """    
        input_file = open(infilename, 'r')

        firstline  = input_file.readline().strip().split()
        secondline = input_file.readline().strip().split()

        # Kind of atom
        self.atom = firstline[6]
        #iterationstep
        self.iterstep = int(firstline[12][:-1]) #[:-1] delete ","
        # Energy of the Fermi level
        self.efermi = float(firstline[15])

        # it keeps just the orbital names
        secondline[0:5] = []
        self.orbitals = secondline 

        lines = input_file.readlines()
   
        eigenvalue = []
        self.occupation = []
        data = []
        self.pdos = []
        for index, line in enumerate(lines):
            data.append(line.split())
            data[index].pop(0)
            eigenvalue.append(float(data[index].pop(0)))
            self.occupation.append(int(float(data[index].pop(0))))
            self.pdos.append([float(i) for i in data[index]])

        self.e = [(x-self.efermi)*27.211384523 for x in eigenvalue] 

        self.tpdos = []
        for i in self.pdos:
           self.tpdos.append(sum(i))

    def __add__(self, other):
        """Return the sum of two PDOS objects"""
        sumtpdos = [i+j for i,j in zip(self.tpdos,other.tpdos)]
        return sumtpdos

    def delta(self,emin,emax,npts,energy,width):
        """Return a delta-function centered at energy
        
        Parameters
        ----------
        emin: float
            minimun eigenvalue
        emax: float
            maximun eigenvalue
        npts: int
            Number of points in the smeared pdos
        energy: float
            energy where the gaussian is centered
        width: float
            dispersion parameter

        Return 
        ------
        delta: numpy array
            array of delta function values

        """
        
        energies = np.linspace(emin, emax, npts)
        x = -((energies - energy) / width)**2
        return np.exp(x) / (sqrt(pi) * width)

    def smearing(self,npts, width,):
        """Return a gaussian smeared DOS"""

        d = np.zeros(npts)
        emin = min(self.e)
        emax = max(self.e)
        for e, pd in zip(self.e,self.tpdos):
            d +=pd*self.delta(emin,emax,npts,e,width)

        return d
    
def sum_tpdos(tpdos1, tpdos2):
    """Return the sum of two PDOS"""
    return [i+j for i,j in zip(tpdos1,tpdos2)]
            

get-smearing-pdos.py

#! /usr/bin/env python
#---------------------------------------------------
# get-smearing-pdos.py: read one or a pair alpha,  
# beta spin files with the cp2k pdos format and
# return a file "smeared.dat" with the Smeared DOS
#---------------------------------------------------
# Usage: ./get-smearing-pdos.py ALPHA.pdos BETA.pdos
#        or
#        ./get-smearing-pdos.py file.pdos 
#
# Output: 
#         smeared.dat: smeared DOS
#---------------------------------------------------
# Todo:
# - Atomatic name generation of output file
# - Move the algorithm to the module pdos
# - Implement printing of d orbitals
# - ...
#---------------------------------------------------
# Author: Juan Garcia e-mail: jcgarcia [at] wpi.edu
# Date:   11-12-2012
#---------------------------------------------------

import sys
from pdos import *


if len(sys.argv) == 2:

    infilename = sys.argv[1]

    alpha = pdos(infilename)
    npts = len(alpha.e)
    alpha_smeared = alpha.smearing(npts,0.2)
    eigenvalues = np.linspace(min(alpha.e), max(alpha.e),npts)
    
    g = open('smeared.dat','w')
    for i,j in zip(eigenvalues, alpha_smeared):
        t = str(i).ljust(15) + '     ' + str(j).ljust(15) + '\n'
        g.write(t)

elif len(sys.argv) == 3:

    infilename1 = sys.argv[1]
    infilename2 = sys.argv[2]

    alpha = pdos(infilename1)
    beta = pdos(infilename2)
    npts = len(alpha.e)
    alpha_smeared = alpha.smearing(npts,0.2)
    beta_smeared = beta.smearing(npts,0.2)
    totalDOS = sum_tpdos(alpha_smeared, beta_smeared)
    
    eigenvalues = np.linspace(min(alpha.e), max(alpha.e),npts)

    g = open('smeared.dat','w')
    for i,j in zip(eigenvalues, totalDOS):
        t = str(i).ljust(15) + '     ' + str(j).ljust(15) + '\n'
        g.write(t)

else:
    print '  Wrong number of arguments!'
    print '  usage:'
    print '  ./get-smearing-pdos.py ALPHA.pdos'
    print '  ./get-smearing-pdos.py ALPHA.pdos BETA.pdos'





或者,您也可以使用Tiziano开发的Python脚本,解析器错误已经修复。

cp2k_pdos.py

#!/usr/bin/env python
# vim: set fileencoding=utf-8 ts=8 sw=4 tw=0 :

"""
Convert the discrete CP2K PDOS points to a smoothed curve using convoluted gaussians.

Also shifts the energies by the Fermi energy (so the Fermi energy will afterwards be at 0),
and normalizes by the number of atoms of this kind.
"""

# Copyright (c) 2017 Tiziano Müller <tiziano.mueller@chem.uzh.ch>,
# based on a Fortran tool written by Marcella Iannuzzi <marcella.iannuzzi@chem.uzh.ch>
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.


from __future__ import print_function

import sys
import re
import argparse
import contextlib

import numpy as np


HEADER_MATCH = re.compile(
    r'\# Projected DOS for atomic kind (?P<element>\w+) at iteration step i = \d+, E\(Fermi\) = [ \t]* (?P<Efermi>[^\t ]+) a\.u\.')

# Column indexes, starting from 0
EIGENVALUE_COLUMN = 1
DENSITY_COLUMN = 3


# from https://stackoverflow.com/a/17603000/1400465
@contextlib.contextmanager
def smart_open(filename=None):
    if filename and filename != '-':
        fhandle = open(filename, 'w')
    else:
        fhandle = sys.stdout

    try:
        yield fhandle
    finally:
        if fhandle is not sys.stdout:
            fhandle.close()


def main():
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument('pdosfilenames', metavar='<PDOS-file>', type=str, nargs='+',
                        help="the PDOS file generated by CP2K, specify two (alpha/beta) files for a common grid")
    parser.add_argument('--sigma', '-s', type=float, default=0.02,
                        help="sigma for the gaussian distribution (default: 0.02)")
    parser.add_argument('--de', '-d', type=float, default=0.001,
                        help="integration step size (default: 0.001)")
    parser.add_argument('--scale', '-c', type=float, default=1,
                        help="scale the density by this factor (default: 1)")
    parser.add_argument('--total-sum', action='store_true',
                        help="calculate and print the total sum for each orbital (default: no)")
    parser.add_argument('--no-header', action='store_true', default=False,
                        help="do not print a header (default: print header)")
    parser.add_argument('--output', '-o', type=str, default=None,
                        help="write output to specified file (default: write to standard output)")
    args = parser.parse_args()

    alldata = []
    orb_headers = []

    for pdosfilename in args.pdosfilenames:
        with open(pdosfilename, 'r') as fhandle:
            match = HEADER_MATCH.match(fhandle.readline().rstrip())
            if not match:
                print(("The file '{}' does not look like a CP2K PDOS output.\n"
                       "If it is indeed a correct output file, please report an issue at\n"
                       "    https://github.com/dev-zero/cp2k-tools/issues").format(pdosfilename))
                sys.exit(1)

            efermi = float(match.group('Efermi'))
            # header is originally: ['#', 'MO', 'Eigenvalue', '[a.u.]', 'Occupation', 's', 'py', ...]
            header = fhandle.readline().rstrip().split()[1:]  # remove the comment directly
            header[1:3] = [' '.join(header[1:3])]  # rejoin "Eigenvalue" and its unit
            data = np.loadtxt(fhandle)  # load the rest directly with numpy

        alldata.append(data)

        orb_headers += header[DENSITY_COLUMN:]

    # take the boundaries over all energy eigenvalues (not guaranteed to be the same)
    # add a margin to not cut-off Gaussians at the borders
    margin = 10 * args.sigma
    emin = min(np.min(data[:, EIGENVALUE_COLUMN]) for data in alldata) - margin
    emax = max(np.max(data[:, EIGENVALUE_COLUMN]) for data in alldata) + margin
    ncols = sum(data.shape[1] - DENSITY_COLUMN for data in alldata)
    nmesh = int((emax-emin)/args.de)+1  # calculate manually instead of using np.arange to ensure emax inside the mesh
    xmesh = np.linspace(emin, emax, nmesh)
    ymesh = np.zeros((nmesh, ncols))

    # printing to stderr makes it possible to simply redirect the stdout to a file
    print("Integration step:  {:14.5f}".format(args.de), file=sys.stderr)
    print("Sigma:             {:14.5f}".format(args.sigma), file=sys.stderr)
    print("Minimum energy:    {:14.5f} - {:.5f}".format(emin+margin, margin), file=sys.stderr)
    print("Maximum energy:    {:14.5f} + {:.5f}".format(emax-margin, margin), file=sys.stderr)
    print("Nr of mesh points: {:14d}".format(nmesh), file=sys.stderr)

    fact = args.de/(args.sigma*np.sqrt(2.0*np.pi))

    coloffset = 0
    for fname, data in zip(args.pdosfilenames, alldata):
        print("Nr of lines:       {:14d} in {}".format(data.shape[0], fname), file=sys.stderr)
        ncol = data.shape[1] - DENSITY_COLUMN

        for idx in range(nmesh):
            func = np.exp(-(xmesh[idx]-data[:, EIGENVALUE_COLUMN])**2/(2.0*args.sigma**2))*fact
            ymesh[idx, coloffset:(coloffset+ncol)] = func.dot(data[:, DENSITY_COLUMN:])

        coloffset += ncol

    if args.total_sum:
        finalsum = np.sum(ymesh, 0)*args.de
        print("Sum over all meshpoints, per orbital:", file=sys.stderr)
        print(("{:16.8f}"*ncols).format(*finalsum), file=sys.stderr)

    xmesh -= efermi  # put the Fermi energy at 0
    xmesh *= 27.211384  # convert to eV
    ymesh *= args.scale  # scale

    with smart_open(args.output) as fhandle:
        if not args.no_header:
            print(("{:>16}" + " {:>16}"*ncols).format("Energy_[eV]", *orb_headers), file=fhandle)
        for idx in range(nmesh):
            print(("{:16.8f}" + " {:16.8f}"*ncols).format(xmesh[idx], *ymesh[idx, :]), file=fhandle)


if __name__ == '__main__':
    main()

#  vim: set ts=4 sw=4 tw=0 :

不同的σ值会产生不同的卷积,这意味着线型不同。需要合理的σ值才能获得良好的PDOS图。当可视化PDO时,只有接近费米能级的能量区域是有趣的。需要适当调整xrange。

还请注意能量单位,单位为Eh。在查看DOS图时,您可能希望将其转换为Electronvolt。在卷积程序中,这已在代码中完成。虽然有助于收敛的一些新选项是数值性质的,但smearing不是。

  • 对不同的多个多胞3x3x3、4x4x4重复上述计算

  • 获得O2p和W5d轨道的总DOS和PDOS,并与文献值进行比较。

  • 你知道为什么有必要进行单胞扩增吗?

  • WO3带隙的值是多少?比较3x3x3和4x4x4的图。

  • .. 哪个状态(ss,pxpx,…)是主要原因?

  • 改变卷积程序中的σ值,并确定PDO图的合理值

2. 获得WO3晶格的能带结构

为了获得WO3的能带结构,与用于计算PDOS的前一示例相比,只需要一些改变:

WO3-bs.inp

&GLOBAL
   PROJECT WO3-kp-bs
   RUN_TYPE ENERGY
   PRINT_LEVEL MEDIUM
&END GLOBAL

&FORCE_EVAL
   METHOD Quickstep
   &DFT
      BASIS_SET_FILE_NAME  BASIS_MOLOPT
      POTENTIAL_FILE_NAME  POTENTIAL

      &POISSON
         PERIODIC XYZ
      &END POISSON
      &QS
         EXTRAPOLATION USE_GUESS ! required for K-Point sampling
      &END QS
      &SCF
         SCF_GUESS ATOMIC
         EPS_SCF 1.0E-6
         MAX_SCF 300

         ADDED_MOS 2
         &DIAGONALIZATION
            ALGORITHM STANDARD
            EPS_ADAPT 0.01
         &END DIAGONALIZATION
         &SMEAR  ON
            METHOD FERMI_DIRAC
            ELECTRONIC_TEMPERATURE [K] 300
         &END SMEAR

         &MIXING
            METHOD BROYDEN_MIXING
            ALPHA 0.2
            BETA 1.5
            NBROYDEN 8
         &END MIXING

      &END SCF
      &XC
         &XC_FUNCTIONAL PBE
         &END XC_FUNCTIONAL
      &END XC
      &KPOINTS
         SCHEME MONKHORST-PACK 3 3 1
         WAVEFUNCTIONS REAL
         SYMMETRY .FALSE.
         FULL_GRID .FALSE.
         PARALLEL_GROUP_SIZE -1
      &END KPOINTS
      &PRINT
         &BAND_STRUCTURE
            ADDED_MOS 2
            FILE_NAME WO3.bs
            &KPOINT_SET
               UNITS B_VECTOR
               SPECIAL_POINT ???   #GAMA
               SPECIAL_POINT ???   #X
               SPECIAL_POINT ???   #M
               SPECIAL_POINT ???   #GAMA
               SPECIAL_POINT ???   #R
               SPECIAL_POINT ???   #M
               NPOINTS ???
            &END
         &END BAND_STRUCTURE
      &END PRINT
   &END DFT

   &SUBSYS
      &CELL
         ABC [angstrom] 3.810000 3.810000 3.810000
         PERIODIC XYZ
         MULTIPLE_UNIT_CELL 1 1 1
      &END CELL
      &TOPOLOGY
         MULTIPLE_UNIT_CELL 1 1 1
      &END TOPOLOGY
      &COORD
         SCALED
         W 0.0 0.0 0.0
         O 0.5 0.0 0.0
         O 0.0 0.5 0.0
         O 0.0 0.0 0.5
      &END
      &KIND W
         ELEMENT W
         BASIS_SET DZVP-MOLOPT-SR-GTH
         POTENTIAL GTH-PBE
      &END KIND
      &KIND O
         ELEMENT O
         BASIS_SET DZVP-MOLOPT-SR-GTH
         POTENTIAL GTH-PBE
      &END KIND
   &END SUBSYS

&END FORCE_EVAL

目前,在进行K点计算时,不可能获得投影态密度。特殊点应根据b-vectors给出。

关于输入文件的一些注释:

  • 通过指定KPOINT小节,可以启用K点计算。

  • 虽然您可以直接指定K点,但我们使用Monkhorst包方案来生成它们。关键字MONKHORST-PACK后面的数字指定布里渊区的tiling。

  • ​在基本计算之后,CP2K计算沿某些线的能量,表示为KPOINT_SET(当您检查文档时,您会注意到本节可以被重复)。

  • 关键字NPOINTS指定在两个特殊点之间应采样多少点(除了起点之外)。

  • SPECIAL_POINT关键字用于指定直线的起点、中点和终点。这些点通常表示倒晶格/晶胞中的特殊点,如Γ、M或K点。您可以在本文的附录部分找到这些的定义。也可以多次指定该关键字,从而可以直接获得完整路径的带结构。

当通过K点进行采样时,建议您使用搜索路径工具,以获得具有倒数空间中重要路径的CP2K骨架输入文件。将以下搜索路径作为XYZ文件,并指定一个简单的立方体单元,其晶格常数a如下所示:

WO3-cubic.xyz

4
WO3; a=3.810000
 W    0.000000    0.000000    0.000000
 O    1.905000    0.000000    0.000000
 O    0.000000    1.905000    0.000000
 O    0.000000    0.000000    1.905000

现在,当您运行这个输入文件时,除了输出文件之外,您还将得到一个名为WO3.bs的文件,该文件看起来类似于:

 SET:       1                 TOTAL POINTS:      26
   POINT   1                     ********    ********    ********
   POINT   2                     ********    ********    ********
   POINT   3                     ********    ********    ********
   POINT   4                     ********    ********    ********
   POINT   5                     ********    ********    ********
   POINT   6                     ********    ********    ********
       Nr.    1    Spin 1        K-Point  0.00000000  0.00000000  0.00000000
               20
           -73.66652408    -38.53370023    -37.80464132    -37.79327769
           -16.71308703    -16.11075946    -16.02553853     -1.43495530
            -1.34739188     -1.33357408      0.37912017      0.38948689
             0.39582882      0.40030859      0.46965212      0.47418816
             2.60728842      2.62105342      3.16044140      6.99806305
       Nr.    2    Spin 1        K-Point  0.00000000  0.10000000  0.00000000
               20
           -73.66647294    -38.53337818    -37.80859042    -37.79536623
           -16.67479677    -16.09554462    -15.96731960     -1.68492873
            -1.44087258     -1.34318045      0.09257368      0.13769271
             0.21643888      0.38447849      0.44179002      0.45757924
             2.61768501      3.02368022      3.51828287      7.06644645

[...]

对于每个集,有一个名为SET的块,其中特殊点列为POINT,然后是每个K点的子块,包含每个MO的能量。

Your tasks:

  • ​查找上述论文中的Γ、X、M、R点的特殊点(确保选择正确的晶格)。计算并绘制WO3晶格沿Γ、X、M、Γ、R、M的能带结构(您可以自由决定是否使用多个K点集或单个集合中的多个特殊点)。标记特殊点。选择适当数量的插值点以获得平滑绘图。

  • 将你的绘图与文献中的绘图进行比较。有什么不同?

  • 你得到了多少轨道能量,为什么?尝试更改输入以获得更多未占用的轨道。

要将能带结构文件转换为可以直接绘制的文件,您可以使用下面的脚本 cp2k_bs2csv.py,当将能带结构文件 WO3.bs 作为参数传递时,该脚本将写入文件 WO3.bs-set-1。 csv 对于每组包含 K 点坐标和能量在一行中。

要绘制 WO3.bs-set-1.csv 文件,您可以将其加载到 MATLAB 或使用 GNUPLOTGNUPLOT 命令行。

gnuplot>set yrange [-8:14]
gnuplot>plot for [i=4:23] "WO3.bs.set-1.csv" u 0:i w l t ""

cp2k_bs2csv.py

#!/usr/bin/env python
"""
Convert the CP2K band structure output to CSV files
"""
 
import re
import argparse
 
SET_MATCH = re.compile(r'''
[ ]*
  SET: [ ]* (?P<setnr>\d+) [ ]*
  TOTAL [ ] POINTS: [ ]* (?P<totalpoints>\d+) [ ]*
  \n
(?P<content>
  [\s\S]*?(?=\n.*?[ ] SET|$)  # match everything until next 'SET' or EOL
)
''', re.VERBOSE)
 
SPOINTS_MATCH = re.compile(r'''
[ ]*
  POINT [ ]+ (?P<pointnr>\d+) [ ]+ (?P<a>\S+) [ ]+ (?P<b>\S+) [ ]+ (?P<c>\S+)
''', re.VERBOSE)
 
POINTS_MATCH = re.compile(r'''
[ ]*
  Nr\. [ ]+ (?P<nr>\d+) [ ]+
  Spin [ ]+ (?P<spin>\d+) [ ]+
  K-Point [ ]+ (?P<a>\S+) [ ]+ (?P<b>\S+) [ ]+ (?P<c>\S+) [ ]*
  \n
[ ]* (?P<npoints>\d+) [ ]* \n
(?P<values>
  [\s\S]*?(?=\n.*?[ ] Nr|$)  # match everything until next 'Nr.' or EOL
)
''', re.VERBOSE)
 
if __name__ == '__main__':
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument('bsfilename', metavar='bandstructure-file', type=str,
                        help="the band structure file generated by CP2K")
 
    args = parser.parse_args()
 
    with open(args.bsfilename, 'r') as fhandle:
        for kpoint_set in SET_MATCH.finditer(fhandle.read()):
            filename = "{}.set-{}.csv".format(args.bsfilename,
                                              kpoint_set.group('setnr'))
            set_content = kpoint_set.group('content')
 
            with open(filename, 'w') as csvout:
                print(("writing point set {}"
                       " (total number of k-points: {totalpoints})"
                       .format(filename, **kpoint_set.groupdict())))
 
                print("  with the following special points:")
                for point in SPOINTS_MATCH.finditer(set_content):
                    print("  {pointnr}: {a}/{b}/{c}".format(
                        **point.groupdict()))
 
                for point in POINTS_MATCH.finditer(set_content):
                    results = point.groupdict()
                    results['values'] = " ".join(results['values'].split())
                    csvout.write("{a} {b} {c} {values}\n".format(**results))

参考资料

exercises:2017_uzh_cmest:pdos [CP2K Open Source Molecular Dynamics ]

https://pubs.acs.org/doi/full/10.1021/cm3032225

https://www.sciencedirect.com/science/article/pii/S0927025610002697

Dsp2p Dsp2c是针对个人与企业之间的融资借贷模式推出的免费p2c开源项目,基于Thinkphp架构的行业系统,全自动安装模式,快速搭建P2C网站。 P2C网贷简介 P2C网贷(Person to company)模式是一种个人与企业之间的融资借贷模式,通过线下开发优质的中小企业客户,并引进实力的融资性担保机构对项目进行担保,在线上通过互联网平台寻找普通投资者,是一种安全、平等、高效、透明的互联网金融创新模式,是一种线上线下(O2O)相结合的全新概念。 Dsp2p发展历史: 2014年6月 Dswjcms P2C网贷系统发布 2015年12月 Dswjcms P2C3.2发布,系统架构已同步PHP版本,不支持之前版本直接升级 2016年5月 Dswjcms P2C改名为dsp2c,插件模板市场正在筹建中,模板市场推出后,将代表着dsp2c又一新的起点 2016年6月 Dsp2c 3.3对安装界面进行了优化,后台增加了开发环境一键切换,并对使用dswjcms项目的站点免费提供售后技术支持服务。 Dsp2c使用要求 1、无任何技术要求,我们的项目全是模块化的,初始安装无任何附加插件,所有插件或模板都可以通过官方提供的插件模板市场进行下载,插件模板的安装也都提供教程,只需几步即可,无需修改源代码,也不需要了解项目目录结构;简单的说只要PHP环境已配置好,会上传文件解压文件,会浏览器网站,那么你就可以利用我们的项目搭建属于自己的平台。 2、如你想对我们的项目进行次开发,要求就比较高了:从事PHP研发工作2年以上,至少参与功能的研发与扩展,需熟练掌握ThinkPHP3.1框架的使用。 Dsp2c建站流程 1、下载dsp2c源码 2、上传到服务器,运行域名/install.php 3、一路下一步,填写数据库用户名、密码、表名 4、项目搭建完成 5、通过官方插件模板市场选择喜欢的模板点击下载 6、根据模板安装说明,轻松几步完成模板替换 7、通过官方插件模板市场选择所需的功能插件(如支付接口、短信接口) 8、根据插件安装说明,轻松几步完成插件安装 9、填写插件所需相关配置,完成插件的应用 10、项目对外发布,上线运营。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值