Moderator: moderators

Fangyong
Posts: 9
Joined: Mon Oct 18, 2021 2:28 pm

Dear all,

is it possible to print the electron density at the bond critical point?

Below is a water molecule density cube file information from cp2k:
320 0.094486 0.000000 0.000000
320 0.000000 0.094486 0.000000
320 0.000000 0.000000 0.094486
8 0.000000 -0.000000 -0.000000 -0.147753
1 0.000000 -0.000000 -1.438202 0.996119
1 0.000000 -0.000000 1.438202 0.996119
...... (electron densities on the grid)

and I use your code to calculate the Bader charge and the bond critical points, ./bader *cube -cp -ls, and I found 2 bond critical points, and 3 nuclei critical points, for the bond critical points, and the details of the critical points can be found in a file called CPFU.dat, here is the bond critical point information, but the electron density at the bond critical point is not printed. Could you please help me with this? Because we need the electron density at the bond critical point to show the bond strength for certain bonds.

We could have just searched the bond critical points in the density cube file, however, since the program has already found the bond critical point, it is convenient to just print the electron density at the point.

Thank you very much!!

Fangyong

Below is the bond critical point information in the CPFU.dat,
************************** (no electron density) ***************
Critical point is found at indices
1.00000000000000 309.125531411542 9.00703206661725
Coordinates in cartesian are
9.448600000000000E-002 29.2080349609510 0.851038431846398
Direct coordinates are
3.125000000000000E-003 0.966017285661070 2.814697520817892E-002
-9.968266401243480E-015 2.674152997699508E-002 2.564549130204341E-002
Hessian is
-42939.6955373299 6.093969405902629E-014 -4.519952531199678E-014
6.093969405902629E-014 1121.76335161199 -33512.6781478708
-4.519952531199678E-014 -33512.6781478708 -15272.6724894016
Eigenvalues are
-42939.6955373299 27425.1827915809 -41576.0919293705
Eigenvectors are
1.00000000000000 -1.784920576630989E-018 -5.283525974885847E-019
1.077871804327750E-018 0.786637187664234 -0.617415528622099
1.517659482752596E-018 0.617415528622099 0.786637187664234
This is a bond critical point
Fangyong
Posts: 9
Joined: Mon Oct 18, 2021 2:28 pm

I think we should just use the density cube file to obtain the electron density on the bond critical point, since we have already found the bond critical point.

Thanks!

Fangyong
Fangyong
Posts: 9
Joined: Mon Oct 18, 2021 2:28 pm

Dear Professor Henkelman,

The bond critical point is not located at the grid point, so I cannot get its electron density from the cube file.

Fangyong
graeme
Posts: 2054
Joined: Tue Apr 26, 2005 4:25 am
Contact:

Sorry for the slow response - I have tasked the person who wrote this code to provide the answer (it should). I will ping him again.
Fangyong
Posts: 9
Joined: Mon Oct 18, 2021 2:28 pm

Dear Professor Henkelman,

Thanks for the reply! Please let me know when the code for printing the electron density at the bond critical point is ready. (Bader's definition of bonding really gives one of the best representations of chemical bonds, considering the ambiguity of defining a bond in physical theory.)

Thank you so much!

Fangyong
wc5879
Posts: 17
Joined: Thu May 22, 2014 7:49 pm

If you check out the latest version through svn, the interpolated charge density is printed out under coordinates.
Fangyong
Posts: 9
Joined: Mon Oct 18, 2021 2:28 pm

Dear wc5879,

How do I check out the latest version through svn? Also, does the latest version contain the execution file, or source code, or both?

Thank you very much!

Fangyong
Fangyong
Posts: 9
Joined: Mon Oct 18, 2021 2:28 pm

Dear Ray (wc5879),

Regards the most up to date source code and your comment:

"Now the charge density is written correctly.
To obtain the most up to date source code, use

I downloaded it and compile it, and used the code to calculate bond critical points for a molecule structure, however, it says:

This function is broken right now.

When I used the code released on 02/08/20 - Version 1.04 Released, I am able to calculate the bond critical points, and it says:
CRITICAL POINTS INFO WRITEN TO CPF.dat
CRITICAL POINTS FOUND: 13880
FINDING UNIQUE CRITICAL POINTS...
Found a unique bond critical point
Found a unique bond critical point
Found a unique bond critical point

So I think the most up to date source code has not updated the critical point calculation yet, because it says the "function is broken right now."

So I compare the code before update and after update, for subroutine critpoint_mod.f90, and found where to print the electron density, so for the before-update version, in line 552,
line 550, !PRINT *, SQRT(cpl(i)%force(1)**2 + cpl(i)%force(2)**2 + &
line 551, ! cpl(i)%force(3)**2)
line 552, I add WRITE (98,*) "Charge density is"
WRITE (98,*) cpl(i)%rho/vol

when I complied it, it tells that

critpoint_mod.f90:553:35:

WRITE (98,*) cpl(i)%rho/vol
1
Error: ‘rho’ at (1) is not a member of the ‘cpc’ structure

So, in summary:
1. looks like critical point does not work for the most update code (that's say it is testing version), but works for the code before update (the code released on 02/08/20 - Version 1.04, that is say it is a stable version)

2. is it possible to write a line for the electron density on the stable version, such as
line 552, I add WRITE (98,*) "Charge density is"
WRITE (98,*) cpl(i)%rho/vol

since I am not familiar with your code, I dont know where the electron density is located, but should not be cpl(i)%rho.
Could you please add such a line, so we can use the stable version and calculate the electron density on the bond critical point?

Thanks!

Fangyong
Fangyong
Posts: 9
Joined: Mon Oct 18, 2021 2:28 pm

Dear Ray and Professor Henkelman,

How do you obtain the charge density at the bond critical point using interpolation? What is the equation for the interpolation?

It will be better if you can add line for printing the electron density, however, if you can tell me the equation for interpolation, maybe I can calculate the electron
density myself.

Thanks!

Fangyong
Fangyong
Posts: 9
Joined: Mon Oct 18, 2021 2:28 pm