求助：在NVT系综的AIMD计算中出现SCF不收敛的问题

anrushan

在已经完成了在npt的弛豫的构象内正式计算AIMD时出现了scf不收敛的问题，应该不是构象的问题，几个参数也都用的默认数值。尝试增添DIIS项后依然是无法收敛。请各位老师指教。
这次计算分子为一多氮物。
SCF WAVEFUNCTION OPTIMIZATION

  Step     Update method      Time    Convergence         Total energy    Change
  ------------------------------------------------------------------------------
     1 NoMix/Dav.  0.40E+00    0.5     3.65652626      -185.3385542749 -1.85E+02
     2 Broy./Dav.  0.40E+00    0.5     5.04097556      -178.4647054514  6.87E+00
     3 Broy./Dav.  0.40E+00    0.5     0.65877679      -196.9598045613 -1.85E+01
。。。。。。
   125 Broy./Dav.  0.40E+00    0.5     0.02160860      -196.5311668555 -4.15E-02
   126 Broy./Dav.  0.40E+00    0.5     0.01407403      -196.5717251976 -4.06E-02
   127 Broy./Dav.  0.40E+00    0.5     0.02355139      -196.6118522437 -4.01E-02
   128 Broy./Dav.  0.40E+00    0.5     0.00619748      -196.6601445349 -4.83E-02

  Leaving inner SCF loop after reaching   128 steps.


  Electronic density on regular grids:       -107.9999999999        0.0000000001
  Core density on regular grids:              107.9999999950       -0.0000000050
  Total charge density on r-space grids:       -0.0000000049
  Total charge density g-space grids:          -0.0000000049

  Overlap energy of the core charge distribution:               1.18503093169915
  Self energy of the core charge distribution:               -503.47211874021690
  Core Hamiltonian energy:                                    156.76882654704423
  Hartree energy:                                             203.45762976821723
  Exchange-correlation energy:                                -54.59951304167091

  Total energy:                                              -196.66014453492721

*** WARNING in qs_scf.F:576 :: SCF run NOT converged ***
ENERGY| Total FORCE_EVAL ( QS ) energy [a.u.]:             -196.809196806938473

MD_INI| MD initialization
MD_INI| Potential energy [hartree]                          -0.196809196807E+03
MD_INI| Kinetic energy [hartree]                             0.218510258687E+00
MD_INI| Temperature [K]                                             2000.000000
MD_INI| Pressure [bar]                                       8.912970522431E+05
MD_INI| Cell volume [bohr^3]                                 1.858398232120E+03
MD_INI| Cell volume [ang^3]                                  2.753862022160E+02
MD_INI| Cell lengths [bohr]      1.22945582E+01  1.22945582E+01  1.22945582E+01
MD_INI| Cell lengths [ang]       6.50600000E+00  6.50600000E+00  6.50600000E+00
MD_INI| Cell angles [deg]        9.00000000E+01  9.00000000E+01  9.00000000E+01

Number of electrons:                                                        108
Number of occupied orbitals:                                                 54
Number of molecular orbitals:                                                54

Number of orbital functions:                                                296
Number of independent orbital functions:                                    296

Extrapolation method: ASPC

以下是部分输入文件


&GLOBAL
  PROJECT 01
  PRINT_LEVEL LOW
  RUN_TYPE MD
&END GLOBAL

&FORCE_EVAL
  METHOD Quickstep
  &SUBSYS
    &CELL
      A     6.50600000     0.00000000     0.00000000
      B     0.00000000     6.50600000     0.00000000
      C     0.00000000     0.00000000     6.50600000
      PERIODIC XYZ #Direction of applied PBC (geometry aspect)
    &END CELL
    &COORD
  C         5.7512427511       14.6247805356       -3.2379190427
  C         3.5879314640       13.6057951729       -4.1422006332
  N         6.7446849278       14.0749270321       -2.4919660132
  N         7.7454339226       14.9329172329       -2.7696084834
  H         8.6506119723       14.8338896903       -2.3040760447
    &END COORD
#   &VELOCITY #You can set initial atomic velocities in this section
#   &END VELOCITY
    &KIND C
      ELEMENT C
      BASIS_SET DZVP-MOLOPT-SR-GTH
      POTENTIAL GTH-PBE
    &END KIND
    &KIND N
      ELEMENT N
      BASIS_SET DZVP-MOLOPT-SR-GTH
      POTENTIAL GTH-PBE
    &END KIND
    &KIND H
      ELEMENT H
      BASIS_SET DZVP-MOLOPT-SR-GTH
      POTENTIAL GTH-PBE
    &END KIND
  &END SUBSYS

  &DFT
    BASIS_SET_FILE_NAME  BASIS_MOLOPT
    POTENTIAL_FILE_NAME  POTENTIAL
