[CP2K-user] [CP2K:20404] Polarization calculation

Sue Park sueparkchemist at gmail.com
Wed Jul 3 04:50:58 UTC 2024


Dear CP2K community,

I'm trying to learn how to calculate polarization using CP2K. I'm using 
2x2x2 supercell of BaTiO3 to benchmark the polarization values compared to 
other references and Quantum Espresso. 

I'm using Berry phase method by setting *PERIODIC .TRUE. in MOMENTS* as 
it's periodic system. And the output for the cubic structure looks like:
  Dipole vectors are based on the periodic (Berry phase) operator.
  They are defined modulo integer multiples of the cell matrix [Debye].
  [X] [   38.42563418     0.00000000     0.00000000 ] [i]
  [Y]=[    0.00000000    38.42563418     0.00000000 ]*[j]
  [Z] [    0.00000000     0.00000000    38.65618798 ] [k]
  Dipole moment [Debye]
    X= -0.1123396E-09 Y= -0.1123351E-09 Z= -0.1414722E-09  Total=     
 0.0000000

How I have been calculating the polarization value is: 1) As I'm only 
interested in the polarization in z direction, I combined the two values 
highlighted in green. And 2) I divide the values by the cell volume. 3) I 
convert the unit to uC/cm^2. The issue with this method is that the 
polarization of the cubic phase is not zero. I've been trying to figure 
this out because the polarization value for cubic phase should be zero or 
close to zero. And the zero polarization for cubic phase is reproduced by 
Quantum Espresso:
                             VALUES OF POLARIZATION
                             ~~~~~~~~~~~~~~~~~~~~~~
        The calculation of phases done along the direction of vector 3
        of the reciprocal lattice gives the following contribution to
        the polarization vector (in different units, and being Omega
        the volume of the unit cell):
           P =   0.0000000  (mod  30.0953409)  (e/Omega).bohr
           P =   0.0000000  (mod   0.0088327)  e/bohr^2
           P =   0.0000000  (mod   0.5053585)  C/m^2
        The polarization direction is:  ( 0.00000 , 0.00000 , 1.00000 )


What I'm confused with is that the difference between the values in the 
dipole vectors and the dipole moment. Is either of them electronic 
polarization (from Berry phase) and ionic polarization? Or is the dipole 
moment in the bottom line calculated from the addition of the electronic 
polarization (Berry phase) and the ionic polarization? Is there any 
reference/literature that describes how CP2K operates the Berry phase 
calculation?


Here's my input script:
&GLOBAL
  PROJECT BTO
  RUN_TYPE CELL_OPT
  PRINT_LEVEL MEDIUM
&END GLOBAL

&FORCE_EVAL
  ! This one is used for cell optimization
  STRESS_TENSOR ANALYTICAL
  METHOD  QS
  &SUBSYS
    &KIND Ba
      ELEMENT   Ba
      BASIS_SET DZVP-MOLOPT-PBE-GTH-q10
      POTENTIAL GTH-PBE-q10
    &END KIND
    &KIND Ti
      ELEMENT   Ti
      BASIS_SET DZVP-MOLOPT-PBE-GTH-q12
      POTENTIAL GTH-PBE-q12
    &END KIND
    &KIND O
      ELEMENT   O
      BASIS_SET DZVP-MOLOPT-PBE-GTH-q6
      POTENTIAL GTH-PBE-q6
    &END KIND
    &CELL
      A     8.000000000    0.000000000    0.000000000
      B     0.000000000    8.000000000    0.000000000
      C     0.000000000    0.000000000    8.000000000
      PERIODIC XYZ
    &END CELL
    &TOPOLOGY
      COORD_FILE_FORMAT pdb
      COORD_FILE_NAME BTO.pdb
    &END
  &END SUBSYS
  &DFT
    BASIS_SET_FILE_NAME  BASIS_MOLOPT_UZH
    POTENTIAL_FILE_NAME  GTH_POTENTIALS
    &QS
      METHOD GPW
      EPS_DEFAULT 1.0E-12
    &END QS
    &POISSON
      PERIODIC XYZ
      PSOLVER PERIODIC
    &END
    &XC
      &XC_FUNCTIONAL
         &PBE
         &END PBE
      &END
      &VDW_POTENTIAL
         POTENTIAL_TYPE PAIR_POTENTIAL
         &PAIR_POTENTIAL
            PARAMETER_FILE_NAME ../../dftd3.dat
            TYPE DFTD3(BJ)
            REFERENCE_FUNCTIONAL PBE
         &END
      &END VDW_POTENTIAL
    &END XC
    &MGRID
      NGRIDS 4
      CUTOFF 600
    &END MGRID
    &SCF
      SCF_GUESS ATOMIC
      EPS_SCF 1.0E-8
      MAX_SCF 300
      &MIXING
         METHOD BROYDEN_MIXING
         ALPHA 0.2
         BETA 1.5
         NBROYDEN 8
      &END MIXING
    &END SCF
    &PRINT
      &MO_CUBES
        WRITE_CUBE .FALSE.
        NHOMO 1
        NLUMO 1
      &END MO_CUBES
      &MOMENTS MEDIUM
        PERIODIC T
      &END MOMENTS
    &END PRINT
  &END DFT
  &PRINT
    &FORCES ON
    &END FORCES
  &END PRINT
&END FORCE_EVAL

