diff --git a/src/aobasis/basis_set_types.F b/src/aobasis/basis_set_types.F
index 3820ca5c71..3b2b3555fa 100644
--- a/src/aobasis/basis_set_types.F
+++ b/src/aobasis/basis_set_types.F
@@ -85,7 +85,7 @@ MODULE basis_set_types
m => NULL(), ncgf_set => NULL(), &
npgf => NULL(), nsgf_set => NULL(), nshell => NULL()
REAL(KIND=dp), DIMENSION(:, :), POINTER :: cphi => NULL(), pgf_radius => NULL(), sphi => NULL(), &
- scon => NULL(), zet => NULL()
+ scon => NULL(), zet => NULL(), ccon => NULL()
INTEGER, DIMENSION(:, :), POINTER :: first_cgf => NULL(), first_sgf => NULL(), l => NULL(), &
last_cgf => NULL(), last_sgf => NULL(), n => NULL()
REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: gcc => NULL()
@@ -193,6 +193,7 @@ CONTAINS
IF (ASSOCIATED(gto_basis_set%pgf_radius)) DEALLOCATE (gto_basis_set%pgf_radius)
IF (ASSOCIATED(gto_basis_set%sphi)) DEALLOCATE (gto_basis_set%sphi)
IF (ASSOCIATED(gto_basis_set%scon)) DEALLOCATE (gto_basis_set%scon)
+ IF (ASSOCIATED(gto_basis_set%ccon)) DEALLOCATE (gto_basis_set%ccon)
IF (ASSOCIATED(gto_basis_set%zet)) DEALLOCATE (gto_basis_set%zet)
IF (ASSOCIATED(gto_basis_set%first_cgf)) DEALLOCATE (gto_basis_set%first_cgf)
IF (ASSOCIATED(gto_basis_set%first_sgf)) DEALLOCATE (gto_basis_set%first_sgf)
@@ -254,9 +255,16 @@ CONTAINS
basis_set_out%nsgf_set = basis_set_in%nsgf_set
maxco = SIZE(basis_set_in%cphi, 1)
ALLOCATE (basis_set_out%cphi(maxco, ncgf), basis_set_out%sphi(maxco, nsgf), basis_set_out%scon(maxco, nsgf))
+ ALLOCATE (basis_set_out%ccon(maxco, ncgf))
basis_set_out%cphi = basis_set_in%cphi
basis_set_out%sphi = basis_set_in%sphi
basis_set_out%scon = basis_set_in%scon
+ basis_set_out%ccon = 0.0_dp
+ IF (ASSOCIATED(basis_set_in%ccon)) THEN
+ IF ((SIZE(basis_set_in%ccon, 1) == maxco) .AND. (SIZE(basis_set_in%ccon, 2) == ncgf)) THEN
+ basis_set_out%ccon = basis_set_in%ccon
+ END IF
+ END IF
maxpgf = MAXVAL(basis_set_in%npgf)
ALLOCATE (basis_set_out%pgf_radius(maxpgf, nset), basis_set_out%zet(maxpgf, nset))
basis_set_out%pgf_radius = basis_set_in%pgf_radius
@@ -443,6 +451,8 @@ CONTAINS
pbasis%sphi = 0.0_dp
ALLOCATE (pbasis%scon(maxco, ncgf))
pbasis%scon = 0.0_dp
+ ALLOCATE (pbasis%ccon(maxco, ncgf))
+ pbasis%ccon = 0.0_dp
ALLOCATE (pbasis%set_radius(nset))
ALLOCATE (pbasis%pgf_radius(mpgf, nset))
pbasis%pgf_radius = 0.0_dp
@@ -570,6 +580,7 @@ CONTAINS
CALL reallocate(basis_set%cphi, 1, maxco, 1, ncgf)
CALL reallocate(basis_set%sphi, 1, maxco, 1, nsgf)
CALL reallocate(basis_set%scon, 1, maxco, 1, nsgf)
+ CALL reallocate(basis_set%ccon, 1, maxco, 1, ncgf)
CALL reallocate(basis_set%pgf_radius, 1, maxpgf, 1, nset)
END SUBROUTINE combine_basis_sets
@@ -622,12 +633,13 @@ CONTAINS
!> \param maxder ...
!> \param short_kind_radius ...
!> \param npgf_seg_sum number of primitives in "segmented contraction format"
+!> \param ccon ...
! **************************************************************************************************
SUBROUTINE get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, &
nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, &
m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, &
last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, &
- npgf_sum, nshell_sum, maxder, short_kind_radius, npgf_seg_sum)
+ npgf_sum, nshell_sum, maxder, short_kind_radius, npgf_seg_sum, ccon)
! Get informations about a Gaussian-type orbital (GTO) basis set.
@@ -655,6 +667,7 @@ CONTAINS
INTEGER, INTENT(IN), OPTIONAL :: maxder
REAL(KIND=dp), INTENT(OUT), OPTIONAL :: short_kind_radius
INTEGER, INTENT(OUT), OPTIONAL :: npgf_seg_sum
+ REAL(KIND=dp), DIMENSION(:, :), OPTIONAL, POINTER :: ccon
INTEGER :: iset, nder
@@ -684,6 +697,7 @@ CONTAINS
IF (PRESENT(pgf_radius)) pgf_radius => gto_basis_set%pgf_radius
IF (PRESENT(sphi)) sphi => gto_basis_set%sphi
IF (PRESENT(scon)) scon => gto_basis_set%scon
+ IF (PRESENT(ccon)) ccon => gto_basis_set%ccon
IF (PRESENT(zet)) zet => gto_basis_set%zet
IF (PRESENT(first_cgf)) first_cgf => gto_basis_set%first_cgf
IF (PRESENT(first_sgf)) first_sgf => gto_basis_set%first_sgf
@@ -789,7 +803,7 @@ CONTAINS
END SELECT
! Initialise the transformation matrices "pgf" -> "cgf"
- CALL init_cphi_and_sphi(gto_basis_set)
+ CALL init_cphi_and_sphi(gto_basis_set, .FALSE.)
CALL timestop(handle)
@@ -798,8 +812,9 @@ CONTAINS
! **************************************************************************************************
!> \brief ...
!> \param gto_basis_set ...
+!> \param lccon ...
! **************************************************************************************************
- SUBROUTINE init_cphi_and_sphi(gto_basis_set)
+ SUBROUTINE init_cphi_and_sphi(gto_basis_set, lccon)
! Initialise the matrices for the transformation of primitive Cartesian
! Gaussian-type functions to contracted Cartesian (cphi) and spherical
@@ -808,6 +823,7 @@ CONTAINS
! - Creation (20.09.2000,MK)
TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
+ LOGICAL, INTENT(IN), OPTIONAL :: lccon
CHARACTER(len=*), PARAMETER :: routineN = 'init_cphi_and_sphi'
@@ -815,7 +831,10 @@ CONTAINS
ipgf, iset, ishell, l, last_sgf, lmax, &
lmin, n, n1, n2, ncgf, nn, nn1, nn2, &
npgf, nsgf
+ LOGICAL :: my_lccon
+ my_lccon = .FALSE.
+ IF (PRESENT(lccon)) my_lccon = lccon
! -------------------------------------------------------------------------
! Build the Cartesian transformation matrix "cphi"
@@ -896,6 +915,36 @@ CONTAINS
END DO
END IF
+ IF (my_lccon) THEN
+ IF (.NOT. ASSOCIATED(gto_basis_set%ccon)) THEN
+ CALL reallocate(gto_basis_set%ccon, 1, SIZE(gto_basis_set%cphi, 1), 1, gto_basis_set%ncgf)
+ ELSE IF ((SIZE(gto_basis_set%ccon, 1) /= SIZE(gto_basis_set%cphi, 1)) .OR. &
+ (SIZE(gto_basis_set%ccon, 2) /= gto_basis_set%ncgf)) THEN
+ CALL reallocate(gto_basis_set%ccon, 1, SIZE(gto_basis_set%cphi, 1), 1, gto_basis_set%ncgf)
+ END IF
+ n = SIZE(gto_basis_set%ccon, 1)
+ gto_basis_set%ccon = 0.0_dp
+ IF (n > 0) THEN
+ DO iset = 1, gto_basis_set%nset
+ lmin = gto_basis_set%lmin(iset)
+ lmax = gto_basis_set%lmax(iset)
+ npgf = gto_basis_set%npgf(iset)
+ nn = ncoset(lmax) - ncoset(lmin - 1)
+ DO ishell = 1, gto_basis_set%nshell(iset)
+ first_sgf = gto_basis_set%first_cgf(ishell, iset)
+ last_sgf = gto_basis_set%last_cgf(ishell, iset)
+ DO ipgf = 1, npgf
+ nn1 = (ipgf - 1)*ncoset(lmax) + ncoset(lmin - 1) + 1
+ nn2 = ipgf*ncoset(lmax)
+ n1 = (ipgf - 1)*nn + 1
+ n2 = ipgf*nn
+ gto_basis_set%ccon(n1:n2, first_sgf:last_sgf) = gto_basis_set%cphi(nn1:nn2, first_sgf:last_sgf)
+ END DO
+ END DO
+ END DO
+ END IF
+ END IF
+
CALL timestop(handle)
END SUBROUTINE init_cphi_and_sphi
@@ -1118,7 +1167,7 @@ CONTAINS
! Initialise the transformation matrices "pgf" -> "cgf"
- CALL init_cphi_and_sphi(gto_basis_set)
+ CALL init_cphi_and_sphi(gto_basis_set, .TRUE.)
CALL timestop(handle)
@@ -1255,6 +1304,7 @@ CONTAINS
CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
+ CALL reallocate(gto_basis_set%ccon, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%lx, 1, ncgf)
CALL reallocate(gto_basis_set%ly, 1, ncgf)
CALL reallocate(gto_basis_set%lz, 1, ncgf)
@@ -1410,6 +1460,7 @@ CONTAINS
CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
+ CALL reallocate(gto_basis_set%ccon, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%lx, 1, ncgf)
CALL reallocate(gto_basis_set%ly, 1, ncgf)
CALL reallocate(gto_basis_set%lz, 1, ncgf)
@@ -1561,6 +1612,7 @@ CONTAINS
CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
+ CALL reallocate(gto_basis_set%ccon, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%lx, 1, ncgf)
CALL reallocate(gto_basis_set%ly, 1, ncgf)
CALL reallocate(gto_basis_set%lz, 1, ncgf)
@@ -1713,6 +1765,7 @@ CONTAINS
CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
+ CALL reallocate(gto_basis_set%ccon, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%lx, 1, ncgf)
CALL reallocate(gto_basis_set%ly, 1, ncgf)
CALL reallocate(gto_basis_set%lz, 1, ncgf)
@@ -1793,13 +1846,14 @@ CONTAINS
!> \param n ...
!> \param gcc ...
!> \param short_kind_radius ...
+!> \param ccon ...
!> \author MK
! **************************************************************************************************
SUBROUTINE set_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, &
nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, &
lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, &
cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, &
- last_cgf, last_sgf, n, gcc, short_kind_radius)
+ last_cgf, last_sgf, n, gcc, short_kind_radius, ccon)
TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
CHARACTER(LEN=default_string_length), INTENT(IN), &
@@ -1818,6 +1872,7 @@ CONTAINS
REAL(KIND=dp), DIMENSION(:, :, :), OPTIONAL, &
POINTER :: gcc
REAL(KIND=dp), INTENT(IN), OPTIONAL :: short_kind_radius
+ REAL(KIND=dp), DIMENSION(:, :), OPTIONAL, POINTER :: ccon
IF (PRESENT(name)) gto_basis_set%name = name
IF (PRESENT(aliases)) gto_basis_set%aliases = aliases
@@ -1845,6 +1900,7 @@ CONTAINS
IF (PRESENT(pgf_radius)) gto_basis_set%pgf_radius(:, :) = pgf_radius(:, :)
IF (PRESENT(sphi)) gto_basis_set%sphi(:, :) = sphi(:, :)
IF (PRESENT(scon)) gto_basis_set%scon(:, :) = scon(:, :)
+ IF (PRESENT(ccon)) gto_basis_set%ccon(:, :) = ccon(:, :)
IF (PRESENT(zet)) gto_basis_set%zet(:, :) = zet(:, :)
IF (PRESENT(first_cgf)) gto_basis_set%first_cgf(:, :) = first_cgf(:, :)
IF (PRESENT(first_sgf)) gto_basis_set%first_sgf(:, :) = first_sgf(:, :)
@@ -2050,6 +2106,8 @@ CONTAINS
WRITE (UNIT=output_unit, FMT="(12F10.5)") gto_basis_set%sphi
WRITE (UNIT=output_unit, FMT="(A1)") "SCON"
WRITE (UNIT=output_unit, FMT="(12F10.5)") gto_basis_set%scon
+ WRITE (UNIT=output_unit, FMT="(A1)") "CCON"
+ WRITE (UNIT=output_unit, FMT="(12F10.5)") gto_basis_set%ccon
END IF
@@ -2536,6 +2594,7 @@ CONTAINS
CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
+ CALL reallocate(gto_basis_set%ccon, 1, maxco, 1, ncgf)
CALL reallocate(gto_basis_set%lx, 1, ncgf)
CALL reallocate(gto_basis_set%ly, 1, ncgf)
CALL reallocate(gto_basis_set%lz, 1, ncgf)
@@ -2711,11 +2770,12 @@ CONTAINS
INTEGER, ALLOCATABLE, DIMENSION(:) :: sort_index
INTEGER, ALLOCATABLE, DIMENSION(:, :) :: icgf_set, isgf_set
INTEGER, DIMENSION(:), POINTER :: lx, ly, lz, m, npgf
+ LOGICAL :: ccon_available
REAL(dp), ALLOCATABLE, DIMENSION(:) :: tmp
REAL(dp), DIMENSION(:), POINTER :: set_radius
REAL(dp), DIMENSION(:, :), POINTER :: zet
REAL(KIND=dp), DIMENSION(:), POINTER :: norm_cgf
- REAL(KIND=dp), DIMENSION(:, :), POINTER :: cphi, scon, sphi
+ REAL(KIND=dp), DIMENSION(:, :), POINTER :: ccon, cphi, scon, sphi
NULLIFY (set_radius, zet)
@@ -2798,6 +2858,13 @@ CONTAINS
sphi = 0.0_dp
ALLOCATE (scon(SIZE(basis_set%scon, 1), SIZE(basis_set%scon, 2)))
scon = 0.0_dp
+ ALLOCATE (ccon(SIZE(basis_set%cphi, 1), SIZE(basis_set%cphi, 2)))
+ ccon = 0.0_dp
+ ccon_available = ASSOCIATED(basis_set%ccon)
+ IF (ccon_available) THEN
+ ccon_available = (SIZE(basis_set%ccon, 1) == SIZE(ccon, 1)) .AND. &
+ (SIZE(basis_set%ccon, 2) == SIZE(ccon, 2))
+ END IF
ALLOCATE (sgf_symbol(SIZE(basis_set%sgf_symbol)))
ALLOCATE (m(SIZE(basis_set%m)))
@@ -2814,6 +2881,7 @@ CONTAINS
ly(icgf_new) = basis_set%ly(icgf_old)
lz(icgf_new) = basis_set%lz(icgf_old)
cphi(:, icgf_new) = basis_set%cphi(:, icgf_old)
+ IF (ccon_available) ccon(:, icgf_new) = basis_set%ccon(:, icgf_old)
cgf_symbol(icgf_new) = basis_set%cgf_symbol(icgf_old)
END DO
DO is = 1, is_max
@@ -2843,6 +2911,8 @@ CONTAINS
basis_set%sphi => sphi
DEALLOCATE (basis_set%scon)
basis_set%scon => scon
+ IF (ASSOCIATED(basis_set%ccon)) DEALLOCATE (basis_set%ccon)
+ basis_set%ccon => ccon
DEALLOCATE (basis_set%m)
basis_set%m => m
diff --git a/src/cp_dbcsr_output.F b/src/cp_dbcsr_output.F
index 45f6b5c4bd..2bb78a1201 100644
--- a/src/cp_dbcsr_output.F
+++ b/src/cp_dbcsr_output.F
@@ -35,7 +35,8 @@ MODULE cp_dbcsr_output
USE machine, ONLY: m_flush
USE mathlib, ONLY: symmetrize_matrix
USE message_passing, ONLY: mp_para_env_type
- USE orbital_pointers, ONLY: nso
+ USE orbital_pointers, ONLY: nco,&
+ nso
USE particle_methods, ONLY: get_particle_set
USE particle_types, ONLY: particle_type
USE qs_environment_types, ONLY: get_qs_env,&
@@ -156,10 +157,11 @@ CONTAINS
!> \param scale ...
!> \param output_unit ...
!> \param omit_headers Write only the matrix data, not the row/column headers
+!> \param cartesian_basis Use Cartesian instead of spherical basis labels
! **************************************************************************************************
SUBROUTINE cp_dbcsr_write_sparse_matrix(sparse_matrix, before, after, qs_env, para_env, &
first_row, last_row, first_col, last_col, scale, &
- output_unit, omit_headers)
+ output_unit, omit_headers, cartesian_basis)
TYPE(dbcsr_type) :: sparse_matrix
INTEGER, INTENT(IN) :: before, after
@@ -168,11 +170,12 @@ CONTAINS
INTEGER, INTENT(IN), OPTIONAL :: first_row, last_row, first_col, last_col
REAL(dp), INTENT(IN), OPTIONAL :: scale
INTEGER, INTENT(IN) :: output_unit
- LOGICAL, INTENT(IN), OPTIONAL :: omit_headers
+ LOGICAL, INTENT(IN), OPTIONAL :: omit_headers, cartesian_basis
CHARACTER(LEN=default_string_length) :: matrix_name
INTEGER :: col1, col2, dim_col, dim_row, row1, row2
- LOGICAL :: my_omit_headers, print_sym
+ LOGICAL :: my_cartesian_basis, my_omit_headers, &
+ print_sym
REAL(KIND=dp), DIMENSION(:, :), POINTER :: matrix
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
@@ -233,11 +236,14 @@ CONTAINS
ELSE
my_omit_headers = .FALSE.
END IF
+ my_cartesian_basis = .FALSE.
+ IF (PRESENT(cartesian_basis)) my_cartesian_basis = cartesian_basis
CALL dbcsr_get_info(sparse_matrix, name=matrix_name)
IF (print_sym) THEN
CALL write_matrix_sym(matrix, matrix_name, before, after, qs_env, para_env, &
- row1, row2, col1, col2, output_unit, my_omit_headers)
+ row1, row2, col1, col2, output_unit, my_omit_headers, &
+ cartesian_basis=my_cartesian_basis)
ELSE
CALL write_matrix_gen(matrix, matrix_name, before, after, para_env, &
row1, row2, col1, col2, output_unit, my_omit_headers)
@@ -328,10 +334,11 @@ CONTAINS
!> \param last_col ...
!> \param output_unit ...
!> \param omit_headers Write only the matrix data, not the row/column headers
+!> \param cartesian_basis Use Cartesian instead of spherical basis labels
!> \author Creation (01.07.2003,MK)
! **************************************************************************************************
SUBROUTINE write_matrix_sym(matrix, matrix_name, before, after, qs_env, para_env, &
- first_row, last_row, first_col, last_col, output_unit, omit_headers)
+ first_row, last_row, first_col, last_col, output_unit, omit_headers, cartesian_basis)
REAL(KIND=dp), DIMENSION(:, :), POINTER :: matrix
CHARACTER(LEN=*), INTENT(IN) :: matrix_name
@@ -341,7 +348,9 @@ CONTAINS
INTEGER, INTENT(IN) :: first_row, last_row, first_col, &
last_col, output_unit
LOGICAL, INTENT(IN) :: omit_headers
+ LOGICAL, INTENT(IN), OPTIONAL :: cartesian_basis
+ CHARACTER(LEN=12), DIMENSION(:), POINTER :: cgf_symbol
CHARACTER(LEN=2) :: element_symbol
CHARACTER(LEN=25) :: fmtstr1
CHARACTER(LEN=35) :: fmtstr2
@@ -353,6 +362,7 @@ CONTAINS
INTEGER, ALLOCATABLE, DIMENSION(:) :: first_sgf, last_sgf
INTEGER, DIMENSION(:), POINTER :: nshell
INTEGER, DIMENSION(:, :), POINTER :: lshell
+ LOGICAL :: my_cartesian_basis
TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
TYPE(gto_basis_set_type), POINTER :: orb_basis_set
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
@@ -368,6 +378,9 @@ CONTAINS
natom = SIZE(particle_set)
+ my_cartesian_basis = .FALSE.
+ IF (PRESENT(cartesian_basis)) my_cartesian_basis = cartesian_basis
+
CALL get_qs_kind_set(qs_kind_set=qs_kind_set, nsgf=nsgf)
ALLOCATE (first_sgf(natom))
@@ -425,20 +438,27 @@ CONTAINS
CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
IF (ASSOCIATED(orb_basis_set)) THEN
CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
- nset=nset, nshell=nshell, l=lshell, sgf_symbol=sgf_symbol)
+ nset=nset, nshell=nshell, l=lshell, &
+ cgf_symbol=cgf_symbol, sgf_symbol=sgf_symbol)
isgf = 1
DO iset = 1, nset
DO ishell = 1, nshell(iset)
l = lshell(ishell, iset)
- DO iso = 1, nso(l)
+ DO iso = 1, MERGE(nco(l), nso(l), my_cartesian_basis)
IF ((irow >= first_row) .AND. (irow <= last_row)) THEN
IF (omit_headers) THEN
WRITE (UNIT=output_unit, FMT=fmtstr2) &
(matrix(irow, jcol), jcol=from, to)
ELSE
- WRITE (UNIT=output_unit, FMT=fmtstr2) &
- irow, iatom, element_symbol, sgf_symbol(isgf), &
- (matrix(irow, jcol), jcol=from, to)
+ IF (my_cartesian_basis) THEN
+ WRITE (UNIT=output_unit, FMT=fmtstr2) &
+ irow, iatom, element_symbol, cgf_symbol(isgf), &
+ (matrix(irow, jcol), jcol=from, to)
+ ELSE
+ WRITE (UNIT=output_unit, FMT=fmtstr2) &
+ irow, iatom, element_symbol, sgf_symbol(isgf), &
+ (matrix(irow, jcol), jcol=from, to)
+ END IF
END IF
END IF
isgf = isgf + 1
diff --git a/src/input_cp2k_print_dft.F b/src/input_cp2k_print_dft.F
index 53037bbfa2..c6292557b5 100644
--- a/src/input_cp2k_print_dft.F
+++ b/src/input_cp2k_print_dft.F
@@ -549,6 +549,11 @@ CONTAINS
default_l_val=.FALSE., lone_keyword_l_val=.TRUE.)
CALL section_add_keyword(print_key, keyword)
CALL keyword_release(keyword)
+ CALL keyword_create(keyword, __LOCATION__, name="CARTESIAN_OVERLAP", &
+ description="Print the Cartesian overlap matrix in line with Cartesian MO coefficients.", &
+ default_l_val=.FALSE., lone_keyword_l_val=.TRUE.)
+ CALL section_add_keyword(print_key, keyword)
+ CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, name="ENERGIES", &
variants=s2a("EIGENVALUES", "EIGVALS"), &
description="Print the MO energies (eigenvalues).", &
diff --git a/src/particle_methods.F b/src/particle_methods.F
index 8cdaa9dbe9..60858eec45 100644
--- a/src/particle_methods.F
+++ b/src/particle_methods.F
@@ -108,6 +108,7 @@ CONTAINS
!> \param nsgf ...
!> \param nmao ...
!> \param basis ...
+!> \param ncgf ...
!> \date 14.01.2002
!> \par History
!> - particle type cleaned (13.10.2003,MK)
@@ -116,12 +117,13 @@ CONTAINS
!> \version 1.0
! **************************************************************************************************
SUBROUTINE get_particle_set(particle_set, qs_kind_set, first_sgf, last_sgf, nsgf, &
- nmao, basis)
+ nmao, basis, ncgf)
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
INTEGER, DIMENSION(:), INTENT(INOUT), OPTIONAL :: first_sgf, last_sgf, nsgf, nmao
TYPE(gto_basis_set_p_type), DIMENSION(:), OPTIONAL :: basis
+ INTEGER, DIMENSION(:), INTENT(INOUT), OPTIONAL :: ncgf
INTEGER :: ikind, iparticle, isgf, nparticle, ns
@@ -140,6 +142,9 @@ CONTAINS
IF (PRESENT(nmao)) THEN
CPASSERT(SIZE(nmao) >= nparticle)
END IF
+ IF (PRESENT(ncgf)) THEN
+ CPASSERT(SIZE(ncgf) >= nparticle)
+ END IF
IF (PRESENT(first_sgf) .OR. PRESENT(last_sgf) .OR. PRESENT(nsgf)) THEN
isgf = 0
@@ -161,6 +166,22 @@ CONTAINS
END DO
END IF
+ IF (PRESENT(ncgf)) THEN
+ DO iparticle = 1, nparticle
+ CALL get_atomic_kind(particle_set(iparticle)%atomic_kind, kind_number=ikind)
+ IF (PRESENT(basis)) THEN
+ IF (ASSOCIATED(basis(ikind)%gto_basis_set)) THEN
+ CALL get_gto_basis_set(gto_basis_set=basis(ikind)%gto_basis_set, ncgf=ns)
+ ELSE
+ ns = 0
+ END IF
+ ELSE
+ CALL get_qs_kind(qs_kind_set(ikind), ncgf=ns)
+ END IF
+ ncgf(iparticle) = ns
+ END DO
+ END IF
+
IF (PRESENT(first_sgf)) THEN
IF (SIZE(first_sgf) > nparticle) first_sgf(nparticle + 1) = isgf + 1
END IF
diff --git a/src/qs_mo_io.F b/src/qs_mo_io.F
index e95b5dd2f2..31ab05dcea 100644
--- a/src/qs_mo_io.F
+++ b/src/qs_mo_io.F
@@ -22,6 +22,7 @@ MODULE qs_mo_io
USE atomic_kind_types, ONLY: get_atomic_kind
USE basis_set_types, ONLY: get_gto_basis_set,&
+ gto_basis_set_p_type,&
gto_basis_set_type
USE cp_dbcsr_api, ONLY: dbcsr_binary_write,&
dbcsr_create,&
@@ -30,7 +31,9 @@ MODULE qs_mo_io
dbcsr_type
USE cp_dbcsr_contrib, ONLY: dbcsr_checksum
USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
- copy_fm_to_dbcsr
+ copy_fm_to_dbcsr,&
+ dbcsr_deallocate_matrix_set
+ USE cp_dbcsr_output, ONLY: cp_dbcsr_write_sparse_matrix
USE cp_files, ONLY: close_file,&
open_file
USE cp_fm_types, ONLY: cp_fm_get_info,&
@@ -68,13 +71,20 @@ MODULE qs_mo_io
USE qs_density_matrices, ONLY: calculate_density_matrix
USE qs_dftb_types, ONLY: qs_dftb_atom_type
USE qs_dftb_utils, ONLY: get_dftb_atom_param
+ USE qs_environment_types, ONLY: get_qs_env,&
+ qs_environment_type
USE qs_kind_types, ONLY: get_qs_kind,&
get_qs_kind_set,&
qs_kind_type
+ USE qs_ks_types, ONLY: qs_ks_env_type
USE qs_mo_methods, ONLY: calculate_subspace_eigenvalues
USE qs_mo_occupation, ONLY: set_mo_occupation
USE qs_mo_types, ONLY: get_mo_set,&
mo_set_type
+ USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type,&
+ release_neighbor_list_sets
+ USE qs_neighbor_lists, ONLY: setup_neighbor_list
+ USE qs_overlap, ONLY: build_overlap_matrix_simple
#include "./base/base_uses.f90"
IMPLICIT NONE
@@ -995,6 +1005,7 @@ CONTAINS
!> \param cpart ...
!> \param sim_step ...
!> \param umo_set ...
+!> \param qs_env ...
!> \date 15.05.2001
!> \par History:
!> - Optionally print Cartesian MOs (20.04.2005, MK)
@@ -1007,7 +1018,7 @@ CONTAINS
! **************************************************************************************************
SUBROUTINE write_mo_set_to_output_unit(mo_set, qs_kind_set, particle_set, &
dft_section, before, kpoint, final_mos, spin, &
- solver_method, rtp, cpart, sim_step, umo_set)
+ solver_method, rtp, cpart, sim_step, umo_set, qs_env)
TYPE(mo_set_type), INTENT(IN) :: mo_set
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
@@ -1020,6 +1031,7 @@ CONTAINS
LOGICAL, INTENT(IN), OPTIONAL :: rtp
INTEGER, INTENT(IN), OPTIONAL :: cpart, sim_step
TYPE(mo_set_type), INTENT(IN), OPTIONAL :: umo_set
+ TYPE(qs_environment_type), OPTIONAL, POINTER :: qs_env
CHARACTER(LEN=12) :: symbol
CHARACTER(LEN=12), DIMENSION(:), POINTER :: bcgf_symbol
@@ -1036,22 +1048,27 @@ CONTAINS
CHARACTER(LEN=40) :: fmtstr3
CHARACTER(LEN=6), DIMENSION(:), POINTER :: bsgf_symbol
INTEGER :: after, first_mo, from, homo, iatom, icgf, ico, icol, ikind, imo, irow, iset, &
- isgf, ishell, iso, iw, jcol, last_mo, left, lmax, lshell, nao, natom, ncgf, ncol, nmo, &
- nset, nsgf, numo, right, scf_step, to, width
+ isgf, ishell, iso, iw, jcol, last_mo, left, lmax, lshell, nao, natom, ncgf, ncol, nkind, &
+ nmo, nset, nsgf, numo, right, scf_step, to, width
INTEGER, DIMENSION(:), POINTER :: mo_index_range, nshell
INTEGER, DIMENSION(:, :), POINTER :: l
- LOGICAL :: ionode, my_final, my_rtp, &
- print_cartesian, print_eigvals, &
- print_eigvecs, print_occup, &
- should_output
+ LOGICAL :: ionode, my_final, my_rtp, omit_headers, print_cartesian, print_cartesian_overlap, &
+ print_eigvals, print_eigvecs, print_occup, should_output
REAL(KIND=dp) :: gap, maxocc
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: mo_eigenvalues, mo_occupation_numbers
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: cmatrix, smatrix
REAL(KIND=dp), DIMENSION(:), POINTER :: eigenvalues, occupation_numbers
TYPE(cp_fm_type), POINTER :: mo_coeff, umo_coeff
TYPE(cp_logger_type), POINTER :: logger
- TYPE(gto_basis_set_type), POINTER :: orb_basis_set
+ TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: sro
+ TYPE(gto_basis_set_p_type), DIMENSION(:), POINTER :: orb_basis_set_list
+ TYPE(gto_basis_set_type), POINTER :: orb_basis_set, orbbasis
+ TYPE(mp_para_env_type), POINTER :: para_env
+ TYPE(neighbor_list_set_p_type), DIMENSION(:), &
+ POINTER :: sro_list
TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter
+ TYPE(qs_kind_type), POINTER :: qs_kind
+ TYPE(qs_ks_env_type), POINTER :: ks_env
NULLIFY (bcgf_symbol)
NULLIFY (bsgf_symbol)
@@ -1068,6 +1085,7 @@ CONTAINS
CALL section_vals_val_get(dft_section, "PRINT%MO%CARTESIAN", l_val=print_cartesian)
CALL section_vals_val_get(dft_section, "PRINT%MO%MO_INDEX_RANGE", i_vals=mo_index_range)
CALL section_vals_val_get(dft_section, "PRINT%MO%NDIGITS", i_val=after)
+ CALL section_vals_val_get(dft_section, "PRINT%MO%CARTESIAN_OVERLAP", l_val=print_cartesian_overlap)
after = MIN(MAX(after, 1), 16)
! Do we print the final MO information after SCF convergence is reached (default: no)
@@ -1161,6 +1179,49 @@ CONTAINS
END IF
END IF
+ IF (PRESENT(qs_env)) THEN
+ IF (ASSOCIATED(qs_env) .AND. my_final .AND. print_cartesian_overlap) THEN
+ NULLIFY (qs_kind_set)
+ CALL get_qs_env(qs_env=qs_env, qs_kind_set=qs_kind_set)
+ nkind = SIZE(qs_kind_set)
+
+ IF (BTEST(cp_print_key_should_output(logger%iter_info, &
+ qs_env%input, "DFT%PRINT%AO_MATRICES/OVERLAP"), cp_p_file)) THEN
+ ALLOCATE (orb_basis_set_list(nkind))
+ DO ikind = 1, nkind
+ qs_kind => qs_kind_set(ikind)
+ NULLIFY (orb_basis_set_list(ikind)%gto_basis_set)
+ NULLIFY (orbbasis)
+ CALL get_qs_kind(qs_kind=qs_kind, basis_set=orbbasis, basis_type="ORB")
+ IF (ASSOCIATED(orbbasis)) orb_basis_set_list(ikind)%gto_basis_set => orbbasis
+ END DO
+ NULLIFY (sro_list)
+ CALL setup_neighbor_list(sro_list, orb_basis_set_list, qs_env=qs_env)
+ NULLIFY (sro)
+ NULLIFY (para_env)
+ CALL get_qs_env(qs_env, ks_env=ks_env, para_env=para_env)
+ CALL build_overlap_matrix_simple(ks_env, sro, &
+ orb_basis_set_list, orb_basis_set_list, sro_list, .TRUE.)
+ CALL release_neighbor_list_sets(sro_list)
+
+ iw = cp_print_key_unit_nr(logger, qs_env%input, "DFT%PRINT%AO_MATRICES/OVERLAP", &
+ extension=".Log")
+ CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%OMIT_HEADERS", l_val=omit_headers)
+ CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%NDIGITS", i_val=after)
+ after = MIN(MAX(after, 1), 16)
+ IF (ASSOCIATED(sro)) THEN
+ CALL cp_dbcsr_write_sparse_matrix(sro(1)%matrix, 4, after, qs_env, para_env, &
+ output_unit=iw, omit_headers=omit_headers, &
+ cartesian_basis=.TRUE.)
+ END IF
+ CALL cp_print_key_finished_output(iw, logger, qs_env%input, &
+ "DFT%PRINT%AO_MATRICES/OVERLAP")
+ IF (ASSOCIATED(sro)) CALL dbcsr_deallocate_matrix_set(sro)
+ DEALLOCATE (orb_basis_set_list)
+ END IF
+ END IF
+ END IF
+
iw = cp_print_key_unit_nr(logger, dft_section, "PRINT%MO", &
ignore_should_output=should_output, &
extension=".MOLog")
@@ -1441,10 +1502,10 @@ CONTAINS
! Release work storage
- DEALLOCATE (smatrix)
IF (print_cartesian) THEN
DEALLOCATE (cmatrix)
END IF
+ DEALLOCATE (smatrix)
ELSE IF (print_occup .OR. print_eigvals) THEN
diff --git a/src/qs_overlap.F b/src/qs_overlap.F
index ba80693bf5..3bc52c5249 100644
--- a/src/qs_overlap.F
+++ b/src/qs_overlap.F
@@ -286,11 +286,11 @@ CONTAINS
IF (dokp) THEN
CALL dbcsr_allocate_matrix_set(matrixkp_s, maxder, nimg)
CALL create_sab_matrix(ks_env, matrixkp_s, matrix_name, basis_set_list_a, basis_set_list_b, &
- sab_nl, do_symmetric)
+ sab_nl, do_symmetric, lcart=.FALSE.)
ELSE
CALL dbcsr_allocate_matrix_set(matrix_s, maxder)
CALL create_sab_matrix(ks_env, matrix_s, matrix_name, basis_set_list_a, basis_set_list_b, &
- sab_nl, do_symmetric)
+ sab_nl, do_symmetric, lcart=.FALSE.)
END IF
maxs = maxder
@@ -563,6 +563,7 @@ CONTAINS
!> \param basis_set_list_a basis set list to be used for bra in
!> \param basis_set_list_b basis set list to be used for ket in
!> \param sab_nl pair list (must be consistent with basis sets!)
+!> \param lcart ...
!> \date 11.03.2016
!> \par History
!> Simplified version of build_overlap_matrix
@@ -570,13 +571,14 @@ CONTAINS
!> \version 1.0
! **************************************************************************************************
SUBROUTINE build_overlap_matrix_simple(ks_env, matrix_s, &
- basis_set_list_a, basis_set_list_b, sab_nl)
+ basis_set_list_a, basis_set_list_b, sab_nl, lcart)
TYPE(qs_ks_env_type), POINTER :: ks_env
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
TYPE(gto_basis_set_p_type), DIMENSION(:), POINTER :: basis_set_list_a, basis_set_list_b
TYPE(neighbor_list_set_p_type), DIMENSION(:), &
POINTER :: sab_nl
+ LOGICAL, OPTIONAL :: lcart
CHARACTER(len=*), PARAMETER :: routineN = 'build_overlap_matrix_simple'
@@ -587,13 +589,14 @@ CONTAINS
INTEGER, DIMENSION(:), POINTER :: la_max, la_min, lb_max, lb_min, npgfa, &
npgfb, nsgfa, nsgfb
INTEGER, DIMENSION(:, :), POINTER :: first_sgfa, first_sgfb
- LOGICAL :: do_symmetric, found, trans
+ LOGICAL :: do_symmetric, found, ldocart, trans
REAL(KIND=dp) :: dab, rab2
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: owork
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: oint
REAL(KIND=dp), DIMENSION(3) :: rab
REAL(KIND=dp), DIMENSION(:), POINTER :: set_radius_a, set_radius_b
- REAL(KIND=dp), DIMENSION(:, :), POINTER :: rpgfa, rpgfb, scon_a, scon_b, zeta, zetb
+ REAL(KIND=dp), DIMENSION(:, :), POINTER :: ccona, cconb, rpgfa, rpgfb, scon_a, &
+ scon_b, zeta, zetb
TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
TYPE(block_p_type), ALLOCATABLE, DIMENSION(:) :: sint
TYPE(dft_control_type), POINTER :: dft_control
@@ -606,6 +609,8 @@ CONTAINS
!$ INTEGER(KIND=int_8) :: iatom8
!$ INTEGER, PARAMETER :: nlock = 501
+ ldocart = .FALSE.
+ IF (PRESENT(lcart)) ldocart = lcart
NULLIFY (dft_control)
CALL timeset(routineN, handle)
@@ -624,8 +629,13 @@ CONTAINS
nkind = SIZE(qs_kind_set)
CALL dbcsr_allocate_matrix_set(matrix_s, 1)
- CALL create_sab_matrix(ks_env, matrix_s, "Matrix", basis_set_list_a, basis_set_list_b, &
- sab_nl, do_symmetric)
+ IF (.NOT. ldocart) THEN
+ CALL create_sab_matrix(ks_env, matrix_s, "Matrix", basis_set_list_a, basis_set_list_b, &
+ sab_nl, do_symmetric, lcart=ldocart)
+ ELSE
+ CALL create_sab_matrix(ks_env, matrix_s, "Cartesian Overlap Matrix", basis_set_list_a, basis_set_list_b, &
+ sab_nl, do_symmetric, lcart=ldocart)
+ END IF
ldsab = 0
DO ikind = 1, nkind
@@ -639,12 +649,12 @@ CONTAINS
!$OMP PARALLEL DEFAULT(NONE) &
!$OMP SHARED (ldsab,sab_nl,do_symmetric,ncoset,natom,&
-!$OMP matrix_s,basis_set_list_a,basis_set_list_b,locks) &
+!$OMP matrix_s,basis_set_list_a,basis_set_list_b,locks,ldocart) &
!$OMP PRIVATE (oint,owork,sint,ikind,jkind,iatom,jatom,rab,basis_set_a,basis_set_b,&
!$OMP first_sgfa, la_max, la_min, npgfa, nsgfa, nseta, rpgfa, set_radius_a, ncoa, ncob, &
!$OMP zeta, first_sgfb, lb_max, lb_min, npgfb, nsetb, rpgfb, set_radius_b, nsgfb, dab, &
!$OMP zetb, scon_a, scon_b, irow, icol, found, trans, rab2, n1, n2, sgfa, sgfb, iset, jset, &
-!$OMP slot, lock_num, hash, hash1, hash2, iatom8 )
+!$OMP slot, lock_num, hash, hash1, hash2, iatom8, ccona, cconb )
!$OMP SINGLE
!$ ALLOCATE (locks(nlock))
@@ -684,6 +694,11 @@ CONTAINS
rpgfa => basis_set_a%pgf_radius
set_radius_a => basis_set_a%set_radius
scon_a => basis_set_a%scon
+ IF (ldocart) THEN
+ first_sgfa => basis_set_a%first_cgf
+ nsgfa => basis_set_a%ncgf_set
+ ccona => basis_set_a%ccon
+ END IF
zeta => basis_set_a%zet
! basis jkind
first_sgfb => basis_set_b%first_sgf
@@ -695,6 +710,11 @@ CONTAINS
rpgfb => basis_set_b%pgf_radius
set_radius_b => basis_set_b%set_radius
scon_b => basis_set_b%scon
+ IF (ldocart) THEN
+ first_sgfb => basis_set_b%first_cgf
+ nsgfb => basis_set_b%ncgf_set
+ cconb => basis_set_b%ccon
+ END IF
zetb => basis_set_b%zet
IF (do_symmetric) THEN
@@ -741,8 +761,13 @@ CONTAINS
lb_max(jset), lb_min(jset), npgfb(jset), rpgfb(:, jset), zetb(:, jset), &
rab, sab=oint(:, :, 1))
! Contraction
- CALL contraction(oint(:, :, 1), owork, ca=scon_a(:, sgfa:), na=n1, ma=nsgfa(iset), &
- cb=scon_b(:, sgfb:), nb=n2, mb=nsgfb(jset), fscale=1.0_dp, trans=trans)
+ IF (.NOT. ldocart) THEN
+ CALL contraction(oint(:, :, 1), owork, ca=scon_a(:, sgfa:), na=n1, ma=nsgfa(iset), &
+ cb=scon_b(:, sgfb:), nb=n2, mb=nsgfb(jset), fscale=1.0_dp, trans=trans)
+ ELSE
+ CALL contraction(oint(:, :, 1), owork, ca=ccona(:, sgfa:), na=n1, ma=nsgfa(iset), &
+ cb=cconb(:, sgfb:), nb=n2, mb=nsgfb(jset), fscale=1.0_dp, trans=trans)
+ END IF
!$OMP CRITICAL(blockadd)
CALL block_add("IN", owork, nsgfa(iset), nsgfb(jset), sint(1)%block, &
sgfa, sgfb, trans=trans)
@@ -1028,9 +1053,10 @@ CONTAINS
!> \param basis_set_list_b Basis set used for |b>
!> \param sab_nl Overlap neighbor list
!> \param symmetric Is symmetry used in the neighbor list?
+!> \param lcart ...
! **************************************************************************************************
SUBROUTINE create_sab_matrix_1d(ks_env, matrix_s, matrix_name, &
- basis_set_list_a, basis_set_list_b, sab_nl, symmetric)
+ basis_set_list_a, basis_set_list_b, sab_nl, symmetric, lcart)
TYPE(qs_ks_env_type), POINTER :: ks_env
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
@@ -1039,6 +1065,7 @@ CONTAINS
TYPE(neighbor_list_set_p_type), DIMENSION(:), &
POINTER :: sab_nl
LOGICAL, INTENT(IN) :: symmetric
+ LOGICAL, INTENT(IN), OPTIONAL :: lcart
CHARACTER(LEN=12) :: cgfsym
CHARACTER(LEN=32) :: symmetry_string
@@ -1049,10 +1076,14 @@ CONTAINS
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
+ LOGICAL:: my_lcart
+
CALL get_ks_env(ks_env=ks_env, particle_set=particle_set, &
qs_kind_set=qs_kind_set, dbcsr_dist=dbcsr_dist)
natom = SIZE(particle_set)
+ my_lcart = .FALSE.
+ IF (PRESENT(lcart)) my_lcart = lcart
IF (PRESENT(matrix_name)) THEN
mname = matrix_name
@@ -1064,10 +1095,17 @@ CONTAINS
ALLOCATE (row_blk_sizes(natom), col_blk_sizes(natom))
- CALL get_particle_set(particle_set, qs_kind_set, nsgf=row_blk_sizes, &
- basis=basis_set_list_a)
- CALL get_particle_set(particle_set, qs_kind_set, nsgf=col_blk_sizes, &
- basis=basis_set_list_b)
+ IF (.NOT. my_lcart) THEN
+ CALL get_particle_set(particle_set, qs_kind_set, nsgf=row_blk_sizes, &
+ basis=basis_set_list_a)
+ CALL get_particle_set(particle_set, qs_kind_set, nsgf=col_blk_sizes, &
+ basis=basis_set_list_b)
+ ELSE
+ CALL get_particle_set(particle_set, qs_kind_set, ncgf=row_blk_sizes, &
+ basis=basis_set_list_a)
+ CALL get_particle_set(particle_set, qs_kind_set, ncgf=col_blk_sizes, &
+ basis=basis_set_list_b)
+ END IF
! prepare for allocation
IF (symmetric) THEN
@@ -1118,9 +1156,10 @@ CONTAINS
!> \param basis_set_list_b Basis set used for |b>
!> \param sab_nl Overlap neighbor list
!> \param symmetric Is symmetry used in the neighbor list?
+!> \param lcart ...
! **************************************************************************************************
SUBROUTINE create_sab_matrix_2d(ks_env, matrix_s, matrix_name, &
- basis_set_list_a, basis_set_list_b, sab_nl, symmetric)
+ basis_set_list_a, basis_set_list_b, sab_nl, symmetric, lcart)
TYPE(qs_ks_env_type), POINTER :: ks_env
TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s
@@ -1129,12 +1168,14 @@ CONTAINS
TYPE(neighbor_list_set_p_type), DIMENSION(:), &
POINTER :: sab_nl
LOGICAL, INTENT(IN) :: symmetric
+ LOGICAL, INTENT(in), OPTIONAL :: lcart
CHARACTER(LEN=12) :: cgfsym
CHARACTER(LEN=32) :: symmetry_string
CHARACTER(LEN=default_string_length) :: mname, name
INTEGER :: i1, i2, natom
INTEGER, DIMENSION(:), POINTER :: col_blk_sizes, row_blk_sizes
+ LOGICAL :: my_lcart
TYPE(dbcsr_distribution_type), POINTER :: dbcsr_dist
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
@@ -1143,6 +1184,8 @@ CONTAINS
qs_kind_set=qs_kind_set, dbcsr_dist=dbcsr_dist)
natom = SIZE(particle_set)
+ my_lcart = .FALSE.
+ IF (PRESENT(lcart)) my_lcart = lcart
IF (PRESENT(matrix_name)) THEN
mname = matrix_name
@@ -1152,10 +1195,17 @@ CONTAINS
ALLOCATE (row_blk_sizes(natom), col_blk_sizes(natom))
- CALL get_particle_set(particle_set, qs_kind_set, nsgf=row_blk_sizes, &
- basis=basis_set_list_a)
- CALL get_particle_set(particle_set, qs_kind_set, nsgf=col_blk_sizes, &
- basis=basis_set_list_b)
+ IF (.NOT. my_lcart) THEN
+ CALL get_particle_set(particle_set, qs_kind_set, nsgf=row_blk_sizes, &
+ basis=basis_set_list_a)
+ CALL get_particle_set(particle_set, qs_kind_set, nsgf=col_blk_sizes, &
+ basis=basis_set_list_b)
+ ELSE
+ CALL get_particle_set(particle_set, qs_kind_set, ncgf=row_blk_sizes, &
+ basis=basis_set_list_a)
+ CALL get_particle_set(particle_set, qs_kind_set, ncgf=col_blk_sizes, &
+ basis=basis_set_list_b)
+ END IF
! prepare for allocation
IF (symmetric) THEN
diff --git a/src/qs_scf_output.F b/src/qs_scf_output.F
index 8a4f61fede..292cd37b00 100644
--- a/src/qs_scf_output.F
+++ b/src/qs_scf_output.F
@@ -220,6 +220,7 @@ CONTAINS
TYPE(mp_para_env_type), POINTER :: para_env
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
TYPE(preconditioner_type), POINTER :: local_preconditioner
+ TYPE(qs_environment_type), POINTER :: cart_overlap_qs_env
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
TYPE(scf_control_type), POINTER :: scf_control
TYPE(section_vals_type), POINTER :: dft_section, input
@@ -434,6 +435,8 @@ CONTAINS
END IF ! OT is used
! Print MO information
+ NULLIFY (cart_overlap_qs_env)
+ IF ((ikp == 1) .AND. (ispin == 1)) cart_overlap_qs_env => qs_env
IF (nspin > 1) THEN
SELECT CASE (ispin)
CASE (1)
@@ -446,19 +449,21 @@ CONTAINS
IF (ASSOCIATED(umo_set)) THEN
CALL write_mo_set_to_output_unit(mo_set, qs_kind_set, particle_set, dft_section, 4, kpoint, &
final_mos=final_mos, spin=TRIM(spin), solver_method=solver_method, &
- umo_set=umo_set)
+ umo_set=umo_set, qs_env=cart_overlap_qs_env)
ELSE
CALL write_mo_set_to_output_unit(mo_set, qs_kind_set, particle_set, dft_section, 4, kpoint, &
- final_mos=final_mos, spin=TRIM(spin), solver_method=solver_method)
+ final_mos=final_mos, spin=TRIM(spin), solver_method=solver_method, &
+ qs_env=cart_overlap_qs_env)
END IF
ELSE
IF (ASSOCIATED(umo_set)) THEN
CALL write_mo_set_to_output_unit(mo_set, qs_kind_set, particle_set, dft_section, 4, kpoint, &
final_mos=final_mos, solver_method=solver_method, &
- umo_set=umo_set)
+ umo_set=umo_set, qs_env=cart_overlap_qs_env)
ELSE
CALL write_mo_set_to_output_unit(mo_set, qs_kind_set, particle_set, dft_section, 4, kpoint, &
- final_mos=final_mos, solver_method=solver_method)
+ final_mos=final_mos, solver_method=solver_method, &
+ qs_env=cart_overlap_qs_env)
END IF
END IF
diff --git a/tests/QS/regtest-gpw-1/TEST_FILES.toml b/tests/QS/regtest-gpw-1/TEST_FILES.toml
index be210bb4ce..8356dca72e 100644
--- a/tests/QS/regtest-gpw-1/TEST_FILES.toml
+++ b/tests/QS/regtest-gpw-1/TEST_FILES.toml
@@ -72,4 +72,6 @@
"h2q.inp" = [{matcher="E_total", tol=1.0E-10, ref=-0.76162210587786}]
# Mol Dipole Voronoi
"moldip_voronoi.inp" = [{matcher="E_total", tol=1.0E-10, ref=-41.88176215391840}]
+# Test printing of cartesian overlap matrix with header
+"h2test.inp" = [{matcher="E_total", tol=3e-13, ref=-1.06531398809803}]
#EOF
diff --git a/tests/QS/regtest-gpw-1/h2test.inp b/tests/QS/regtest-gpw-1/h2test.inp
new file mode 100644
index 0000000000..d2802d479c
--- /dev/null
+++ b/tests/QS/regtest-gpw-1/h2test.inp
@@ -0,0 +1,73 @@
+&GLOBAL
+ PRINT_LEVEL MEDIUM
+ PROJECT h2
+ RUN_TYPE ENERGY
+&END GLOBAL
+
+&FORCE_EVAL
+ METHOD Quickstep
+ &DFT
+ BASIS_SET_FILE_NAME BASIS_MOLOPT
+ BASIS_SET_FILE_NAME EMSL_BASIS_SETS
+ BASIS_SET_FILE_NAME BASIS_MOLOPT_UCL
+ &MGRID
+ CUTOFF 400
+ REL_CUTOFF 60
+ &END MGRID
+ &POISSON
+ PERIODIC NONE
+ PSOLVER WAVELET
+ &END POISSON
+ &PRINT
+ &AO_MATRICES
+ OVERLAP
+ &END AO_MATRICES
+ &MO
+ CARTESIAN .TRUE.
+ CARTESIAN_OVERLAP .TRUE.
+ COEFFICIENTS .True.
+ FILENAME ./mos
+ MO_INDEX_RANGE 0 200
+ &END MO
+ &END PRINT
+ &QS
+ EPS_DEFAULT 1.0E-14
+ LMAXN0 4
+ METHOD gapw
+ &END QS
+ &SCF
+ ADDED_MOS -1
+ EPS_SCF 1.0E-6
+ MAX_SCF 300
+ SCF_GUESS RESTART
+ &DIAGONALIZATION
+ &END DIAGONALIZATION
+ &END SCF
+ &XC
+ &XC_FUNCTIONAL pbe
+ &END XC_FUNCTIONAL
+ &END XC
+ &END DFT
+ &SUBSYS
+ &CELL
+ ABC 6. 6. 6.
+ PERIODIC NONE
+ &END CELL
+ &COORD
+ h 0.00000000000000 0.00000000000000 0.58312867712598
+ h 0.00000000000000 0.00000000000000 0.08687132338553
+ &END COORD
+ &KIND H
+ BASIS_SET Ahlrichs-VDZ
+ ELEMENT H
+ HARD_EXP_RADIUS 1.00
+ LEBEDEV_GRID 50
+ POTENTIAL ALL
+ RADIAL_GRID 50
+ &END KIND
+ &TOPOLOGY
+ &CENTER_COORDINATES T
+ &END CENTER_COORDINATES
+ &END TOPOLOGY
+ &END SUBSYS
+&END FORCE_EVAL