Constrained DFT (CDFT) Errors

[Brian Day]

Hi everyone,

I am currently trying to run a CDFT calculation to determine the electronic coupling integral of a charged system in CP2K. I have successfully run the Zn dimer tutorial, but am running into issues when I switch to my system; specifically, I get the following error: 

===================================================================================

=   BAD TERMINATION OF ONE OF YOUR APPLICATION PROCESSES

=   PID 21797 RUNNING AT opa-n68

=   EXIT CODE: 9

=   CLEANING UP REMAINING PROCESSES

=   YOU CAN IGNORE THE BELOW CLEANUP MESSAGES

===================================================================================

   Intel(R) MPI Library troubleshooting guide:

      https://software.intel.com/node/561764

===================================================================================

I have attached all files used to submit and run this job, including the two outputs (Ni3HITP2-Mcent-cdft.out = cluster outputs, Ni3HITP2-Mcent-state1-run.out = cp2k outputs). 

Note about the system: I am studying is a bonded system, and I am localizing charge onto a certain part of the molecule, it is not two distinct molecules (technically, this is a periodic system, but I am using an H-terminated unit to approximate the framework). 

Additional questions in case anyone happens to know the answers:

1. If the net charge on the a region of the neutral framework is non-zero, should the charged state contain one more electron than the neutral-state charge in that region, or simply one extra valence electron total? (Although I do plan to test this once I get the code running)

2. If I am using an atom like copper, and a basis set which only contains 11 valence electrons, should the 'strength' term in CDFT reflect the number of valence electrons in the basis set, or the number of electrons total (i.e. 11 vs 29)?

Thanks!

-Brian

[Nico Holmberg]

Hi Brian, 

I'm not 100 % sure but I think that "EXIT CODE: 9" means that your calculation ran out of memory. Try increasing the number of nodes or initially for testing purposes decrease the amount of memory required by using: smaller cutoff and/or smaller cell (your should be able to get away with a smaller cell with wavelet, but do test this).

By the way, I noticed you're using Diagonalization instead of OT as the solver. If at all possible for your system (converges, produces reasonable results), switch to OT because in my experience CDFT converges much more nicely with OT than diagonalization.

To your questions, 

1) Do you want the charge to be positive or negative? If in the neutral case the CDFT constraint group has X electrons, then X+1 electrons will give it a negative -1 charge, and X-1 electrons will give it a positive +1 charge. The total number of electrons must be conserved in both states if you subsequently plan on calculating the electronic coupling. You can also calculate the energy difference between CDFT states with a different total net charge as long as you take into account the changes in charge properly (but not any properties called with MIXED_CDFT section). 

2) I guess you mean the "TARGET" keyword not "STRENGTH". Use the number of valence electrons as determined by the pseudo and basis set if you're doing GPW calculations (not all-electron GAPW). So a neutral Cu atom would have 11 valence electrons (target 11). 

BR, 

Nico

maanantai 22. lokakuuta 2018 17.44.59 UTC+3 Brian Day kirjoitti:

Hi everyone,

I am currently trying to run a CDFT calculation to determine the electronic coupling integral of a charged system in CP2K. I have successfully run the Zn dimer tutorial, but am running into issues when I switch to my system; specifically, I get the following error: 

===================================================================================

=   BAD TERMINATION OF ONE OF YOUR APPLICATION PROCESSES

=   PID 21797 RUNNING AT opa-n68

=   EXIT CODE: 9

=   CLEANING UP REMAINING PROCESSES

=   YOU CAN IGNORE THE BELOW CLEANUP MESSAGES

===================================================================================

   Intel(R) MPI Library troubleshooting guide:

      https://software.intel.com/node/561764

===================================================================================

I have attached all files used to submit and run this job, including the two outputs (Ni3HITP2-Mcent-cdft.out = cluster outputs, Ni3HITP2-Mcent-state1-run.inp = cp2k outputs). 

[Brian Day]

Hi Nico,

Thanks for your reply.

I was able to get the calculation running by decreasing the box size (I did not initially realize that the memory scaled that quickly with box size, I thought it was mostly impacted by number of atoms). I do get some warnings though, which I was hoping you might be able to explain.

When running the CDFT calculation for one of the states, I get the following:

 *** WARNING in pw/pw_pool_types.F:470 :: hit max_cache ***

When I try to run a calulcation for the electronic coupling integral, I get the following set or errors:

 *** WARNING in mixed_cdft_methods.F:1192 :: The usage of SVD based matrix ***

 *** inversions with fractionally occupied orbitals is strongly            ***

 *** recommended to screen nearly orthogonal states.                       ***

 *** WARNING in mixed_cdft_methods.F:1241 :: The number of alpha electrons ***

 *** is not constant across all CDFT states. Calculation proceeds but the  ***

 *** results will likely be meaningless.                                   ***

 *** WARNING in mixed_cdft_methods.F:1245 :: The number of beta electrons ***

 *** is not constant across all CDFT states. Calculation proceeds but the ***

 *** results will likely be meaningless.                                  ***

Lastly, this is an error I have gotten with almost all of my charged state calculations that I do not know how to handle:

WARNING: S**2 computation does not yet treat fractional occupied orbitals

The files from the CDFT calculation set can be found here: https://drive.google.com/open?id=1k3D-cNtOUmJSxW4lzVrC_d7Zu9uTtn_O (Also, I just noticed that I didn't call the restart files correctly, so I will be fixing that and repeating the run, in case that is the cause of any of the above issues).

Thanks again!

-Brian

[Nico Holmberg]

Hi,

