(git:e5b1968)
Loading...
Searching...
No Matches
qs_mo_methods.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief collects routines that perform operations directly related to MOs
10!> \note
11!> first version : most routines imported
12!> \author Joost VandeVondele (2003-08)
13! **************************************************************************************************
15 USE admm_types, ONLY: admm_type
20 USE cp_dbcsr_api, ONLY: dbcsr_copy,&
27 dbcsr_type_no_symmetry
28 USE cp_dbcsr_diag, ONLY: cp_dbcsr_syevd,&
43 USE cp_fm_types, ONLY: cp_fm_create,&
49 USE kinds, ONLY: dp
52 USE physcon, ONLY: evolt
54 USE qs_mo_types, ONLY: get_mo_set,&
57#include "./base/base_uses.f90"
58
59 IMPLICIT NONE
60 PRIVATE
61 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_mo_methods'
62
66
68 MODULE PROCEDURE subspace_eigenvalues_ks_fm
69 MODULE PROCEDURE subspace_eigenvalues_ks_dbcsr
70 END INTERFACE
71
72 INTERFACE make_basis_sv
73 MODULE PROCEDURE make_basis_sv_fm
74 MODULE PROCEDURE make_basis_sv_dbcsr
75 END INTERFACE
76
77CONTAINS
78
79! **************************************************************************************************
80!> \brief returns an S-orthonormal basis v (v^T S v ==1)
81!> \param vmatrix ...
82!> \param ncol ...
83!> \param matrix_s ...
84!> \par History
85!> 03.2006 created [Joost VandeVondele]
86! **************************************************************************************************
87 SUBROUTINE make_basis_sm(vmatrix, ncol, matrix_s)
88 TYPE(cp_fm_type), INTENT(IN) :: vmatrix
89 INTEGER, INTENT(IN) :: ncol
90 TYPE(dbcsr_type), POINTER :: matrix_s
91
92 CHARACTER(LEN=*), PARAMETER :: routinen = 'make_basis_sm'
93 REAL(kind=dp), PARAMETER :: rone = 1.0_dp, rzero = 0.0_dp
94
95 INTEGER :: handle, n, ncol_global
96 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
97 TYPE(cp_fm_type) :: overlap_vv, svmatrix
98
99 IF (ncol .EQ. 0) RETURN
100
101 CALL timeset(routinen, handle)
102
103 CALL cp_fm_get_info(matrix=vmatrix, nrow_global=n, ncol_global=ncol_global)
104 IF (ncol .GT. ncol_global) cpabort("Wrong ncol value")
105
106 CALL cp_fm_create(svmatrix, vmatrix%matrix_struct, "SV")
107 CALL cp_dbcsr_sm_fm_multiply(matrix_s, vmatrix, svmatrix, ncol)
108
109 NULLIFY (fm_struct_tmp)
110 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=ncol, ncol_global=ncol, &
111 para_env=vmatrix%matrix_struct%para_env, &
112 context=vmatrix%matrix_struct%context)
113 CALL cp_fm_create(overlap_vv, fm_struct_tmp, "overlap_vv")
114 CALL cp_fm_struct_release(fm_struct_tmp)
115
116 CALL parallel_gemm('T', 'N', ncol, ncol, n, rone, vmatrix, svmatrix, rzero, overlap_vv)
117 CALL cp_fm_cholesky_decompose(overlap_vv)
118 CALL cp_fm_triangular_multiply(overlap_vv, vmatrix, n_cols=ncol, side='R', invert_tr=.true.)
119
120 CALL cp_fm_release(overlap_vv)
121 CALL cp_fm_release(svmatrix)
122
123 CALL timestop(handle)
124
125 END SUBROUTINE make_basis_sm
126
127! **************************************************************************************************
128!> \brief returns an S-orthonormal basis v and the corresponding matrix S*v as well
129!> \param vmatrix ...
130!> \param ncol ...
131!> \param svmatrix ...
132!> \par History
133!> 03.2006 created [Joost VandeVondele]
134! **************************************************************************************************
135 SUBROUTINE make_basis_sv_fm(vmatrix, ncol, svmatrix)
136
137 TYPE(cp_fm_type), INTENT(IN) :: vmatrix
138 INTEGER, INTENT(IN) :: ncol
139 TYPE(cp_fm_type), INTENT(IN) :: svmatrix
140
141 CHARACTER(LEN=*), PARAMETER :: routineN = 'make_basis_sv_fm'
142 REAL(KIND=dp), PARAMETER :: rone = 1.0_dp, rzero = 0.0_dp
143
144 INTEGER :: handle, n, ncol_global
145 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
146 TYPE(cp_fm_type) :: overlap_vv
147
148 IF (ncol .EQ. 0) RETURN
149
150 CALL timeset(routinen, handle)
151 NULLIFY (fm_struct_tmp)
152
153 CALL cp_fm_get_info(matrix=vmatrix, nrow_global=n, ncol_global=ncol_global)
154 IF (ncol .GT. ncol_global) cpabort("Wrong ncol value")
155
156 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=ncol, ncol_global=ncol, &
157 para_env=vmatrix%matrix_struct%para_env, &
158 context=vmatrix%matrix_struct%context)
159 CALL cp_fm_create(overlap_vv, fm_struct_tmp, "overlap_vv")
160 CALL cp_fm_struct_release(fm_struct_tmp)
161
162 CALL parallel_gemm('T', 'N', ncol, ncol, n, rone, vmatrix, svmatrix, rzero, overlap_vv)
163 CALL cp_fm_cholesky_decompose(overlap_vv)
164 CALL cp_fm_triangular_multiply(overlap_vv, vmatrix, n_cols=ncol, side='R', invert_tr=.true.)
165 CALL cp_fm_triangular_multiply(overlap_vv, svmatrix, n_cols=ncol, side='R', invert_tr=.true.)
166
167 CALL cp_fm_release(overlap_vv)
168
169 CALL timestop(handle)
170
171 END SUBROUTINE make_basis_sv_fm
172
173! **************************************************************************************************
174!> \brief ...
175!> \param vmatrix ...
176!> \param ncol ...
177!> \param svmatrix ...
178!> \param para_env ...
179!> \param blacs_env ...
180! **************************************************************************************************
181 SUBROUTINE make_basis_sv_dbcsr(vmatrix, ncol, svmatrix, para_env, blacs_env)
182
183 TYPE(dbcsr_type) :: vmatrix
184 INTEGER, INTENT(IN) :: ncol
185 TYPE(dbcsr_type) :: svmatrix
186 TYPE(mp_para_env_type), POINTER :: para_env
187 TYPE(cp_blacs_env_type), POINTER :: blacs_env
188
189 CHARACTER(LEN=*), PARAMETER :: routineN = 'make_basis_sv_dbcsr'
190 REAL(KIND=dp), PARAMETER :: rone = 1.0_dp, rzero = 0.0_dp
191
192 INTEGER :: handle, n, ncol_global
193 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
194 TYPE(cp_fm_type) :: fm_svmatrix, fm_vmatrix, overlap_vv
195
196 IF (ncol .EQ. 0) RETURN
197
198 CALL timeset(routinen, handle)
199
200 !CALL cp_fm_get_info(matrix=vmatrix,nrow_global=n,ncol_global=ncol_global)
201 CALL dbcsr_get_info(vmatrix, nfullrows_total=n, nfullcols_total=ncol_global)
202 IF (ncol .GT. ncol_global) cpabort("Wrong ncol value")
203
204 CALL cp_fm_struct_create(fm_struct_tmp, context=blacs_env, nrow_global=ncol, &
205 ncol_global=ncol, para_env=para_env)
206 CALL cp_fm_create(overlap_vv, fm_struct_tmp, name="fm_overlap_vv")
207 CALL cp_fm_struct_release(fm_struct_tmp)
208
209 CALL cp_fm_struct_create(fm_struct_tmp, context=blacs_env, nrow_global=n, &
210 ncol_global=ncol_global, para_env=para_env)
211 CALL cp_fm_create(fm_vmatrix, fm_struct_tmp, name="fm_vmatrix")
212 CALL cp_fm_create(fm_svmatrix, fm_struct_tmp, name="fm_svmatrix")
213 CALL cp_fm_struct_release(fm_struct_tmp)
214
215 CALL copy_dbcsr_to_fm(vmatrix, fm_vmatrix)
216 CALL copy_dbcsr_to_fm(svmatrix, fm_svmatrix)
217
218 CALL parallel_gemm('T', 'N', ncol, ncol, n, rone, fm_vmatrix, fm_svmatrix, rzero, overlap_vv)
219 CALL cp_fm_cholesky_decompose(overlap_vv)
220 CALL cp_fm_triangular_multiply(overlap_vv, fm_vmatrix, n_cols=ncol, side='R', invert_tr=.true.)
221 CALL cp_fm_triangular_multiply(overlap_vv, fm_svmatrix, n_cols=ncol, side='R', invert_tr=.true.)
222
223 CALL copy_fm_to_dbcsr(fm_vmatrix, vmatrix)
224 CALL copy_fm_to_dbcsr(fm_svmatrix, svmatrix)
225
226 CALL cp_fm_release(overlap_vv)
227 CALL cp_fm_release(fm_vmatrix)
228 CALL cp_fm_release(fm_svmatrix)
229
230 CALL timestop(handle)
231
232 END SUBROUTINE make_basis_sv_dbcsr
233
234! **************************************************************************************************
235!> \brief return a set of S orthonormal vectors (C^T S C == 1) where
236!> the cholesky decomposed form of S is passed as an argument
237!> \param vmatrix ...
238!> \param ncol ...
239!> \param ortho cholesky decomposed S matrix
240!> \par History
241!> 03.2006 created [Joost VandeVondele]
242!> \note
243!> if the cholesky decomposed S matrix is not available
244!> use make_basis_sm since this is much faster than computing the
245!> cholesky decomposition of S
246! **************************************************************************************************
247 SUBROUTINE make_basis_cholesky(vmatrix, ncol, ortho)
248
249 TYPE(cp_fm_type), INTENT(IN) :: vmatrix
250 INTEGER, INTENT(IN) :: ncol
251 TYPE(cp_fm_type), INTENT(IN) :: ortho
252
253 CHARACTER(LEN=*), PARAMETER :: routinen = 'make_basis_cholesky'
254 REAL(kind=dp), PARAMETER :: rone = 1.0_dp, rzero = 0.0_dp
255
256 INTEGER :: handle, n, ncol_global
257 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
258 TYPE(cp_fm_type) :: overlap_vv
259
260 IF (ncol .EQ. 0) RETURN
261
262 CALL timeset(routinen, handle)
263 NULLIFY (fm_struct_tmp)
264
265 CALL cp_fm_get_info(matrix=vmatrix, nrow_global=n, ncol_global=ncol_global)
266 IF (ncol .GT. ncol_global) cpabort("Wrong ncol value")
267
268 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=ncol, ncol_global=ncol, &
269 para_env=vmatrix%matrix_struct%para_env, &
270 context=vmatrix%matrix_struct%context)
271 CALL cp_fm_create(overlap_vv, fm_struct_tmp, "overlap_vv")
272 CALL cp_fm_struct_release(fm_struct_tmp)
273
274 CALL cp_fm_triangular_multiply(ortho, vmatrix, n_cols=ncol)
275 CALL cp_fm_syrk('U', 'T', n, rone, vmatrix, 1, 1, rzero, overlap_vv)
276 CALL cp_fm_cholesky_decompose(overlap_vv)
277 CALL cp_fm_triangular_multiply(overlap_vv, vmatrix, n_cols=ncol, side='R', invert_tr=.true.)
278 CALL cp_fm_triangular_multiply(ortho, vmatrix, n_cols=ncol, invert_tr=.true.)
279
280 CALL cp_fm_release(overlap_vv)
281
282 CALL timestop(handle)
283
284 END SUBROUTINE make_basis_cholesky
285
286! **************************************************************************************************
287!> \brief return a set of S orthonormal vectors (C^T S C == 1) where
288!> a Loedwin transformation is applied to keep the rotated vectors as close
289!> as possible to the original ones
290!> \param vmatrix ...
291!> \param ncol ...
292!> \param matrix_s ...
293!> \param
294!> \par History
295!> 05.2009 created [MI]
296!> \note
297! **************************************************************************************************
298 SUBROUTINE make_basis_lowdin(vmatrix, ncol, matrix_s)
299
300 TYPE(cp_fm_type), INTENT(IN) :: vmatrix
301 INTEGER, INTENT(IN) :: ncol
302 TYPE(dbcsr_type) :: matrix_s
303
304 CHARACTER(LEN=*), PARAMETER :: routinen = 'make_basis_lowdin'
305 REAL(kind=dp), PARAMETER :: rone = 1.0_dp, rzero = 0.0_dp
306
307 INTEGER :: handle, n, ncol_global, ndep
308 REAL(dp) :: threshold
309 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
310 TYPE(cp_fm_type) :: csc, sc, work
311
312 IF (ncol .EQ. 0) RETURN
313
314 CALL timeset(routinen, handle)
315 NULLIFY (fm_struct_tmp)
316 threshold = 1.0e-7_dp
317 CALL cp_fm_get_info(matrix=vmatrix, nrow_global=n, ncol_global=ncol_global)
318 IF (ncol .GT. ncol_global) cpabort("Wrong ncol value")
319
320 CALL cp_fm_create(sc, vmatrix%matrix_struct, "SC")
321 CALL cp_dbcsr_sm_fm_multiply(matrix_s, vmatrix, sc, ncol)
322
323 NULLIFY (fm_struct_tmp)
324 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=ncol, ncol_global=ncol, &
325 para_env=vmatrix%matrix_struct%para_env, &
326 context=vmatrix%matrix_struct%context)
327 CALL cp_fm_create(csc, fm_struct_tmp, "csc")
328 CALL cp_fm_create(work, fm_struct_tmp, "work")
329 CALL cp_fm_struct_release(fm_struct_tmp)
330
331 CALL parallel_gemm('T', 'N', ncol, ncol, n, rone, vmatrix, sc, rzero, csc)
332 CALL cp_fm_power(csc, work, -0.5_dp, threshold, ndep)
333 CALL parallel_gemm('N', 'N', n, ncol, ncol, rone, vmatrix, csc, rzero, sc)
334 CALL cp_fm_to_fm(sc, vmatrix, ncol, 1, 1)
335
336 CALL cp_fm_release(csc)
337 CALL cp_fm_release(sc)
338 CALL cp_fm_release(work)
339
340 CALL timestop(handle)
341
342 END SUBROUTINE make_basis_lowdin
343
344! **************************************************************************************************
345!> \brief given a set of vectors, return an orthogonal (C^T C == 1) set
346!> spanning the same space (notice, only for cases where S==1)
347!> \param vmatrix ...
348!> \param ncol ...
349!> \par History
350!> 03.2006 created [Joost VandeVondele]
351! **************************************************************************************************
352 SUBROUTINE make_basis_simple(vmatrix, ncol)
353
354 TYPE(cp_fm_type), INTENT(IN) :: vmatrix
355 INTEGER, INTENT(IN) :: ncol
356
357 CHARACTER(LEN=*), PARAMETER :: routinen = 'make_basis_simple'
358 REAL(kind=dp), PARAMETER :: rone = 1.0_dp, rzero = 0.0_dp
359
360 INTEGER :: handle, n, ncol_global
361 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
362 TYPE(cp_fm_type) :: overlap_vv
363
364 IF (ncol .EQ. 0) RETURN
365
366 CALL timeset(routinen, handle)
367
368 NULLIFY (fm_struct_tmp)
369
370 CALL cp_fm_get_info(matrix=vmatrix, nrow_global=n, ncol_global=ncol_global)
371 IF (ncol .GT. ncol_global) cpabort("Wrong ncol value")
372
373 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=ncol, ncol_global=ncol, &
374 para_env=vmatrix%matrix_struct%para_env, &
375 context=vmatrix%matrix_struct%context)
376 CALL cp_fm_create(overlap_vv, fm_struct_tmp, "overlap_vv")
377 CALL cp_fm_struct_release(fm_struct_tmp)
378
379 CALL parallel_gemm('T', 'N', ncol, ncol, n, rone, vmatrix, vmatrix, rzero, overlap_vv)
380 CALL cp_fm_cholesky_decompose(overlap_vv)
381 CALL cp_fm_triangular_multiply(overlap_vv, vmatrix, n_cols=ncol, side='R', invert_tr=.true.)
382
383 CALL cp_fm_release(overlap_vv)
384
385 CALL timestop(handle)
386
387 END SUBROUTINE make_basis_simple
388
389! **************************************************************************************************
390!> \brief computes ritz values of a set of orbitals given a ks_matrix
391!> rotates the orbitals into eigenstates depending on do_rotation
392!> writes the evals to the screen depending on ionode/scr
393!> \param orbitals S-orthonormal orbitals
394!> \param ks_matrix Kohn-Sham matrix
395!> \param evals_arg optional, filled with the evals
396!> \param ionode , scr if present write to unit scr where ionode
397!> \param scr ...
398!> \param do_rotation optional rotate orbitals (default=.TRUE.)
399!> note that rotating the orbitals is slower
400!> \param co_rotate an optional set of orbitals rotated by the same rotation matrix
401!> \param co_rotate_dbcsr ...
402!> \par History
403!> 08.2004 documented and added do_rotation [Joost VandeVondele]
404!> 09.2008 only compute eigenvalues if rotation is not needed
405! **************************************************************************************************
406 SUBROUTINE subspace_eigenvalues_ks_fm(orbitals, ks_matrix, evals_arg, ionode, scr, &
407 do_rotation, co_rotate, co_rotate_dbcsr)
408
409 TYPE(cp_fm_type), INTENT(IN) :: orbitals
410 TYPE(dbcsr_type), POINTER :: ks_matrix
411 REAL(KIND=dp), DIMENSION(:), OPTIONAL :: evals_arg
412 LOGICAL, INTENT(IN), OPTIONAL :: ionode
413 INTEGER, INTENT(IN), OPTIONAL :: scr
414 LOGICAL, INTENT(IN), OPTIONAL :: do_rotation
415 TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: co_rotate
416 TYPE(dbcsr_type), OPTIONAL, POINTER :: co_rotate_dbcsr
417
418 CHARACTER(len=*), PARAMETER :: routineN = 'subspace_eigenvalues_ks_fm'
419
420 INTEGER :: handle, i, j, n, ncol_global, nrow_global
421 LOGICAL :: compute_evecs, do_rotation_local
422 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: evals
423 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
424 TYPE(cp_fm_type) :: e_vectors, h_block, weighted_vectors, &
425 weighted_vectors2
426
427 CALL timeset(routinen, handle)
428
429 do_rotation_local = .true.
430 IF (PRESENT(do_rotation)) do_rotation_local = do_rotation
431
432 NULLIFY (fm_struct_tmp)
433 CALL cp_fm_get_info(matrix=orbitals, &
434 ncol_global=ncol_global, &
435 nrow_global=nrow_global)
436
437 IF (do_rotation_local) THEN
438 compute_evecs = .true.
439 ELSE
440 ! this would be the logical choice if syevx computing only evals were faster than syevd computing evecs and evals.
441 compute_evecs = .false.
442 ! this is not the case, so lets compute evecs always
443 compute_evecs = .true.
444 END IF
445
446 IF (ncol_global .GT. 0) THEN
447
448 ALLOCATE (evals(ncol_global))
449
450 CALL cp_fm_create(weighted_vectors, orbitals%matrix_struct, "weighted_vectors")
451 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=ncol_global, ncol_global=ncol_global, &
452 para_env=orbitals%matrix_struct%para_env, &
453 context=orbitals%matrix_struct%context)
454 CALL cp_fm_create(h_block, fm_struct_tmp, name="h block")
455 IF (compute_evecs) THEN
456 CALL cp_fm_create(e_vectors, fm_struct_tmp, name="e vectors")
457 END IF
458 CALL cp_fm_struct_release(fm_struct_tmp)
459
460 ! h subblock and diag
461 CALL cp_dbcsr_sm_fm_multiply(ks_matrix, orbitals, weighted_vectors, ncol_global)
462
463 CALL parallel_gemm('T', 'N', ncol_global, ncol_global, nrow_global, 1.0_dp, &
464 orbitals, weighted_vectors, 0.0_dp, h_block)
465
466 ! if eigenvectors are required, go for syevd, otherwise only compute eigenvalues
467 IF (compute_evecs) THEN
468 CALL choose_eigv_solver(h_block, e_vectors, evals)
469 ELSE
470 CALL cp_fm_syevx(h_block, eigenvalues=evals)
471 END IF
472
473 ! rotate the orbitals
474 IF (do_rotation_local) THEN
475 CALL parallel_gemm('N', 'N', nrow_global, ncol_global, ncol_global, 1.0_dp, &
476 orbitals, e_vectors, 0.0_dp, weighted_vectors)
477 CALL cp_fm_to_fm(weighted_vectors, orbitals)
478 IF (PRESENT(co_rotate)) THEN
479 CALL parallel_gemm('N', 'N', nrow_global, ncol_global, ncol_global, 1.0_dp, &
480 co_rotate, e_vectors, 0.0_dp, weighted_vectors)
481 CALL cp_fm_to_fm(weighted_vectors, co_rotate)
482 END IF
483 IF (PRESENT(co_rotate_dbcsr)) THEN
484 IF (ASSOCIATED(co_rotate_dbcsr)) THEN
485 CALL cp_fm_create(weighted_vectors2, orbitals%matrix_struct, "weighted_vectors")
486 CALL copy_dbcsr_to_fm(co_rotate_dbcsr, weighted_vectors2)
487 CALL parallel_gemm('N', 'N', nrow_global, ncol_global, ncol_global, 1.0_dp, &
488 weighted_vectors2, e_vectors, 0.0_dp, weighted_vectors)
489 CALL copy_fm_to_dbcsr(weighted_vectors, co_rotate_dbcsr)
490 CALL cp_fm_release(weighted_vectors2)
491 END IF
492 END IF
493 END IF
494
495 ! give output
496 IF (PRESENT(evals_arg)) THEN
497 n = min(SIZE(evals_arg), SIZE(evals))
498 evals_arg(1:n) = evals(1:n)
499 END IF
500
501 IF (PRESENT(ionode) .OR. PRESENT(scr)) THEN
502 IF (.NOT. PRESENT(ionode)) cpabort("IONODE?")
503 IF (.NOT. PRESENT(scr)) cpabort("SCR?")
504 IF (ionode) THEN
505 DO i = 1, ncol_global, 4
506 j = min(3, ncol_global - i)
507 SELECT CASE (j)
508 CASE (3)
509 WRITE (scr, '(1X,4F16.8)') evals(i:i + j)
510 CASE (2)
511 WRITE (scr, '(1X,3F16.8)') evals(i:i + j)
512 CASE (1)
513 WRITE (scr, '(1X,2F16.8)') evals(i:i + j)
514 CASE (0)
515 WRITE (scr, '(1X,1F16.8)') evals(i:i + j)
516 END SELECT
517 END DO
518 END IF
519 END IF
520
521 CALL cp_fm_release(weighted_vectors)
522 CALL cp_fm_release(h_block)
523 IF (compute_evecs) THEN
524 CALL cp_fm_release(e_vectors)
525 END IF
526
527 DEALLOCATE (evals)
528
529 END IF
530
531 CALL timestop(handle)
532
533 END SUBROUTINE subspace_eigenvalues_ks_fm
534
535! **************************************************************************************************
536!> \brief ...
537!> \param orbitals ...
538!> \param ks_matrix ...
539!> \param evals_arg ...
540!> \param ionode ...
541!> \param scr ...
542!> \param do_rotation ...
543!> \param co_rotate ...
544!> \param para_env ...
545!> \param blacs_env ...
546! **************************************************************************************************
547 SUBROUTINE subspace_eigenvalues_ks_dbcsr(orbitals, ks_matrix, evals_arg, ionode, scr, &
548 do_rotation, co_rotate, para_env, blacs_env)
549
550 TYPE(dbcsr_type), POINTER :: orbitals, ks_matrix
551 REAL(KIND=dp), DIMENSION(:), OPTIONAL :: evals_arg
552 LOGICAL, INTENT(IN), OPTIONAL :: ionode
553 INTEGER, INTENT(IN), OPTIONAL :: scr
554 LOGICAL, INTENT(IN), OPTIONAL :: do_rotation
555 TYPE(dbcsr_type), OPTIONAL, POINTER :: co_rotate
556 TYPE(mp_para_env_type), POINTER :: para_env
557 TYPE(cp_blacs_env_type), POINTER :: blacs_env
558
559 CHARACTER(len=*), PARAMETER :: routineN = 'subspace_eigenvalues_ks_dbcsr'
560
561 INTEGER :: handle, i, j, ncol_global, nrow_global
562 LOGICAL :: compute_evecs, do_rotation_local
563 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: evals
564 TYPE(dbcsr_type), POINTER :: e_vectors, h_block, weighted_vectors
565
566 CALL timeset(routinen, handle)
567
568 do_rotation_local = .true.
569 IF (PRESENT(do_rotation)) do_rotation_local = do_rotation
570
571 NULLIFY (e_vectors, h_block, weighted_vectors)
572
573 CALL dbcsr_get_info(matrix=orbitals, &
574 nfullcols_total=ncol_global, &
575 nfullrows_total=nrow_global)
576
577 IF (do_rotation_local) THEN
578 compute_evecs = .true.
579 ELSE
580 ! this would be the logical choice if syevx computing only evals were faster than syevd computing evecs and evals.
581 compute_evecs = .false.
582 ! this is not the case, so lets compute evecs always
583 compute_evecs = .true.
584 END IF
585
586 IF (ncol_global .GT. 0) THEN
587
588 ALLOCATE (evals(ncol_global))
589
590 CALL dbcsr_init_p(weighted_vectors)
591 CALL dbcsr_copy(weighted_vectors, orbitals, name="weighted_vectors")
592
593 CALL dbcsr_init_p(h_block)
594 CALL cp_dbcsr_m_by_n_from_template(h_block, template=orbitals, m=ncol_global, n=ncol_global, &
595 sym=dbcsr_type_no_symmetry)
596
597!!!!!!!!!!!!!!XXXXXXXXXXXXXXXXXXX!!!!!!!!!!!!!
598 IF (compute_evecs) THEN
599 CALL dbcsr_init_p(e_vectors)
600 CALL cp_dbcsr_m_by_n_from_template(e_vectors, template=orbitals, m=ncol_global, n=ncol_global, &
601 sym=dbcsr_type_no_symmetry)
602 END IF
603
604 ! h subblock and diag
605 CALL dbcsr_multiply('N', 'N', 1.0_dp, ks_matrix, orbitals, &
606 0.0_dp, weighted_vectors)
607 !CALL cp_dbcsr_sm_fm_multiply(ks_matrix,orbitals,weighted_vectors, ncol_global)
608
609 CALL dbcsr_multiply('T', 'N', 1.0_dp, orbitals, weighted_vectors, 0.0_dp, h_block)
610 !CALL parallel_gemm('T','N',ncol_global,ncol_global,nrow_global,1.0_dp, &
611 ! orbitals,weighted_vectors,0.0_dp,h_block)
612
613 ! if eigenvectors are required, go for syevd, otherwise only compute eigenvalues
614 IF (compute_evecs) THEN
615 CALL cp_dbcsr_syevd(h_block, e_vectors, evals, para_env=para_env, blacs_env=blacs_env)
616 ELSE
617 CALL cp_dbcsr_syevx(h_block, eigenvalues=evals, para_env=para_env, blacs_env=blacs_env)
618 END IF
619
620 ! rotate the orbitals
621 IF (do_rotation_local) THEN
622 CALL dbcsr_multiply('N', 'N', 1.0_dp, orbitals, e_vectors, 0.0_dp, weighted_vectors)
623 !CALL parallel_gemm('N','N',nrow_global,ncol_global,ncol_global,1.0_dp, &
624 ! orbitals,e_vectors,0.0_dp,weighted_vectors)
625 CALL dbcsr_copy(orbitals, weighted_vectors)
626 !CALL cp_fm_to_fm(weighted_vectors,orbitals)
627 IF (PRESENT(co_rotate)) THEN
628 IF (ASSOCIATED(co_rotate)) THEN
629 CALL dbcsr_multiply('N', 'N', 1.0_dp, co_rotate, e_vectors, 0.0_dp, weighted_vectors)
630 !CALL parallel_gemm('N','N',nrow_global,ncol_global,ncol_global,1.0_dp, &
631 ! co_rotate,e_vectors,0.0_dp,weighted_vectors)
632 CALL dbcsr_copy(co_rotate, weighted_vectors)
633 !CALL cp_fm_to_fm(weighted_vectors,co_rotate)
634 END IF
635 END IF
636 END IF
637
638 ! give output
639 IF (PRESENT(evals_arg)) THEN
640 evals_arg(:) = evals(:)
641 END IF
642
643 IF (PRESENT(ionode) .OR. PRESENT(scr)) THEN
644 IF (.NOT. PRESENT(ionode)) cpabort("IONODE?")
645 IF (.NOT. PRESENT(scr)) cpabort("SCR?")
646 IF (ionode) THEN
647 DO i = 1, ncol_global, 4
648 j = min(3, ncol_global - i)
649 SELECT CASE (j)
650 CASE (3)
651 WRITE (scr, '(1X,4F16.8)') evals(i:i + j)
652 CASE (2)
653 WRITE (scr, '(1X,3F16.8)') evals(i:i + j)
654 CASE (1)
655 WRITE (scr, '(1X,2F16.8)') evals(i:i + j)
656 CASE (0)
657 WRITE (scr, '(1X,1F16.8)') evals(i:i + j)
658 END SELECT
659 END DO
660 END IF
661 END IF
662
663 CALL dbcsr_release_p(weighted_vectors)
664 CALL dbcsr_release_p(h_block)
665 IF (compute_evecs) THEN
666 CALL dbcsr_release_p(e_vectors)
667 END IF
668
669 DEALLOCATE (evals)
670
671 END IF
672
673 CALL timestop(handle)
674
675 END SUBROUTINE subspace_eigenvalues_ks_dbcsr
676
677! computes the effective orthonormality of a set of mos given an s-matrix
678! orthonormality is the max deviation from unity of the C^T S C
679! **************************************************************************************************
680!> \brief ...
681!> \param orthonormality ...
682!> \param mo_array ...
683!> \param matrix_s ...
684! **************************************************************************************************
685 SUBROUTINE calculate_orthonormality(orthonormality, mo_array, matrix_s)
686 REAL(kind=dp) :: orthonormality
687 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mo_array
688 TYPE(dbcsr_type), OPTIONAL, POINTER :: matrix_s
689
690 CHARACTER(len=*), PARAMETER :: routinen = 'calculate_orthonormality'
691
692 INTEGER :: handle, i, ispin, j, k, n, ncol_local, &
693 nrow_local, nspin
694 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
695 REAL(kind=dp) :: alpha, max_alpha
696 TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
697 TYPE(cp_fm_type) :: overlap, svec
698
699 NULLIFY (tmp_fm_struct)
700
701 CALL timeset(routinen, handle)
702
703 nspin = SIZE(mo_array)
704 max_alpha = 0.0_dp
705
706 DO ispin = 1, nspin
707 IF (PRESENT(matrix_s)) THEN
708 ! get S*C
709 CALL cp_fm_create(svec, mo_array(ispin)%mo_coeff%matrix_struct)
710 CALL cp_fm_get_info(mo_array(ispin)%mo_coeff, &
711 nrow_global=n, ncol_global=k)
712 CALL cp_dbcsr_sm_fm_multiply(matrix_s, mo_array(ispin)%mo_coeff, &
713 svec, k)
714 ! get C^T (S*C)
715 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=k, ncol_global=k, &
716 para_env=mo_array(ispin)%mo_coeff%matrix_struct%para_env, &
717 context=mo_array(ispin)%mo_coeff%matrix_struct%context)
718 CALL cp_fm_create(overlap, tmp_fm_struct)
719 CALL cp_fm_struct_release(tmp_fm_struct)
720 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, mo_array(ispin)%mo_coeff, &
721 svec, 0.0_dp, overlap)
722 CALL cp_fm_release(svec)
723 ELSE
724 ! orthogonal basis C^T C
725 CALL cp_fm_get_info(mo_array(ispin)%mo_coeff, &
726 nrow_global=n, ncol_global=k)
727 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=k, ncol_global=k, &
728 para_env=mo_array(ispin)%mo_coeff%matrix_struct%para_env, &
729 context=mo_array(ispin)%mo_coeff%matrix_struct%context)
730 CALL cp_fm_create(overlap, tmp_fm_struct)
731 CALL cp_fm_struct_release(tmp_fm_struct)
732 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, mo_array(ispin)%mo_coeff, &
733 mo_array(ispin)%mo_coeff, 0.0_dp, overlap)
734 END IF
735 CALL cp_fm_get_info(overlap, nrow_local=nrow_local, ncol_local=ncol_local, &
736 row_indices=row_indices, col_indices=col_indices)
737 DO i = 1, nrow_local
738 DO j = 1, ncol_local
739 alpha = overlap%local_data(i, j)
740 IF (row_indices(i) .EQ. col_indices(j)) alpha = alpha - 1.0_dp
741 max_alpha = max(max_alpha, abs(alpha))
742 END DO
743 END DO
744 CALL cp_fm_release(overlap)
745 END DO
746 CALL mo_array(1)%mo_coeff%matrix_struct%para_env%max(max_alpha)
747 orthonormality = max_alpha
748 ! write(*,*) "max deviation from orthonormalization ",orthonormality
749
750 CALL timestop(handle)
751
752 END SUBROUTINE calculate_orthonormality
753
754! computes the minimum/maximum magnitudes of C^T C. This could be useful
755! to detect problems in the case of nearly singular overlap matrices.
756! in this case, we expect the ratio of min/max to be large
757! this routine is only similar to mo_orthonormality if S==1
758! **************************************************************************************************
759!> \brief ...
760!> \param mo_array ...
761!> \param mo_mag_min ...
762!> \param mo_mag_max ...
763! **************************************************************************************************
764 SUBROUTINE calculate_magnitude(mo_array, mo_mag_min, mo_mag_max)
765 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mo_array
766 REAL(kind=dp) :: mo_mag_min, mo_mag_max
767
768 CHARACTER(len=*), PARAMETER :: routinen = 'calculate_magnitude'
769
770 INTEGER :: handle, ispin, k, n, nspin
771 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: evals
772 TYPE(cp_fm_struct_type), POINTER :: tmp_fm_struct
773 TYPE(cp_fm_type) :: evecs, overlap
774
775 NULLIFY (tmp_fm_struct)
776
777 CALL timeset(routinen, handle)
778
779 nspin = SIZE(mo_array)
780 mo_mag_min = huge(0.0_dp)
781 mo_mag_max = -huge(0.0_dp)
782 DO ispin = 1, nspin
783 CALL cp_fm_get_info(mo_array(ispin)%mo_coeff, &
784 nrow_global=n, ncol_global=k)
785 ALLOCATE (evals(k))
786 CALL cp_fm_struct_create(tmp_fm_struct, nrow_global=k, ncol_global=k, &
787 para_env=mo_array(ispin)%mo_coeff%matrix_struct%para_env, &
788 context=mo_array(ispin)%mo_coeff%matrix_struct%context)
789 CALL cp_fm_create(overlap, tmp_fm_struct)
790 CALL cp_fm_create(evecs, tmp_fm_struct)
791 CALL cp_fm_struct_release(tmp_fm_struct)
792 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, mo_array(ispin)%mo_coeff, &
793 mo_array(ispin)%mo_coeff, 0.0_dp, overlap)
794 CALL choose_eigv_solver(overlap, evecs, evals)
795 mo_mag_min = min(minval(evals), mo_mag_min)
796 mo_mag_max = max(maxval(evals), mo_mag_max)
797 CALL cp_fm_release(overlap)
798 CALL cp_fm_release(evecs)
799 DEALLOCATE (evals)
800 END DO
801 CALL timestop(handle)
802
803 END SUBROUTINE calculate_magnitude
804
805! **************************************************************************************************
806!> \brief Calculate KS eigenvalues starting from OF MOS
807!> \param mos ...
808!> \param nspins ...
809!> \param ks_rmpv ...
810!> \param scf_control ...
811!> \param mo_derivs ...
812!> \param admm_env ...
813!> \param hairy_probes ...
814!> \param probe ...
815!> \par History
816!> 02.2013 moved from qs_scf_post_gpw
817!>
818! **************************************************************************************************
819 SUBROUTINE make_mo_eig(mos, nspins, ks_rmpv, scf_control, mo_derivs, admm_env, &
820 hairy_probes, probe)
821
822 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
823 INTEGER, INTENT(IN) :: nspins
824 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: ks_rmpv
825 TYPE(scf_control_type), POINTER :: scf_control
826 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: mo_derivs
827 TYPE(admm_type), OPTIONAL, POINTER :: admm_env
828 LOGICAL, OPTIONAL :: hairy_probes
829 TYPE(hairy_probes_type), DIMENSION(:), OPTIONAL, &
830 POINTER :: probe
831
832 CHARACTER(len=*), PARAMETER :: routinen = 'make_mo_eig'
833
834 INTEGER :: handle, homo, ispin, nmo, output_unit
835 REAL(kind=dp), DIMENSION(:), POINTER :: mo_eigenvalues
836 TYPE(cp_fm_type), POINTER :: mo_coeff
837 TYPE(dbcsr_type), POINTER :: mo_coeff_deriv
838
839 CALL timeset(routinen, handle)
840
841 NULLIFY (mo_coeff_deriv, mo_coeff, mo_eigenvalues)
842 output_unit = cp_logger_get_default_io_unit()
843
844 DO ispin = 1, nspins
845 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff, &
846 eigenvalues=mo_eigenvalues, homo=homo, nmo=nmo)
847 IF (output_unit > 0) WRITE (output_unit, *) " "
848 IF (output_unit > 0) WRITE (output_unit, *) " Eigenvalues of the occupied subspace spin ", ispin
849 ! IF (homo .NE. nmo) THEN
850 ! IF (output_unit>0) WRITE(output_unit,*)" and ",nmo-homo," added MO eigenvalues"
851 ! END IF
852 IF (output_unit > 0) WRITE (output_unit, *) "---------------------------------------------"
853
854 IF (scf_control%use_ot) THEN
855 ! Also rotate the OT derivs, since they are needed for force calculations
856 IF (ASSOCIATED(mo_derivs)) THEN
857 mo_coeff_deriv => mo_derivs(ispin)%matrix
858 ELSE
859 mo_coeff_deriv => null()
860 END IF
861
862 ! ** If we do ADMM, we add have to modify the kohn-sham matrix
863 IF (PRESENT(admm_env)) THEN
864 CALL admm_correct_for_eigenvalues(ispin, admm_env, ks_rmpv(ispin)%matrix)
865 END IF
866
867 CALL calculate_subspace_eigenvalues(mo_coeff, ks_rmpv(ispin)%matrix, mo_eigenvalues, &
868 scr=output_unit, ionode=output_unit > 0, do_rotation=.true., &
869 co_rotate_dbcsr=mo_coeff_deriv)
870
871 ! ** If we do ADMM, we restore the original kohn-sham matrix
872 IF (PRESENT(admm_env)) THEN
873 CALL admm_uncorrect_for_eigenvalues(ispin, admm_env, ks_rmpv(ispin)%matrix)
874 END IF
875 ELSE
876 IF (output_unit > 0) WRITE (output_unit, '(4(1X,1F16.8))') mo_eigenvalues(1:homo)
877 END IF
878 IF (.NOT. scf_control%diagonalization%mom) THEN
879 IF (PRESENT(hairy_probes) .EQV. .true.) THEN
880 scf_control%smear%do_smear = .false.
881 CALL set_mo_occupation(mo_set=mos(ispin), &
882 smear=scf_control%smear, &
883 probe=probe)
884 ELSE
885 CALL set_mo_occupation(mo_set=mos(ispin), smear=scf_control%smear)
886 END IF
887 END IF
888 IF (output_unit > 0) WRITE (output_unit, '(T2,A,F12.6)') &
889 "Fermi Energy [eV] :", mos(ispin)%mu*evolt
890 END DO
891
892 CALL timestop(handle)
893
894 END SUBROUTINE make_mo_eig
895
896END MODULE qs_mo_methods
Types and set/get functions for auxiliary density matrix methods.
Definition admm_types.F:15
Contains methods used in the context of density fitting.
Definition admm_utils.F:15
subroutine, public admm_uncorrect_for_eigenvalues(ispin, admm_env, ks_matrix)
...
Definition admm_utils.F:127
subroutine, public admm_correct_for_eigenvalues(ispin, admm_env, ks_matrix)
...
Definition admm_utils.F:53
methods related to the blacs parallel environment
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_release_p(matrix)
...
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
subroutine, public dbcsr_multiply(transa, transb, alpha, matrix_a, matrix_b, beta, matrix_c, first_row, last_row, first_column, last_column, first_k, last_k, retain_sparsity, filter_eps, flop)
...
subroutine, public dbcsr_get_info(matrix, nblkrows_total, nblkcols_total, nfullrows_total, nfullcols_total, nblkrows_local, nblkcols_local, nfullrows_local, nfullcols_local, my_prow, my_pcol, local_rows, local_cols, proc_row_dist, proc_col_dist, row_blk_size, col_blk_size, row_blk_offset, col_blk_offset, distribution, name, matrix_type, group)
...
subroutine, public dbcsr_init_p(matrix)
...
Interface to (sca)lapack for the Cholesky based procedures.
subroutine, public cp_dbcsr_syevx(matrix, eigenvectors, eigenvalues, neig, work_syevx, para_env, blacs_env)
compute eigenvalues and optionally eigenvectors of a real symmetric matrix using scalapack....
subroutine, public cp_dbcsr_syevd(matrix, eigenvectors, eigenvalues, para_env, blacs_env)
...
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.
subroutine, public cp_dbcsr_m_by_n_from_template(matrix, template, m, n, sym)
Utility function to create an arbitrary shaped dbcsr matrix with the same processor grid as the templ...
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
Basic linear algebra operations for full matrices.
subroutine, public cp_fm_syrk(uplo, trans, k, alpha, matrix_a, ia, ja, beta, matrix_c)
performs a rank-k update of a symmetric matrix_c matrix_c = beta * matrix_c + alpha * matrix_a * tran...
subroutine, public cp_fm_triangular_multiply(triangular_matrix, matrix_b, side, transpose_tr, invert_tr, uplo_tr, unit_diag_tr, n_rows, n_cols, alpha)
multiplies in place by a triangular matrix: matrix_b = alpha op(triangular_matrix) matrix_b or (if si...
various cholesky decomposition related routines
subroutine, public cp_fm_cholesky_decompose(matrix, n, info_out)
used to replace a symmetric positive def. matrix M with its cholesky decomposition U: M = U^T * U,...
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition cp_fm_diag.F:17
subroutine, public cp_fm_power(matrix, work, exponent, threshold, n_dependent, verbose, eigvals)
...
subroutine, public choose_eigv_solver(matrix, eigenvectors, eigenvalues, info)
Choose the Eigensolver depending on which library is available ELPA seems to be unstable for small sy...
Definition cp_fm_diag.F:238
subroutine, public cp_fm_syevx(matrix, eigenvectors, eigenvalues, neig, work_syevx)
compute eigenvalues and optionally eigenvectors of a real symmetric matrix using scalapack....
Definition cp_fm_diag.F:687
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
Definition cp_fm_types.F:15
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
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Interface to the message passing library MPI.
basic linear algebra operations for full matrixes
Definition of physical constants:
Definition physcon.F:68
real(kind=dp), parameter, public evolt
Definition physcon.F:183
collects routines that perform operations directly related to MOs
subroutine, public make_basis_simple(vmatrix, ncol)
given a set of vectors, return an orthogonal (C^T C == 1) set spanning the same space (notice,...
subroutine, public make_basis_lowdin(vmatrix, ncol, matrix_s)
return a set of S orthonormal vectors (C^T S C == 1) where a Loedwin transformation is applied to kee...
subroutine, public make_basis_cholesky(vmatrix, ncol, ortho)
return a set of S orthonormal vectors (C^T S C == 1) where the cholesky decomposed form of S is passe...
subroutine, public calculate_magnitude(mo_array, mo_mag_min, mo_mag_max)
...
subroutine, public make_basis_sm(vmatrix, ncol, matrix_s)
returns an S-orthonormal basis v (v^T S v ==1)
subroutine, public make_mo_eig(mos, nspins, ks_rmpv, scf_control, mo_derivs, admm_env, hairy_probes, probe)
Calculate KS eigenvalues starting from OF MOS.
subroutine, public calculate_orthonormality(orthonormality, mo_array, matrix_s)
...
Set occupation of molecular orbitals.
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
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.
parameters that control an scf iteration
stores some data used in wavefunction fitting
Definition admm_types.F:120
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
keeps the information about the structure of a full matrix
represent a full matrix
stores all the informations relevant to an mpi environment