26 USE dbcsr_api,
ONLY: dbcsr_create,&
27 dbcsr_distribution_type,&
33 dbcsr_type_no_symmetry
56#include "./base/base_uses.f90"
63 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'optbas_opt_utils'
80 fval, energy, S_cond_number)
81 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(IN) :: mos, mos_aux_fit
82 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_ks
83 TYPE(dbcsr_type),
POINTER :: q, snew
85 REAL(kind=
dp) :: fval, energy, s_cond_number
87 CHARACTER(len=*),
PARAMETER :: routinen =
'evaluate_optvals'
89 INTEGER :: handle, ispin, iunit, naux, nmo, norb, &
91 INTEGER,
DIMENSION(:),
POINTER :: col_blk_sizes, row_blk_sizes
92 REAL(kind=
dp) :: tmp_energy, trace
93 REAL(kind=
dp),
DIMENSION(2) :: condnum
96 TYPE(
cp_fm_type),
POINTER :: mo_coeff, mo_coeff_aux_fit
97 TYPE(dbcsr_distribution_type) :: dbcsr_dist
98 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: smat
99 TYPE(dbcsr_type) :: qt
101 CALL timeset(routinen, handle)
105 NULLIFY (col_blk_sizes, row_blk_sizes)
106 CALL dbcsr_get_info(q, distribution=dbcsr_dist, &
107 nfullrows_total=naux, nfullcols_total=norb, &
108 row_blk_size=row_blk_sizes, col_blk_size=col_blk_sizes)
109 CALL dbcsr_create(matrix=qt, name=
"Qt", &
110 dist=dbcsr_dist, matrix_type=dbcsr_type_no_symmetry, &
111 row_blk_size=col_blk_sizes, col_blk_size=row_blk_sizes, &
113 CALL dbcsr_transposed(qt, q)
118 CALL get_mo_set(mos(ispin), mo_coeff=mo_coeff)
119 CALL get_mo_set(mos_aux_fit(ispin), mo_coeff=mo_coeff_aux_fit)
121 CALL cp_fm_create(tmp1, matrix_struct=mo_coeff%matrix_struct)
124 fval = fval - 2.0_dp*trace + 2.0_dp*nmo
126 CALL cp_fm_create(tmp2, matrix_struct=mo_coeff%matrix_struct)
127 CALL parallel_gemm(
'N',
'N', norb, nmo, norb, 1.0_dp, s_inv_orb, tmp1, 0.0_dp, tmp2)
130 energy = energy + tmp_energy*(3.0_dp - real(nspins, kind=
dp))
134 CALL dbcsr_release(qt)
136 ALLOCATE (smat(1, 1))
137 smat(1, 1)%matrix => snew
140 CALL overlap_condnum(smat, condnum, iunit, .false., .true., .false., blacs_env)
141 s_cond_number = condnum(2)
144 CALL timestop(handle)
156 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: saux, sauxorb
157 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(IN) :: mos, mosaux
159 CHARACTER(len=*),
PARAMETER :: routinen =
'fit_mo_coeffs'
160 REAL(kind=
dp),
PARAMETER :: threshold = 1.e-12_dp
162 INTEGER :: handle, ispin, naux, ndep, nmo, norb, &
165 TYPE(
cp_fm_type) :: fm_s, fm_sinv, tmat, tmp1, tmp2, work
166 TYPE(
cp_fm_type),
POINTER :: mo_coeff, mo_coeff_aux
168 CALL timeset(routinen, handle)
170 CALL dbcsr_get_info(saux(1)%matrix, nfullrows_total=naux)
171 CALL dbcsr_get_info(sauxorb(1)%matrix, nfullcols_total=norb)
175 context=mo_coeff%matrix_struct%context, &
176 para_env=mo_coeff%matrix_struct%para_env)
185 CALL get_mo_set(mos(ispin), mo_coeff=mo_coeff)
186 CALL get_mo_set(mosaux(ispin), mo_coeff=mo_coeff_aux)
188 CALL cp_fm_create(tmp1, matrix_struct=mo_coeff_aux%matrix_struct)
189 CALL cp_fm_create(tmp2, matrix_struct=mo_coeff_aux%matrix_struct)
191 context=mo_coeff%matrix_struct%context, &
192 para_env=mo_coeff%matrix_struct%para_env)
198 CALL parallel_gemm(
'N',
'N', naux, nmo, naux, 1.0_dp, fm_sinv, tmp1, 0.0_dp, tmp2)
199 CALL parallel_gemm(
'T',
'N', nmo, nmo, naux, 1.0_dp, tmp1, tmp2, 0.0_dp, tmat)
200 CALL cp_fm_power(tmat, work, -0.5_dp, threshold, ndep)
201 CALL parallel_gemm(
'N',
'N', naux, nmo, nmo, 1.0_dp, tmp2, tmat, 0.0_dp, mo_coeff_aux)
210 CALL timestop(handle)
226 POINTER :: sab_aux, sab_aux_orb
227 CHARACTER(*) :: basis_type
229 CHARACTER(LEN=*),
PARAMETER :: routinen =
'optbas_build_neighborlist'
231 INTEGER :: handle, ikind, nkind
232 LOGICAL :: mic, molecule_only
233 LOGICAL,
ALLOCATABLE,
DIMENSION(:) :: aux_fit_present, orb_present
235 REAL(
dp),
ALLOCATABLE,
DIMENSION(:) :: aux_fit_radius, orb_radius
236 REAL(
dp),
ALLOCATABLE,
DIMENSION(:, :) :: pair_radius
246 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
249 CALL timeset(routinen, handle)
253 molecule_only = .false.
258 atomic_kind_set=atomic_kind_set, &
259 qs_kind_set=qs_kind_set, &
261 distribution_2d=distribution_2d, &
262 molecule_set=molecule_set, &
263 local_particles=distribution_1d, &
264 particle_set=particle_set, &
270 nkind =
SIZE(atomic_kind_set)
271 ALLOCATE (orb_radius(nkind), aux_fit_radius(nkind))
272 orb_radius(:) = 0.0_dp
273 aux_fit_radius(:) = 0.0_dp
274 ALLOCATE (orb_present(nkind), aux_fit_present(nkind))
275 ALLOCATE (pair_radius(nkind, nkind))
276 ALLOCATE (atom2d(nkind))
278 CALL atom2d_build(atom2d, distribution_1d, distribution_2d, atomic_kind_set, &
279 molecule_set, molecule_only, particle_set=particle_set)
282 CALL get_atomic_kind(atomic_kind_set(ikind), atom_list=atom2d(ikind)%list)
283 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type=
"ORB")
284 IF (
ASSOCIATED(orb_basis_set))
THEN
285 orb_present(ikind) = .true.
286 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, kind_radius=orb_radius(ikind))
288 orb_present(ikind) = .false.
289 orb_radius(ikind) = 0.0_dp
291 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_fit_basis_set, basis_type=basis_type)
292 IF (
ASSOCIATED(aux_fit_basis_set))
THEN
293 aux_fit_present(ikind) = .true.
294 CALL get_gto_basis_set(gto_basis_set=aux_fit_basis_set, kind_radius=aux_fit_radius(ikind))
296 aux_fit_present(ikind) = .false.
297 aux_fit_radius(ikind) = 0.0_dp
301 CALL pair_radius_setup(aux_fit_present, aux_fit_present, aux_fit_radius, aux_fit_radius, pair_radius)
303 mic=mic, molecular=molecule_only, subcells=subcells, nlname=
"sab_aux")
304 CALL pair_radius_setup(aux_fit_present, orb_present, aux_fit_radius, orb_radius, pair_radius)
306 mic=mic, symmetric=.false., molecular=molecule_only, subcells=subcells, &
307 nlname=
"sab_aux_orb")
312 DEALLOCATE (orb_present, aux_fit_present)
313 DEALLOCATE (orb_radius, aux_fit_radius)
314 DEALLOCATE (pair_radius)
316 CALL timestop(handle)
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
subroutine, public 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)
...
Handles all functions related to the CELL.
methods related to the blacs parallel environment
DBCSR operations in CP2K.
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
basic linear algebra operations for full matrices
subroutine, public cp_fm_invert(matrix_a, matrix_inverse, det_a, eps_svd, eigval)
Inverts a cp_fm_type matrix, optionally returning the determinant of the input matrix.
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
subroutine, public cp_fm_power(matrix, work, exponent, threshold, n_dependent, verbose, eigvals)
...
represent the structure of a full matrix
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
represent a full matrix distributed on many processors
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
stores a mapping of 2D info (e.g. matrix) on a 2D processor distribution (i.e. blacs grid) where cpus...
Defines the basic variable types.
integer, parameter, public dp
Interface to the message passing library MPI.
Define the data structure for the molecule information.
subroutine, public optbas_build_neighborlist(qs_env, sab_aux, sab_aux_orb, basis_type)
rebuilds neighborlist for absis sets
subroutine, public evaluate_optvals(mos, mos_aux_fit, matrix_ks, q, snew, s_inv_orb, fval, energy, s_cond_number)
...
subroutine, public fit_mo_coeffs(saux, sauxorb, mos, mosaux)
...
basic linear algebra operations for full matrixes
Define the data structure for the particle information.
Calculation of overlap matrix condition numbers.
subroutine, public overlap_condnum(matrixkp_s, condnum, iunit, norml1, norml2, use_arnoldi, blacs_env)
Calculation of the overlap matrix Condition Number.
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, ecoul_1c, rho0_s_rs, rho0_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Get the QUICKSTEP environment.
Define the quickstep kind type and their sub types.
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_r3d_rs_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition and initialisation of the mo data type.
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count)
Get the components of a MO set data structure.
Define the neighbor list data types and the corresponding functionality.
Generate the atomic neighbor lists.
subroutine, public atom2d_cleanup(atom2d)
free the internals of atom2d
subroutine, public pair_radius_setup(present_a, present_b, radius_a, radius_b, pair_radius, prmin)
...
subroutine, public build_neighbor_lists(ab_list, particle_set, atom, cell, pair_radius, subcells, mic, symmetric, molecular, subset_of_mol, current_subset, operator_type, nlname, atomb_to_keep)
Build simple pair neighbor lists.
subroutine, public atom2d_build(atom2d, distribution_1d, distribution_2d, atomic_kind_set, molecule_set, molecule_only, particle_set)
Build some distribution structure of atoms, refactored from build_qs_neighbor_lists.
Provides all information about an atomic kind.
Type defining parameters related to the simulation cell.
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
keeps the information about the structure of a full matrix
structure to store local (to a processor) ordered lists of integers.
distributes pairs on a 2d grid of processors
stores all the informations relevant to an mpi environment
Provides all information about a quickstep kind.
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...