&MOTION
  &CELL_OPT
    TYPE DIRECT_CELL_OPT
    MAX_ITER 300
    OPTIMIZER BFGS    #Most efficient minimizer, but only for 'small' 
systems
  &END CELL_OPT
  &GEO_OPT
    MAX_DR 0.003
    MAX_FORCE 0.00045
    RMS_DR 0.0015
    RMS_FORCE 0.0003
    MAX_ITER 200
    OPTIMIZER BFGS
    &BFGS
      TRUST_RADIUS 0.25
    &END BFGS
  &END GEO_OPT
&END MOTION

Here's the input structure:
CRYST1    8.000    8.000    8.000  90.00  90.00  90.00 P 1           1
ATOM      1  Ba  BTO A   1       0.000   0.000   0.000  1.00  0.00         
 Ba
ATOM      2  Ti  BTO A   1       2.000   2.000   2.000  1.00  0.00         
 Ti
ATOM      3  O   BTO A   1       0.000   2.000   2.000  1.00  0.00         
  O
ATOM      4  O   BTO A   1       2.000   0.000   2.000  1.00  0.00         
  O
ATOM      5  O   BTO A   1       2.000   2.000   0.000  1.00  0.00         
  O
ATOM      6  Ba  BTO A   1       0.000   0.000   4.000  1.00  0.00         
 Ba
ATOM      7  Ti  BTO A   1       2.000   2.000   6.000  1.00  0.00         
 Ti
ATOM      8  O   BTO A   1       0.000   2.000   6.000  1.00  0.00         
  O
ATOM      9  O   BTO A   1       2.000   0.000   6.000  1.00  0.00         
  O
ATOM     10  O   BTO A   1       2.000   2.000   4.000  1.00  0.00         
  O
ATOM     11  Ba  BTO A   1       0.000   4.000   0.000  1.00  0.00         
 Ba
ATOM     12  Ti  BTO A   1       2.000   6.000   2.000  1.00  0.00         
 Ti
ATOM     13  O   BTO A   1       0.000   6.000   2.000  1.00  0.00         
  O
ATOM     14  O   BTO A   1       2.000   4.000   2.000  1.00  0.00         
  O
ATOM     15  O   BTO A   1       2.000   6.000   0.000  1.00  0.00         
  O
ATOM     16  Ba  BTO A   1       0.000   4.000   4.000  1.00  0.00         
 Ba
ATOM     17  Ti  BTO A   1       2.000   6.000   6.000  1.00  0.00         
 Ti
ATOM     18  O   BTO A   1       0.000   6.000   6.000  1.00  0.00         
  O
ATOM     19  O   BTO A   1       2.000   4.000   6.000  1.00  0.00         
  O
ATOM     20  O   BTO A   1       2.000   6.000   4.000  1.00  0.00         
  O
ATOM     21  Ba  BTO A   1       4.000   0.000   0.000  1.00  0.00         
 Ba
ATOM     22  Ti  BTO A   1       6.000   2.000   2.000  1.00  0.00         
 Ti
ATOM     23  O   BTO A   1       4.000   2.000   2.000  1.00  0.00         
  O
ATOM     24  O   BTO A   1       6.000   0.000   2.000  1.00  0.00         
  O
ATOM     25  O   BTO A   1       6.000   2.000   0.000  1.00  0.00         
  O
ATOM     26  Ba  BTO A   1       4.000   0.000   4.000  1.00  0.00         
 Ba
ATOM     27  Ti  BTO A   1       6.000   2.000   6.000  1.00  0.00         
 Ti
ATOM     28  O   BTO A   1       4.000   2.000   6.000  1.00  0.00         
  O
ATOM     29  O   BTO A   1       6.000   0.000   6.000  1.00  0.00         
  O
ATOM     30  O   BTO A   1       6.000   2.000   4.000  1.00  0.00         
  O
ATOM     31  Ba  BTO A   1       4.000   4.000   0.000  1.00  0.00         
 Ba
ATOM     32  Ti  BTO A   1       6.000   6.000   2.000  1.00  0.00         
 Ti
ATOM     33  O   BTO A   1       4.000   6.000   2.000  1.00  0.00         
  O
ATOM     34  O   BTO A   1       6.000   4.000   2.000  1.00  0.00         
  O
ATOM     35  O   BTO A   1       6.000   6.000   0.000  1.00  0.00         
  O
ATOM     36  Ba  BTO A   1       4.000   4.000   4.000  1.00  0.00         
 Ba
ATOM     37  Ti  BTO A   1       6.000   6.000   6.000  1.00  0.00         
 Ti
ATOM     38  O   BTO A   1       4.000   6.000   6.000  1.00  0.00         
  O
ATOM     39  O   BTO A   1       6.000   4.000   6.000  1.00  0.00         
  O
ATOM     40  O   BTO A   1       6.000   6.000   4.000  1.00  0.00         
  O
END


Any of your inputs would be greatly appreciated. Thank you for your time 
and consideration in advance.

Best, 
Sue

-- 
You received this message because you are subscribed to the Google Groups "cp2k" group.
To unsubscribe from this group and stop receiving emails from it, send an email to cp2k+unsubscribe at googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/cp2k/d69b786c-2b3e-47be-bada-caef6b123f7dn%40googlegroups.com.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20240702/38792979/attachment-0001.htm>


More information about the CP2K-user mailing list