The first error is harmless. When you compute CDFT charges, the atomic weight functions for each atom are stored as separate grids. This warning merely states that your simulation allocated (and later released) more of these grids from the grid pool than the hardcoded "limit". I don't know the history behind this choice, but I would say that its safe to ignore this warning. The last error is also harmless: you're using diagonalization with orbital occupation smearing, which cannot be used to compute the spin contamination metric S**2.

The 2nd-4th errors are more serious. These are related to the use of Diagonalization + orbital smearing. As I mentioned in my previous post, the CDFT coupling is well-defined only when computed between states with the same total number of electrons. This also applies to spins: the states should be in the same spin state, i.e., the number of alpha/beta electrons should be the same in both states. These errors indicate that the number of alpha/beta electrons is not the same in both states, which are computed as the sum of orbital occupations (you can print these if you want to check for yourself). If you want to keep using diagonalization, you need to ensure that the calculation preserves the number of alpha/beta electrons. You can use the keyword FIXED_MAGNETIC_MOMENT to fix the number spin electrons. I encourage you to test this keyword first without CDFT. 

You also have another choice if the difference between the number of alpha/beta electrons is small (ideally <0.1). The CDFT calculation considers all MOs with occupation larger than EPS_OCCUPIED as occupied orbitals (default value 1e-6). It's possible that you could remove some of the fractionally occupied orbitals by increasing EPS_OCCUPIED, which could then ensure that both states are in the same spin state. In my experience, fractionally occupied orbitals can exhibit linear dependencies causing the coupling to vanish. You can use a singular value decomposition (keyword EPS_SVD ) to screen linearly-dependent MOs by setting EPS_SVD to some small value e.g. 1e-3. In this case, only those MO pais whose singular value is above 1e-3 is included in the CDFT calculation. You can use the print keyword MO_OVERLAP_EIGENVALUES  to see which MO pairs get eliminated if you vary EPS_SVD.

As you can see using diagonalization with CDFT can be quite complex. If possible I would recommend switching to OT.

- Nico

[Brian Day]

Hi Nico,

I'm pleased to report that after switching to OT, the only warning that remains is  *** WARNING in pw/pw_pool_types.F:470 :: hit max_cache ***.

If you wouldn't mind answering one last question regarding the output, that would be appreciated. Specifically, what are the differences between "Overlap between states I and J", "Charge transfer energy (J-I)", and "Diabatic electronic coupling" as listed below. I'm am running these calculations to calculate electron transfer rates using Marcus Theory, so presumably I should use the electronic coupling value, but wanted to check. 

  ------------------------- CDFT coupling information --------------------------

  Information at step (fs):                                                 0.00

  ############################################

  ###### CDFT states I =  1 and J =   2 ######

  ############################################

  Atomic group:                                                                1

  Strength of constraint I:                                      -0.254713549611

  Strength of constraint J:                                      -0.256763104414

  Final value of constraint I:                                  119.000059959930

  Final value of constraint J:                                  119.000008498861

  Overlap between states I and J:                                 0.069002784505

  Charge transfer energy (J-I) (Hartree):                         0.001330331665

  Diabatic electronic coupling (Lowdin, mHartree):               19.740562532228

  ------------------------------------------------------------------------------

 NO FORCE_EVAL section calculated the dipole

 ENERGY| Total FORCE_EVAL ( MIXED ) energy (a.u.):          -512.280992021731890

Thanks!

-Brian

[Nico Holmberg]

Hi Brian,

I am glad to hear that switching to OT resolved your issue. The definitions for the terms you asked about are the following

• Charge transfer energy (J-I): This is the (vertical) energy difference between the two CDFT states. You can compute the solvent reorganization energy and the reaction (free) energy terms from this quantity, which both appear in Marcus theory. In practice, you would need to sample this quantity over a CDFT MD simulation, see e.g. the first CP2K CDFT implementation paper. 
• Diabatic electronic coupling: This is the effective coupling between the CDFT states after orthogonalization, i.e., the term you would use for H_ab in Marcus' rate equation for ET (see e.g. the CP2K CDFT tutorial for the actual equation).
• Overlap between states: This is the value of the overlap integral between CDFT states, which are represented by Kohn-Sham Slater determinants (not true wavefunctions). This quantity is needed to compute the electronic coupling, see the CP2K CDFT tutorial for the relevant equations. This quantity measures how similar (value close to 1) or dissimilar (close to 0) the involved CDFT states are.

BR,

Nico  

[Brian Day]

Hi Nico,

I have some follow up questions about CDFT if you don't mind. I figured I'd post them here so that many CDFT questions were addressed in the same thread. 

1. In a few mixed CDFT calculations, none of the following files were printed: *.inverseJacobian, *-JacobianInfo.out, and -LineSearch.out. Fortunately, there were also no warnings in the *mixed-run.out, but I was curious as to why these files never got printed. In other runs, where there are errors such as LS not converging, these files all get printed. Do these files only get printed when there are errors?

2. I am frequently having issues with the LineSearch not converging, but even when I increase to 50 iterations for the line search, I still get that issue? Is that a serious problem, and do you know of any ways to solve it? 

3. With regards to the coupling reliability metric, it seems that the paper implies that roughly anything less than 0.10 can be considered reliable - Would you agree that this in an acceptable cutoff? Additionally, how serious are more moderate coupling metrics, such as 0.15-0.25? Lastly, is there anything to your knowledge that can be changed in the code to reduce this value towards 0, or it it just system dependent? (I realize there is probably not a clean-cut answer to any of these questions, but any extra insight would be a huge help.)

Best,

     Brian