#   WFN_RESTART_FILE_NAME to_be_specified
#Net charge
    MULTIPLICITY  1 #Spin multiplicity
    &QS
      EPS_DEFAULT 1E-10 #This is default. Set all EPS_xxx to values such that the energy will be correct up to this value
      EXTRAPOLATION ASPC #Extrapolation for wavefunction during e.g. MD. ASPC is default, PS also be used
      EXTRAPOLATION_ORDER 3 #Order for PS or ASPC extrapolation. 3 is default
    &END QS
    &POISSON
      PERIODIC XYZ #Direction(s) of PBC for calculating electrostatics
      PSOLVER PERIODIC #The way to solve Poisson equation
    &END POISSON
    &XC
      &XC_FUNCTIONAL PBE
      &END XC_FUNCTIONAL
    &END XC
    &MGRID
      CUTOFF 300
      REL_CUTOFF 40
    &END MGRID
    &SCF
      MAX_SCF 128
      EPS_SCF 1E-5 #Convergence threshold of density matrix during SCF process
#     SCF_GUESS RESTART #Use wavefunction from WFN_RESTART_FILE_NAME file as initial guess
      &DIAGONALIZATION
        ALGORITHM  DAVIDSON
        #Algorithm for diagonalization. DAVIDSON is faster for large systems
      &END DIAGONALIZATION
      &MIXING #How to mix old and new density matrices
        METHOD BROYDEN_MIXING
        ALPHA 0.4 #Default. Mixing 40% of new density matrix with the old one
        NBROYDEN 8 #Default is 4. Number of previous steps stored for the actual mixing scheme
      &END MIXING
      &PRINT
        &RESTART #Use "&RESTART OFF" can prevent generating wfn file
          BACKUP_COPIES 0 #Maximum number of backup copies of wfn file
        &END RESTART
      &END PRINT
    &END SCF
  &END DFT
     STRESS_TENSOR ANALYTICAL
  &END FORCE_EVAL

&MOTION
   &MD
     ENSEMBLE NVT                #constant temperature and pressure using an isotropic cell
     STEPS 14000
     TIMESTEP 1
     TEMPERATURE 2000
  &THERMOSTAT
    &NOSE
      TIMECON 1000
    &END NOSE
    &END THERMOSTAT
   &END MD
  &PRINT
    &TRAJECTORY
      &EACH
        MD     1 #Output frequency of geometry
      &END EACH
      FORMAT xyz
    &END TRAJECTORY
    &VELOCITIES
      &EACH
        MD     1 #Output frequency of velocity
      &END EACH
    &END VELOCITIES
    &RESTART
      BACKUP_COPIES 0 #Maximum number of backup copies of MD restart file
      &EACH
        MD 10 #Frequency of updating restart file
      &END EACH
    &END RESTART
  &END PRINT
&END MOTION

wzkchem5

你确定初始构象合理？光是弛豫过不能保证合理，还要看过弛豫后的结构，确定是你想要的结构才能保证合理

anrushan

wzkchem5 发表于 2021-3-7 18:30
你确定初始构象合理？光是弛豫过不能保证合理，还要看过弛豫后的结构，确定是你想要的结构才能保证合理

确定合理，上面的in文件坐标部分不全，在弛豫过后检查过输出的XYZ轨迹，并没有发生结构变化，原子间的距离也合理，我曾经用初始结构跑了一下上面in文件里的nvt收敛的可以，不明白为什么跑了弛豫之后会出现不收敛的状况，谢谢老师的提醒

anrushan

升高了收敛的上限，观察输出文件发现结构已经散了，是方法有问题。

sobereva

没你的体系细节信息没法说，起码给个当前体系结构图，并且显示出盒子
如果是单个分子，用ORCA甭用CP2K
使用ORCA做从头算动力学(AIMD)的简单例子
http://sobereva.com/576（http://bbs.keinsci.com/thread-20800-1-1.html）

anrushan

sobereva 发表于 2021-3-7 23:53
没你的体系细节信息没法说，起码给个当前体系结构图，并且显示出盒子
如果是单个分子，用ORCA甭用CP2K
使 ...

卢老师，如果只考虑分子在动力学中的轨迹问题是不是不用关注SCF的收敛（分子在高温条件下发生的分解反应是否也会影响SCF方程的收敛呢）？

sobereva

anrushan 发表于 2021-3-8 20:50
卢老师，如果只考虑分子在动力学中的轨迹问题是不是不用关注SCF的收敛（分子在高温条件下发生的分解反应 ...

需要关注
如果某一步没收敛，而且最后一轮SCF的时候收敛程度很差，那么这一步的受力就会算得很糟糕，轻则不准确，重则造成动力学不稳定

收敛成功几率和这一步结构与上一步结构差异有关，差异越小收敛成功几率越大。高温下如果发现中途难收敛，尝试降低步长

anrushan

sobereva 发表于 2021-3-9 04:19
需要关注
如果某一步没收敛，而且最后一轮SCF的时候收敛程度很差，那么这一步的受力就会算得很糟糕，轻 ...

谢谢卢老师，注意休息

LLYYTT

我也遇到了这个问题，请问应该怎么解决呢