Octopus
test.F90
Go to the documentation of this file.
1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
2!!
3!! This program is free software; you can redistribute it and/or modify
4!! it under the terms of the GNU General Public License as published by
5!! the Free Software Foundation; either version 2, or (at your option)
6!! any later version.
7!!
8!! This program is distributed in the hope that it will be useful,
9!! but WITHOUT ANY WARRANTY; without even the implied warranty of
10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11!! GNU General Public License for more details.
12!!
13!! You should have received a copy of the GNU General Public License
14!! along with this program; if not, write to the Free Software
15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
16!! 02110-1301, USA.
17!!
18
20#include "global.h"
21
22module test_oct_m
23 use accel_oct_m
24 use batch_oct_m
25 use batch_ops_oct_m
26 use boundaries_oct_m
27 use box_oct_m
28 use box_factory_oct_m
29 use calc_mode_par_oct_m
30 use cgal_polyhedra_oct_m
31 use chebyshev_filter_bounds_oct_m
32 use clock_oct_m
33 use current_oct_m
34 use debug_oct_m
35 use density_oct_m
36 use derivatives_oct_m
37 use derivatives_test_oct_m
38 use eigen_chebyshev_oct_m
39 use electron_space_oct_m
40 use electrons_oct_m
41 use exponential_oct_m
42 use global_oct_m
43 use grid_oct_m
44 use hamiltonian_elec_oct_m
45 use hamiltonian_elec_base_oct_m
46 use hardware_oct_m
47 use helmholtz_decomposition_m
48 use helmholtz_decomposition_test_m
49 use iihash_oct_m
50 use ion_interaction_oct_m
51 use iso_c_binding
52 use io_csv_oct_m
53 use io_oct_m
54 use io_function_oct_m
55 use, intrinsic :: iso_fortran_env
56 use iteration_counter_oct_m
57 use lalg_basic_oct_m
58 use lalg_adv_oct_m
59 use lalg_adv_tests_oct_m
60 use lapl_operator_oct_m
61 use mesh_oct_m
62 use mesh_batch_oct_m
63 use mesh_function_oct_m
64 use mesh_interpolation_oct_m
65 use messages_oct_m
66 use minimizer_tests_oct_m
67 use maxwell_oct_m
68 use mix_tests_oct_m
69 use mpi_oct_m
70 use mpi_test_oct_m
71 use multicomm_oct_m
72 use multigrid_oct_m
73 use namespace_oct_m
74 use orbitalbasis_oct_m
75 use orbitalset_oct_m
76 use parser_oct_m
77 use poisson_test_oct_m
78 use preconditioners_oct_m
79 use profiling_oct_m
80 use projector_oct_m
81 use regridding_oct_m
82 use sihash_oct_m
83 use solvers_oct_m
84 use space_oct_m
85 use sphash_oct_m
86 use states_abst_oct_m
87 use states_elec_oct_m
88 use states_elec_calc_oct_m
89 use states_elec_dim_oct_m
90 use states_mxll_oct_m
91 use submesh_oct_m
92 use subspace_oct_m
93 use test_parameters_oct_m
94 use types_oct_m
95 use v_ks_oct_m
96 use wfs_elec_oct_m
97 use unit_system_oct_m
98 use io_function_oct_m
99
100 ! TODO(Alex) Issue 951. Migrate the test implementations to separate modules in $OCTOPUS_ROOT/test_integrations/
101 ! Test implementation modules
102 use test_kmeans_clustering_oct_m
103 use test_isdf_oct_m
104 use test_poisson_fft_batch_oct_m
105
106 implicit none
107
108 !Auxiliary quantities needed by the linear solver
109 real(real64), target :: shift_aux
110 type(preconditioner_t) :: prec_aux
111
112 public :: test_run
113
114contains
115
119 subroutine test_run(namespace, grp)
120 type(namespace_t), intent(in) :: namespace
121 type(mpi_grp_t), intent(in) :: grp
122
123 type(test_parameters_t) :: param
124 integer :: test_mode
125
126 push_sub(test_run)
127
128 call messages_obsolete_variable(namespace, 'WhichTest', 'TestMode')
129
130 call param%init_from_file(namespace)
131
132 !%Variable TestMode
133 !%Type integer
134 !%Default hartree
135 !%Section Calculation Modes::Test
136 !%Description
137 !% Decides what kind of test should be performed.
138 !%Option hartree 1
139 !% Tests the Poisson solvers used to calculate the Hartree potential.
140 !%Option derivatives 2
141 !% Tests and benchmarks the implementation of the finite-difference operators, used to calculate derivatives.
142 !%Option orthogonalization 3
143 !% Tests the implementation of the orthogonalization routines.
144 !%Option interpolation 4
145 !% Test the interpolation routines.
146 !%Option ion_interaction 5
147 !% Tests the ion-ion interaction routines.
148 !%Option projector 6
149 !% Tests the code that applies the nonlocal part of the pseudopotentials
150 !% in case of spin-orbit coupling
151 !%Option dft_u 7
152 !% Tests the DFT+U part of the code for projections on the basis.
153 !%Option hamiltonian_apply 8
154 !% Tests the application of the Hamiltonian, or a part of it
155 !%Option density_calc 9
156 !% Calculation of the density.
157 !%Option exp_apply 10
158 !% Tests the exponential of the Hamiltonian
159 !%Option boundaries 11
160 !% Tests the boundaries conditions
161 !%Option subspace_diag 12
162 !% Tests the subspace diagonalization
163 !%Option batch_ops 13
164 !% Tests the batch operations
165 !%Option clock 18
166 !% Tests for clock
167 !%Option linear_solver 19
168 !% Tests the linear solvers
169 !%Option cgal 20
170 !% Tests for cgal interface
171 !%Option dense_eigensolver 21
172 !% Tests for dense eigensolvers (especially parallel ones)
173 !%Option grid_interpolation 22
174 !% Tests for grid interpolation and multigrid methods.
175 !%Option iihash 23
176 !% Tests for the integer-integer hash table.
177 !%Option sihash 24
178 !% Tests for the string-integer hash table.
179 !%Option sphash 25
180 !% Tests for the string-polymorphic hash table.
181 !%Option mpiwrappers 26
182 !% Tests for the MPI wrappers with large integer displacements.
183 !%Option regridding 27
184 !% Tests the regridding between two different grids.
185 !%Option helmholtz_decomposition 28
186 !% Test for the Helmholtz decomposition subroutines
187 !%Option vecpot_analytical 29
188 !% Tests analytically the vector potential from B field.
189 !%Option current_density 30
190 !% Tests the different contributions to the total electronic current density
191 !%Option mixing_tests 31
192 !% Unit test for the mixing
193 !%Option optimizers 32
194 !% Tests for the optimizers, using standard test functions
195 !%Option weighted_kmeans 33
196 !% Tests weighted kmeans algorithm, for assignment of centroids at clustered
197 !% around regions of high weight.
198 !%Option csv_input 34
199 !% Tests the reading of csv files in various dimensions
200 !%Option composition_chebyshev 35
201 !% Tests the chebyshev filtering and the related composition rule
202 !%Option lalg_adv 36
203 !% Tests for the advanced linear algebra routines
204 !%Option isdf_serial 37
205 !% Tests construction of interoperable separable density fit vectors using the
206 !% serial reference implementation.
207 !%Option isdf 38
208 !% Tests construction of interoperable separable density fit vectors.
209 !%Option centroids_t 39
210 !% Test centroids type methods
211 !%Option mesh_generation
212 !% Compares mesh operations for two electronic systems under sibling namespaces.
213 !%Option poisson_fft_batch
214 !% Tests batched FFT Poisson solves for periodic pair-product densities.
215 !%End
216
217 call parse_variable(namespace, 'TestMode', option__testmode__hartree, test_mode)
218
219 call messages_print_with_emphasis(msg="Test mode", namespace=namespace)
220 call messages_print_var_option("TestMode", test_mode, namespace=namespace)
221 call messages_print_var_option("TestType", param%type, namespace=namespace)
222 call messages_print_var_value("TestRepetitions", param%repetitions, namespace=namespace)
223 call messages_print_var_value("TestMinBlockSize", param%min_blocksize, namespace=namespace)
224 call messages_print_var_value("TestMaxBlockSize", param%max_blocksize, namespace=namespace)
225 call messages_print_with_emphasis(namespace=namespace)
226
227 select case (test_mode)
228 case (option__testmode__hartree)
229 call test_hartree(param, namespace, grp)
230 case (option__testmode__derivatives)
231 call test_derivatives(param, namespace, grp)
232 case (option__testmode__orthogonalization)
233 call test_orthogonalization(param, namespace, grp)
234 case (option__testmode__interpolation)
235 call test_interpolation(param, namespace, grp)
236 case (option__testmode__ion_interaction)
237 call test_ion_interaction(namespace, grp)
238 case (option__testmode__projector)
239 call test_projector(param, namespace, grp)
240 case (option__testmode__dft_u)
241 call test_dft_u(param, namespace, grp)
242 case (option__testmode__hamiltonian_apply)
243 call test_hamiltonian(param, namespace, grp)
244 case (option__testmode__density_calc)
245 call test_density_calc(param, namespace, grp)
246 case (option__testmode__exp_apply)
247 call test_exponential(param, namespace, grp)
248 case (option__testmode__boundaries)
249 call test_boundaries(param, namespace, grp)
250 case (option__testmode__subspace_diag)
251 call test_subspace_diagonalization(param, namespace, grp)
252 case (option__testmode__batch_ops)
253 call test_batch_ops(param, namespace, grp)
254 case (option__testmode__clock)
255 call test_clock()
256 case (option__testmode__linear_solver)
257 call test_linear_solver(namespace, grp)
258 case (option__testmode__cgal)
259 call test_cgal()
260 case (option__testmode__dense_eigensolver)
261 call test_dense_eigensolver()
262 case (option__testmode__grid_interpolation)
263 call test_grid_interpolation(grp)
264 case (option__testmode__iihash)
265 call test_iihash()
266 case (option__testmode__sihash)
267 call test_sihash()
268 case (option__testmode__sphash)
269 call test_sphash(namespace)
270 case (option__testmode__mpiwrappers)
271 call test_mpiwrappers()
272 case (option__testmode__regridding)
273 call test_regridding(namespace, grp)
274 case (option__testmode__helmholtz_decomposition)
275 call test_helmholtz_decomposition(namespace, grp)
276 case (option__testmode__vecpot_analytical)
277 call test_vecpot_analytical(namespace, grp)
278 case (option__testmode__current_density)
279 call test_current_density(namespace, grp)
280 case (option__testmode__mixing_tests)
281 call mix_tests_run()
282 case (option__testmode__optimizers)
283 call test_optimizers(namespace)
284 case (option__testmode__weighted_kmeans)
285 call test_weighted_kmeans(namespace)
286 case (option__testmode__csv_input)
287 call test_csv_input(namespace)
288 case (option__testmode__composition_chebyshev)
289 call test_composition_chebyshev(namespace, grp)
290 case (option__testmode__lalg_adv)
291 call test_exponential_matrix(namespace)
292 case (option__testmode__isdf_serial)
293 call test_isdf(namespace, serial=.true.)
294 case (option__testmode__isdf)
295 call test_isdf(namespace)
296 case (option__testmode__mesh_generation)
297 call test_mesh_generation(namespace, grp)
298 case (option__testmode__poisson_fft_batch)
299 call test_poisson_fft_batch(namespace)
300 end select
301
302 pop_sub(test_run)
303 end subroutine test_run
304
305 ! ---------------------------------------------------------
306 subroutine test_hartree(param, namespace, grp)
307 type(test_parameters_t), intent(in) :: param
308 type(namespace_t), intent(in) :: namespace
309 type(mpi_grp_t), intent(in) :: grp
310
311 type(electrons_t), pointer :: sys
312
313 push_sub(test_hartree)
314
315 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
316
317 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
318 call poisson_test(sys%hm%psolver, sys%space, sys%gr, sys%ions%latt, namespace, param%repetitions)
319 safe_deallocate_p(sys)
320
321 pop_sub(test_hartree)
322 end subroutine test_hartree
323
324 ! ---------------------------------------------------------
325 subroutine test_helmholtz_decomposition(namespace, grp)
326 type(namespace_t), intent(in) :: namespace
327 type(mpi_grp_t), intent(in) :: grp
328
329 type(electrons_t), pointer :: sys
330 type(helmholtz_decomposition_t) :: helmholtz
331
332 push_sub(test_helmholtz_decomposition)
333
334 ! First of all we have to initialize the grid and the poisson solver
335 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
336
337 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
338
339 call helmholtz%init(namespace, sys%gr, sys%mc, sys%space)
340
341 ! Then we have to initialize the exact fields
342 write(message(1),'(a)') "Helmholtz decomposition: beginning Hertzian dipole test"
343 call messages_info(1, namespace=namespace)
344 call hertzian_dipole_test(helmholtz, sys%gr, sys%namespace, sys%space)
345
346 write(message(1),'(a)') "Helmholtz decomposition: beginning Gaussian test"
347 call messages_info(1, namespace=namespace)
348 call gaussian_test(helmholtz, sys%gr, sys%namespace, sys%space)
349
350 safe_deallocate_p(sys)
351
352 pop_sub(test_helmholtz_decomposition)
353 end subroutine test_helmholtz_decomposition
354
355 ! ---------------------------------------------------------
356 subroutine test_linear_solver(namespace, grp)
357 use, intrinsic :: iso_c_binding, only: c_ptr
358 type(namespace_t), intent(in) :: namespace
359 type(mpi_grp_t), intent(in) :: grp
360
361 type(electrons_t), pointer :: sys
362 real(real64), allocatable :: rho(:), x(:), center(:)
363 real(real64) :: rr, alpha, beta, res
364 integer :: ip, iter
365 integer, target :: prefactor = -1
366
367 real(real64), parameter :: threshold = 1.e-7_real64
368
369 push_sub(test_linear_solver)
370
371 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
372
373 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
374
375 ! We need to set up some auxiliary quantities called by the linear solver
376 call mesh_init_mesh_aux(sys%gr)
377 ! Shift of the linear equation
378 shift_aux = 0.25_real64
379 ! Preconditioner used for the QMR algorithm
380 call preconditioner_init(prec_aux, namespace, sys%gr, sys%mc, sys%space)
381
382 ! Here we put a Gaussian as the right-hand side of the linear solver
383 ! Values are taken from the poisson_test routine
384 alpha = m_four * sys%gr%spacing(1)
385 beta = m_one / (alpha**sys%space%dim * sqrt(m_pi)**sys%space%dim)
386 ! The Gaussian is centered around the origin
387 safe_allocate(center(1:sys%space%dim))
388 center = m_zero
389
390 safe_allocate(rho(1:sys%gr%np))
391 rho = m_zero
392 do ip = 1, sys%gr%np
393 call mesh_r(sys%gr, ip, rr, origin = center(:))
394 rho(ip) = beta*exp(-(rr/alpha)**2)
395 end do
396
397 safe_allocate(x(1:sys%gr%np))
398
399 !Test the CG linear solver
400 x = m_zero
401 iter = 1000
402
403 call dconjugate_gradients(sys%gr%np, x, rho, dshifted_laplacian_op, dmf_dotp_aux, iter, res, &
404 threshold, userdata=[c_loc(sys%gr%der),c_loc(shift_aux),c_loc(prefactor)])
405 write(message(1),'(a,i6,a)') "Info: CG converged with ", iter, " iterations."
406 write(message(2),'(a,e14.6)') "Info: The residue is ", res
407 write(message(3),'(a,e14.6)') "Info: Norm solution CG ", dmf_nrm2(sys%gr, x)
408 call messages_info(3, namespace=namespace)
409
410 call preconditioner_end(prec_aux)
411 safe_deallocate_a(x)
412 safe_deallocate_a(rho)
413 safe_deallocate_p(sys)
414
415 pop_sub(test_linear_solver)
416 end subroutine test_linear_solver
417
418 ! ---------------------------------------------------------
419 subroutine test_projector(param, namespace, grp)
420 type(test_parameters_t), intent(in) :: param
421 type(namespace_t), intent(in) :: namespace
422 type(mpi_grp_t), intent(in) :: grp
423
424 type(electrons_t), pointer :: sys
425 type(wfs_elec_t) :: epsib
426 integer :: itime
427 complex(real64), allocatable :: psi(:, :)
428
429 push_sub(test_projector)
430
431 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
432
433 call messages_write('Info: Testing the nonlocal part of the pseudopotential with SOC')
434 call messages_new_line()
435 call messages_new_line()
436 call messages_info(namespace=namespace)
437
438 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
439
440 call states_elec_allocate_wfns(sys%st, sys%gr, wfs_type = type_cmplx)
441 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
442
443 ! Initialize external potential
444 call hamiltonian_elec_epot_generate(sys%hm, sys%namespace, sys%space, sys%gr, sys%ions, sys%ext_partners, sys%st)
445
446 call sys%st%group%psib(1, 1)%copy_to(epsib)
447
448 call batch_set_zero(epsib)
449
450 do itime = 1, param%repetitions
451 call zproject_psi_batch(sys%gr, sys%gr%der%boundaries, sys%hm%ep%proj, &
452 sys%hm%ep%natoms, 2, sys%st%group%psib(1, 1), epsib)
453 end do
454
455 safe_allocate(psi(1:sys%gr%np, 1:sys%st%d%dim))
456 do itime = 1, epsib%nst
457 call batch_get_state(epsib, itime, sys%gr%np, psi)
458 write(message(1),'(a,i1,3x, f12.6)') "Norm state ", itime, zmf_nrm2(sys%gr, 2, psi)
459 call messages_info(1, namespace=sys%namespace)
460 end do
461 safe_deallocate_a(psi)
462
463 call epsib%end()
464 call states_elec_deallocate_wfns(sys%st)
465 safe_deallocate_p(sys)
466
467 pop_sub(test_projector)
468 end subroutine test_projector
469
470 ! ---------------------------------------------------------
471 subroutine test_dft_u(param, namespace, grp)
472 type(test_parameters_t), intent(in) :: param
473 type(namespace_t), intent(in) :: namespace
474 type(mpi_grp_t), intent(in) :: grp
475
476 type(electrons_t), pointer :: sys
477 type(wfs_elec_t) :: epsib, epsib2
478 integer :: itime, ist
479 type(orbitalbasis_t) :: basis
480 real(real64), allocatable :: ddot(:,:,:), dweight(:,:)
481 complex(real64), allocatable :: zdot(:,:,:), zweight(:,:)
482 logical :: skipSOrbitals, useAllOrbitals
483
484 push_sub(test_dft_u)
485
486 call calc_mode_par%unset_parallelization(p_strategy_states)
487 call calc_mode_par%unset_parallelization(p_strategy_kpoints)
488 call calc_mode_par%set_parallelization(p_strategy_domains, default = .true.)
489
490 call messages_write('Info: Testing some DFT+U routines')
491 call messages_new_line()
492 call messages_new_line()
493 call messages_info(namespace=namespace)
494
495 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
496
497 call states_elec_allocate_wfns(sys%st, sys%gr)
498 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
499 if (sys%st%pack_states) call sys%st%pack()
500
501 call sys%st%group%psib(1, 1)%copy_to(epsib2, copy_data = .true.)
502
503 ! We set the phase of the batch if needed
504 if (.not. sys%hm%phase%is_allocated()) then
505 call sys%st%group%psib(1, 1)%copy_to(epsib, copy_data = .true.)
506 else
507 call sys%st%group%psib(1, 1)%copy_to(epsib)
508 call sys%hm%phase%apply_to(sys%gr, sys%gr%np, &
509 .false., epsib, src=sys%st%group%psib(1, 1))
510 epsib2%has_phase = .true.
511 end if
512
513 ! Initialize the orbital basis
514 call orbitalbasis_init(basis, sys%namespace, sys%space%periodic_dim)
515
516 ! Parsed in src/hamiltonian/lda_u.F90
517 call parse_variable(namespace, 'UseAllAtomicOrbitals', .false., useallorbitals)
518 call parse_variable(namespace, 'SkipSOrbitals', .true.,skipsorbitals)
519
520 if (states_are_real(sys%st)) then
521 call dorbitalbasis_build(basis, sys%namespace, sys%ions, sys%gr, sys%st%d%kpt, sys%st%d%dim, &
522 skipsorbitals, useallorbitals, verbose=.false.)
523 safe_allocate(dweight(1:basis%orbsets(1)%norbs, 1:epsib%nst_linear))
524 safe_allocate(ddot(1:sys%st%d%dim, 1:basis%orbsets(1)%norbs, 1:epsib%nst))
525 else
526 call zorbitalbasis_build(basis, sys%namespace, sys%ions, sys%gr, sys%st%d%kpt, sys%st%d%dim, &
527 skipsorbitals, useallorbitals, verbose=.false.)
528 call orbitalset_update_phase(basis%orbsets(1), sys%space%dim, sys%st%d%kpt, sys%kpoints, (sys%st%d%ispin==spin_polarized))
529 safe_allocate(zweight(1:basis%orbsets(1)%norbs, 1:epsib%nst_linear))
530 safe_allocate(zdot(1:sys%st%d%dim, 1:basis%orbsets(1)%norbs, 1:epsib%nst))
531
532 !We set the phase of the orbitals if needed
533 if (sys%hm%phase%is_allocated()) then
534 call orbitalset_update_phase(basis%orbsets(1), sys%space%dim, sys%st%d%kpt, sys%kpoints, &
535 (sys%st%d%ispin==spin_polarized))
536 end if
537 end if
538
539 do itime = 1, param%repetitions
540 call batch_set_zero(epsib2)
541 if (states_are_real(sys%st)) then
542 dweight = m_one
543 ddot = m_zero
544 call dorbitalset_get_coeff_batch(basis%orbsets(1), sys%st%d%dim, sys%st%group%psib(1, 1), ddot)
545 call dorbitalset_add_to_batch(basis%orbsets(1), sys%st%d%dim, epsib2, dweight)
546 else
547 zweight = m_one
548 zdot = m_zero
549 call zorbitalset_get_coeff_batch(basis%orbsets(1), sys%st%d%dim, epsib, zdot)
550 call zorbitalset_add_to_batch(basis%orbsets(1), sys%st%d%dim, epsib2, zweight)
551 end if
552 end do
553
554 if (epsib%is_packed()) then
555 call epsib%do_unpack(force = .true.)
556 end if
557
558 do ist = 1, epsib%nst
559 if (states_are_real(sys%st)) then
560 write(message(1),'(a,i2,3x,e13.6)') "Dotp state ", ist, ddot(1,1,ist)
561 else
562 write(message(1),'(a,i2,2(3x,e13.6))') "Dotp state ", ist, zdot(1,1,ist)
563 end if
564 call messages_info(1)
565 end do
566
567
568 call test_prints_info_batch(sys%st, sys%gr, epsib2)
569
570 safe_deallocate_a(dweight)
571 safe_deallocate_a(zweight)
572 safe_deallocate_a(ddot)
573 safe_deallocate_a(zdot)
574
575 call epsib%end()
576 call epsib2%end()
577 call orbitalbasis_end(basis)
578 call states_elec_deallocate_wfns(sys%st)
579 safe_deallocate_p(sys)
580
581 pop_sub(test_dft_u)
582 end subroutine test_dft_u
583
584 ! ---------------------------------------------------------
585 subroutine test_hamiltonian(param, namespace, grp)
586 type(test_parameters_t), intent(in) :: param
587 type(namespace_t), intent(in) :: namespace
588 type(mpi_grp_t), intent(in) :: grp
589
590 type(electrons_t), pointer :: sys
591 type(hamiltonian_elec_t) :: hm_copy
592 type(wfs_elec_t) :: hpsib, hpsib_copy, hpsib_diff
593 integer :: itime, terms
594 integer :: ist
595 logical :: copy_before_apply
596 real(real64), allocatable :: norm_hpsib(:), norm_diff(:)
597 real(real64) :: max_rel_diff
598 real(real64), parameter :: hamiltonian_copy_rel_tol = 1.0e-12_real64
599
600 push_sub(test_hamiltonian)
601
602 !%Variable TestHamiltonianApply
603 !%Type integer
604 !%Default term_all
605 !%Section Calculation Modes::Test
606 !%Description
607 !% Decides which part of the Hamiltonian is applied.
608 !%Option term_all 0
609 !% Apply the full Hamiltonian.
610 !%Option term_kinetic 1
611 !% Apply only the kinetic operator
612 !%Option term_local_potential 2
613 !% Apply only the local potential.
614 !%Option term_non_local_potential 4
615 !% Apply only the non_local potential.
616 !%End
617 call parse_variable(namespace, 'TestHamiltonianApply', option__testhamiltonianapply__term_all, terms)
618 if (terms == 0) terms = huge(1)
619
620 !%Variable TestHamiltonianCopyBeforeApply
621 !%Type logical
622 !%Default no
623 !%Section Calculation Modes::Test
624 !%Description
625 !% If set to yes, the Hamiltonian is copied with hamiltonian_elec_copy
626 !% and the copied Hamiltonian is used for application.
627 !%End
628 call parse_variable(namespace, 'TestHamiltonianCopyBeforeApply', .false., copy_before_apply)
629 call messages_print_var_value('TestHamiltonianCopyBeforeApply', copy_before_apply, namespace=namespace)
630
631
632 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
633
634 call messages_write('Info: Testing the application of the Hamiltonian')
635 call messages_new_line()
636 call messages_new_line()
637 call messages_info(namespace=namespace)
638
639 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
640
641 call states_elec_allocate_wfns(sys%st, sys%gr)
642 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
643
644 ! Initialize external potential
645 if (sys%st%pack_states .and. sys%hm%apply_packed()) call sys%st%pack()
646 call hamiltonian_elec_epot_generate(sys%hm, sys%namespace, sys%space, sys%gr, sys%ions, sys%ext_partners, sys%st)
647 call density_calc(sys%st, sys%gr, sys%st%rho)
648 call v_ks_calc(sys%ks, sys%namespace, sys%space, sys%hm, sys%st, sys%ions, sys%ext_partners)
649
650 if (copy_before_apply) then
651 call hamiltonian_elec_copy(hm_copy, sys%hm)
652 end if
653
654 call boundaries_set(sys%gr%der%boundaries, sys%gr, sys%st%group%psib(1, 1))
655
656 call sys%st%group%psib(1, 1)%copy_to(hpsib)
657 if (copy_before_apply) call sys%st%group%psib(1, 1)%copy_to(hpsib_copy)
658
659 if (sys%hm%apply_packed()) then
660 call sys%st%group%psib(1, 1)%do_pack()
661 call hpsib%do_pack(copy = .false.)
662 if (copy_before_apply) call hpsib_copy%do_pack(copy = .false.)
663 end if
664
665 do itime = 1, param%repetitions
666 if (states_are_real(sys%st)) then
667 call dhamiltonian_elec_apply_batch(sys%hm, sys%namespace, sys%gr, sys%st%group%psib(1, 1), hpsib, &
668 terms = terms, set_bc = .false.)
669 if (copy_before_apply) then
670 call dhamiltonian_elec_apply_batch(hm_copy, sys%namespace, sys%gr, sys%st%group%psib(1, 1), hpsib_copy, &
671 terms = terms, set_bc = .false.)
672 end if
673 else
674 call zhamiltonian_elec_apply_batch(sys%hm, sys%namespace, sys%gr, sys%st%group%psib(1, 1), hpsib, &
675 terms = terms, set_bc = .false.)
676 if (copy_before_apply) then
677 call zhamiltonian_elec_apply_batch(hm_copy, sys%namespace, sys%gr, sys%st%group%psib(1, 1), hpsib_copy, &
678 terms = terms, set_bc = .false.)
679 end if
680 end if
681 end do
682
683 if (hpsib%is_packed()) then
684 call hpsib%do_unpack(force = .true.)
685 end if
686 if (copy_before_apply .and. hpsib_copy%is_packed()) then
687 call hpsib_copy%do_unpack(force = .true.)
688 end if
689
690 if (copy_before_apply) then
691 call hpsib_copy%copy_to(hpsib_diff, copy_data = .true.)
692 call batch_axpy(sys%gr%np, -1.0_real64, hpsib, hpsib_diff)
693
694 safe_allocate(norm_hpsib(1:hpsib%nst))
695 safe_allocate(norm_diff(1:hpsib_diff%nst))
696 call mesh_batch_nrm2(sys%gr, hpsib, norm_hpsib)
697 call mesh_batch_nrm2(sys%gr, hpsib_diff, norm_diff)
698
699 max_rel_diff = maxval(norm_diff / max(norm_hpsib, tiny(1.0_real64)))
700 write(message(1),'(a,e23.16)') 'Hamiltonian copy max relative difference ', max_rel_diff
701 call messages_info(1, namespace=namespace)
702
703 do ist = 1, hpsib%nst
704 if (norm_diff(ist) > hamiltonian_copy_rel_tol*max(norm_hpsib(ist), tiny(1.0_real64))) then
705 write(message(1),'(a,i6)') 'Hamiltonian copy mismatch for state ', ist
706 call messages_fatal(1, namespace=namespace)
707 end if
708 end do
709
710 safe_deallocate_a(norm_hpsib)
711 safe_deallocate_a(norm_diff)
712 call hpsib_diff%end(copy = .false.)
713 end if
714
715 call test_prints_info_batch(sys%st, sys%gr, hpsib)
716
717 if (copy_before_apply) then
718 call hamiltonian_elec_end(hm_copy)
719 call hpsib_copy%end(copy = .false.)
720 ! check that original hm still works after teardown of copy
721 if (sys%hm%apply_packed()) then
722 call hpsib%do_pack(copy = .false.)
723 end if
724 if (states_are_real(sys%st)) then
725 call dhamiltonian_elec_apply_batch(sys%hm, sys%namespace, sys%gr, sys%st%group%psib(1, 1), hpsib, &
726 terms = terms, set_bc = .false.)
727 else
728 call zhamiltonian_elec_apply_batch(sys%hm, sys%namespace, sys%gr, sys%st%group%psib(1, 1), hpsib, &
729 terms = terms, set_bc = .false.)
730 end if
731 end if
732 call hpsib%end(copy = .false.)
733 call states_elec_deallocate_wfns(sys%st)
734 safe_deallocate_p(sys)
735
736 pop_sub(test_hamiltonian)
737 end subroutine test_hamiltonian
738
739
740 ! ---------------------------------------------------------
741 subroutine test_density_calc(param, namespace, grp)
742 type(test_parameters_t), intent(in) :: param
743 type(namespace_t), intent(in) :: namespace
744 type(mpi_grp_t), intent(in) :: grp
745
746 type(electrons_t), pointer :: sys
747 integer :: itime
748
749 push_sub(test_density_calc)
750
751 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
752
753 call messages_write('Info: Testing density calculation')
754 call messages_new_line()
755 call messages_new_line()
756 call messages_info(namespace=namespace)
757
758 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
759
760 call states_elec_allocate_wfns(sys%st, sys%gr)
761 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
762 if (sys%st%pack_states) call sys%st%pack()
763
764 do itime = 1, param%repetitions
765 call density_calc(sys%st, sys%gr, sys%st%rho)
766 end do
767
768 write(message(1),'(a,3x, f12.6)') "Norm density ", dmf_nrm2(sys%gr, sys%st%rho(:,1))
769 call messages_info(1, namespace=namespace)
770
771 call states_elec_deallocate_wfns(sys%st)
772 safe_deallocate_p(sys)
773
774 pop_sub(test_density_calc)
775 end subroutine test_density_calc
776
777
778 ! ---------------------------------------------------------
779 subroutine test_boundaries(param, namespace, grp)
780 type(test_parameters_t), intent(in) :: param
781 type(namespace_t), intent(in) :: namespace
782 type(mpi_grp_t), intent(in) :: grp
783
784 type(electrons_t), pointer :: sys
785 integer :: itime
786
787 push_sub(test_boundaries)
788
789 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
790
791 call messages_write('Info: Testing boundary conditions')
792 call messages_new_line()
793 call messages_new_line()
794 call messages_info(namespace=namespace)
795
796 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
797
798 call states_elec_allocate_wfns(sys%st, sys%gr)
799 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
800 if (sys%st%pack_states) call sys%st%pack()
801
802 do itime = 1, param%repetitions
803 call boundaries_set(sys%gr%der%boundaries, sys%gr, sys%st%group%psib(1, 1))
804 end do
805
806 call test_prints_info_batch(sys%st, sys%gr, sys%st%group%psib(1, 1))
807
808 call states_elec_deallocate_wfns(sys%st)
809 safe_deallocate_p(sys)
810
811 pop_sub(test_boundaries)
812 end subroutine test_boundaries
813
814 ! ---------------------------------------------------------
818 subroutine test_composition_chebyshev(namespace, grp)
819 type(namespace_t), intent(in) :: namespace
820 type(mpi_grp_t), intent(in) :: grp
821
822 type(electrons_t), pointer :: sys
823 type(chebyshev_filter_bounds_t), pointer :: bounds
824 integer, allocatable :: degree(:)
825 integer :: degree_n
826 type(wfs_elec_t) :: psib
827
828 push_sub(test_composition_chebyshev)
829
830 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
831
832 call messages_write('Info: Testing Chebyshev filtering and its composition rule')
833 call messages_new_line()
834 call messages_new_line()
835 call messages_info(namespace=namespace)
836
837 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
838
839 call states_elec_allocate_wfns(sys%st, sys%gr)
840
841 if (sys%st%pack_states .and. sys%hm%apply_packed()) call sys%st%pack()
842 call hamiltonian_elec_epot_generate(sys%hm, sys%namespace, sys%space, sys%gr, sys%ions, sys%ext_partners, sys%st)
843 call density_calc(sys%st, sys%gr, sys%st%rho)
844 call v_ks_calc(sys%ks, sys%namespace, sys%space, sys%hm, sys%st, sys%ions, sys%ext_partners)
845
846 bounds => chebyshev_filter_bounds_t(-1.0_real64, 1.0_real64)
847 safe_allocate(degree(1:sys%st%group%nblocks))
848
849 !%Variable TestCompositionOrder
850 !%Type integer
851 !%Section Calculation Modes::Test
852 !%Description
853 !% Specifies the value of m for the test of the composition rule
854 !% $T_{2n}(x) = 2*T_n(x)^2-1$.
855 !%
856 !%End
857 call parse_variable(namespace, 'TestCompositionOrder', 2, degree_n)
858
859 !Initialization
860 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
861 call messages_new_line()
862 call messages_write('Info: Result for the original vector')
863 call messages_info(namespace=namespace)
864 call test_prints_info_batch(sys%st, sys%gr, sys%st%group%psib(1, 1))
865
866 ! Direct evaluation of mxn-th order Chebyshev filtering
867 degree = 2*degree_n
868 if (states_are_real(sys%st)) then
869 call dchebyshev_filter(namespace, sys%gr, sys%st, sys%hm, degree, bounds, 1, normalize=.false.)
870 else
871 call zchebyshev_filter(namespace, sys%gr, sys%st, sys%hm, degree, bounds, 1, normalize=.false.)
872 end if
873
874 call messages_new_line()
875 call messages_write('Info: Result after calling 2n-th order filtering')
876 call messages_info(namespace=namespace)
877 call test_prints_info_batch(sys%st, sys%gr, sys%st%group%psib(1, 1))
878
879 ! Composition rule: m-th order Chebyshev filtering followed by n-th order filtering
880 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
881 call sys%st%group%psib(1, 1)%copy_to(psib, copy_data=.true.)
882 degree = degree_n
883 if (states_are_real(sys%st)) then
884 call dchebyshev_filter(namespace, sys%gr, sys%st, sys%hm, degree, bounds, 1, normalize=.false.)
885 call dchebyshev_filter(namespace, sys%gr, sys%st, sys%hm, degree, bounds, 1, normalize=.false.)
886 else
887 call zchebyshev_filter(namespace, sys%gr, sys%st, sys%hm, degree, bounds, 1, normalize=.false.)
888 call zchebyshev_filter(namespace, sys%gr, sys%st, sys%hm, degree, bounds, 1, normalize=.false.)
889 end if
890 call batch_scal(sys%gr%np, m_two, sys%st%group%psib(1, 1))
891 call batch_axpy(sys%gr%np, -m_one, psib, sys%st%group%psib(1, 1))
892
893 call messages_new_line()
894 call messages_write('Info: Result using the composition rule')
895 call messages_info(namespace=namespace)
896 call test_prints_info_batch(sys%st, sys%gr, sys%st%group%psib(1, 1))
897
898
899 safe_deallocate_a(degree)
900 safe_deallocate_p(bounds)
901
902 call psib%end()
903 call states_elec_deallocate_wfns(sys%st)
904 safe_deallocate_p(sys)
905
906 pop_sub(test_composition_chebyshev)
907 end subroutine test_composition_chebyshev
908
909
910
911 ! ---------------------------------------------------------
912 subroutine test_exponential(param, namespace, grp)
913 type(test_parameters_t), intent(in) :: param
914 type(namespace_t), intent(in) :: namespace
915 type(mpi_grp_t), intent(in) :: grp
916
917 type(electrons_t), pointer :: sys
918 type(exponential_t) :: te
919 integer :: itime
920
921 push_sub(test_exponential)
922
923 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
924
925 call messages_write('Info: Testing exponential')
926 call messages_new_line()
927 call messages_new_line()
928 call messages_info(namespace=namespace)
929
930 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
931
932 call states_elec_allocate_wfns(sys%st, sys%gr, wfs_type=type_cmplx)
933 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
934
935 ! Initialize external potential
936 if (sys%st%pack_states .and. sys%hm%apply_packed()) call sys%st%pack()
937 call hamiltonian_elec_epot_generate(sys%hm, sys%namespace, sys%space, sys%gr, sys%ions, sys%ext_partners, sys%st)
938 call density_calc(sys%st, sys%gr, sys%st%rho)
939 call v_ks_calc(sys%ks, sys%namespace, sys%space, sys%hm, sys%st, sys%ions, sys%ext_partners)
940
941 call exponential_init(te, namespace)
942
943 if (sys%hm%apply_packed()) then
944 call sys%st%group%psib(1, 1)%do_pack()
945 end if
946
947 do itime = 1, param%repetitions
948 call te%apply_batch(sys%namespace, sys%gr, sys%hm, sys%st%group%psib(1, 1), 1e-3_real64)
949 end do
950
951 call test_prints_info_batch(sys%st, sys%gr, sys%st%group%psib(1, 1))
952
953 call states_elec_deallocate_wfns(sys%st)
954 safe_deallocate_p(sys)
955
956 pop_sub(test_exponential)
957 end subroutine test_exponential
958
959
960 ! ---------------------------------------------------------
961 subroutine test_subspace_diagonalization(param, namespace, grp)
962 type(test_parameters_t), intent(in) :: param
963 type(namespace_t), intent(in) :: namespace
964 type(mpi_grp_t), intent(in) :: grp
965
966 type(electrons_t), pointer :: sys
967 integer :: itime
968 type(subspace_t) :: sdiag
969 real(real64), allocatable :: diff(:)
970
971 push_sub(test_subspace_diagonalization)
972
973 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
974
975 call messages_write('Info: Testing subspace diagonalization')
976 call messages_new_line()
977 call messages_new_line()
978 call messages_info(namespace=namespace)
979
980 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
981
982 call states_elec_allocate_wfns(sys%st, sys%gr)
983 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
984
985 if (sys%st%pack_states .and. sys%hm%apply_packed()) call sys%st%pack()
986 call hamiltonian_elec_epot_generate(sys%hm, sys%namespace, sys%space, sys%gr, sys%ions, sys%ext_partners, sys%st)
987 call density_calc(sys%st, sys%gr, sys%st%rho)
988 call v_ks_calc(sys%ks, sys%namespace, sys%space, sys%hm, sys%st, sys%ions, sys%ext_partners)
989
990 call sdiag%init(sys%namespace, sys%st)
991
992 safe_allocate(diff(1:sys%st%nst))
993
994 do itime = 1, param%repetitions
995 if (states_are_real(sys%st)) then
996 call dsubspace_diag(sdiag, sys%namespace, sys%gr, sys%st, sys%hm, 1, sys%st%eigenval(:, 1), diff)
997 else
998 call zsubspace_diag(sdiag, sys%namespace, sys%gr, sys%st, sys%hm, 1, sys%st%eigenval(:, 1), diff)
999 end if
1000 end do
1001
1002 safe_deallocate_a(diff)
1003
1004 call test_prints_info_batch(sys%st, sys%gr, sys%st%group%psib(1, 1))
1005
1006 call states_elec_deallocate_wfns(sys%st)
1007 safe_deallocate_p(sys)
1008
1009 pop_sub(test_subspace_diagonalization)
1010 end subroutine test_subspace_diagonalization
1011
1012
1013 ! ---------------------------------------------------------
1014 subroutine test_batch_ops(param, namespace, grp)
1015 type(test_parameters_t), intent(in) :: param
1016 type(namespace_t), intent(in) :: namespace
1017 type(mpi_grp_t), intent(in) :: grp
1018
1019 type(electrons_t), pointer :: sys
1020 integer :: itime, ops, ops_default, ist, jst, nst, ip
1021 type(wfs_elec_t) :: xx, yy, zz
1022 real(real64), allocatable :: tmp(:)
1023 real(real64), allocatable :: ddotv(:)
1024 complex(real64), allocatable :: zdotv(:)
1025 real(real64), allocatable :: ddot(:,:), df(:,:), dweight(:), dpoints(:,:,:)
1026 complex(real64), allocatable :: zdot(:,:), zf(:,:), zweight(:), zpoints(:,:,:)
1027 real(real64), allocatable :: dff_mf(:)
1028 complex(real64), allocatable :: zff_mf(:)
1029 integer :: sp, block_size, size
1030
1031
1032 push_sub(test_batch_ops)
1033
1034 !%Variable TestBatchOps
1035 !%Type flag
1036 !%Default all
1037 !%Section Calculation Modes::Test
1038 !%Description
1039 !% Decides which part of the Hamiltonian is applied.
1040 !%Option ops_axpy bit(1)
1041 !% Tests batch_axpy operation
1042 !%Option ops_scal bit(2)
1043 !% Tests batch_scal operation
1044 !%Option ops_nrm2 bit(3)
1045 !% Tests batch_nrm2 operation
1046 !%Option ops_dotp_matrix bit(4)
1047 !% Tests X(mesh_batch_dotp_matrix)
1048 !%Option ops_dotp_self bit(5)
1049 !% Tests X(mesh_batch_dotp_self)
1050 !%Option ops_dotp_vector bit(6)
1051 !% Tests X(mesh_batch_dotp_vector)
1052 !%Option ops_ax_function_py bit(7)
1053 !% Tests X(batch_ax_function_py)
1054 !%Option ops_get_points bit(8)
1055 !% Tests batch_get_points and batch_set_points
1056 !%End
1057 ops_default = &
1058 option__testbatchops__ops_axpy + &
1059 option__testbatchops__ops_scal + &
1060 option__testbatchops__ops_nrm2 + &
1061 option__testbatchops__ops_dotp_matrix + &
1062 option__testbatchops__ops_dotp_self + &
1063 option__testbatchops__ops_dotp_vector + &
1064 option__testbatchops__ops_ax_function_py + &
1065 option__testbatchops__ops_get_points
1066
1067 call parse_variable(namespace, 'TestBatchOps', ops_default, ops)
1068
1069 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
1070
1071 call messages_write('Info: Testing batch operations')
1072 call messages_new_line()
1073 call messages_new_line()
1074 call messages_info(namespace=namespace)
1075
1076 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
1077
1078 call states_elec_allocate_wfns(sys%st, sys%gr)
1079 call test_batch_set_gaussian(sys%st%group%psib(1, 1), sys%gr)
1080 if (sys%st%pack_states) call sys%st%pack()
1081
1082 if (bitand(ops, option__testbatchops__ops_axpy) /= 0) then
1083 message(1) = 'Info: Testing axpy'
1084 call messages_info(1, namespace=sys%namespace)
1085
1086 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1087 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .true.)
1088
1089 do itime = 1, param%repetitions
1090 call batch_axpy(sys%gr%np, 0.1_real64, xx, yy)
1091 end do
1092 call test_prints_info_batch(sys%st, sys%gr, yy, string = "axpy")
1093
1094 call xx%end()
1095 call yy%end()
1096
1097 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1098 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .true.)
1099 do itime = 1, param%repetitions
1100 call batch_axpby(sys%gr%np, 0.1_real64, xx, 0.2_real64, yy)
1101 end do
1102 call test_prints_info_batch(sys%st, sys%gr, yy, string = "axpby")
1103
1104 call xx%end()
1105 call yy%end()
1106 end if
1107
1108 if (bitand(ops, option__testbatchops__ops_scal) /= 0) then
1109 message(1) = 'Info: Testing scal'
1110 call messages_info(1, namespace=sys%namespace)
1111
1112 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1113
1114 do itime = 1, param%repetitions
1115 call batch_scal(sys%gr%np, 0.1_real64, xx)
1116 end do
1117 call test_prints_info_batch(sys%st, sys%gr, xx, string="scal")
1118 call xx%end()
1119
1120 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1121 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .false.)
1122 do itime = 1, param%repetitions
1123 call batch_scal2v(sys%gr%np, 0.1_real64, xx, yy)
1124 end do
1125 call test_prints_info_batch(sys%st, sys%gr, yy, string="scal2v")
1126
1127 call xx%end()
1128 call yy%end()
1129 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1130 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .true.)
1131 call sys%st%group%psib(1, 1)%copy_to(zz, copy_data = .false.)
1132
1133 do itime = 1, param%repetitions
1134 call batch_pointwise_mul(sys%gr%np, xx, yy, zz)
1135 end do
1136 call test_prints_info_batch(sys%st, sys%gr, zz, string="hadamard")
1137 call xx%end()
1138 call yy%end()
1139 call zz%end()
1140
1141 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1142 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .false.)
1143 if (states_are_real(sys%st)) then
1144 safe_allocate(dff_mf(1:sys%gr%np))
1145 do ip = 1, sys%gr%np
1146 dff_mf(ip) = 0.01_real64 + 0.001_real64*real(sys%gr%x(1, ip), real64)
1147 end do
1148
1149 do itime = 1, param%repetitions
1150 call batch_mul_mf(sys%gr%np, dff_mf, xx, yy)
1151 end do
1152 else
1153 safe_allocate(zff_mf(1:sys%gr%np))
1154 do ip = 1, sys%gr%np
1155 zff_mf(ip) = cmplx(0.01_real64 + 0.001_real64*real(sys%gr%x(1, ip), real64), 0.02_real64, real64)
1156 end do
1157
1158 do itime = 1, param%repetitions
1159 call batch_mul_mf(sys%gr%np, zff_mf, xx, yy)
1160 end do
1161 end if
1162
1163 call test_prints_info_batch(sys%st, sys%gr, yy, string="mul_mf")
1164
1165
1166 safe_deallocate_a(dff_mf)
1167 safe_deallocate_a(zff_mf)
1168
1169 call xx%end()
1170 call yy%end()
1171 end if
1172
1173 if (bitand(ops, option__testbatchops__ops_nrm2) /= 0) then
1174 message(1) = 'Info: Testing nrm2'
1175 call messages_info(1, namespace=sys%namespace)
1176
1177 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1178 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .true.)
1179
1180 safe_allocate(tmp(1:xx%nst))
1181
1182 do itime = 1, param%repetitions
1183 call mesh_batch_nrm2(sys%gr, yy, tmp)
1184 end do
1185 do itime = 1, xx%nst
1186 write(message(1),'(a,i1,3x,e23.16)') "Nrm2 norm state ", itime, tmp(itime)
1187 call messages_info(1, namespace=sys%namespace)
1188 end do
1189
1190 safe_deallocate_a(tmp)
1191
1192 call xx%end()
1193 call yy%end()
1194 end if
1195
1196 if (bitand(ops, option__testbatchops__ops_dotp_matrix) /= 0) then
1197
1198 message(1) = 'Info: Testing dotp_matrix'
1199 call messages_info(1, namespace=sys%namespace)
1200
1201 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1202 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .true.)
1203
1204 nst = sys%st%group%psib(1, 1)%nst
1205
1206 if (states_are_real(sys%st)) then
1207 safe_allocate(ddot(1:nst, 1:nst))
1208 do itime = 1, param%repetitions
1209 call dmesh_batch_dotp_matrix(sys%gr, xx, yy, ddot)
1210 end do
1211
1212 do ist = 1, nst
1213 do jst = 1, nst
1214 write(message(jst), '(a,2i3,3x,e23.16)') 'Dotp_matrix states', ist, jst, ddot(ist,jst)
1215 end do
1216 call messages_info(nst, namespace=sys%namespace)
1217 end do
1218 safe_deallocate_a(ddot)
1219 else
1220 safe_allocate(zdot(1:nst, 1:nst))
1221 do itime = 1, param%repetitions
1222 call zmesh_batch_dotp_matrix(sys%gr, xx, yy, zdot)
1223 end do
1224
1225 do ist = 1, nst
1226 do jst = 1, nst
1227 write(message(jst), '(a,2i3,3x,2(e23.16,1x))') 'Dotp_matrix states', ist, jst, zdot(ist,jst)
1228 end do
1229 call messages_info(nst, namespace=sys%namespace)
1230 end do
1231 safe_deallocate_a(zdot)
1232 end if
1233
1234 call xx%end()
1235 call yy%end()
1236
1237 message(1) = 'Info: Testing dotp_matrix with distinct batches'
1238 call messages_info(1, namespace=sys%namespace)
1239
1240 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1241 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .false.)
1242
1243 nst = sys%st%group%psib(1, 1)%nst
1244
1245 if (states_are_real(sys%st)) then
1246 safe_allocate(dff_mf(1:sys%gr%np))
1247 do ip = 1, sys%gr%np
1248 dff_mf(ip) = 0.5_real64 + 0.1_real64*real(sys%gr%x(1, ip), real64)
1249 end do
1250 call batch_mul_mf(sys%gr%np, dff_mf, xx, yy)
1251 safe_deallocate_a(dff_mf)
1252
1253 safe_allocate(ddot(1:nst, 1:nst))
1254 call dmesh_batch_dotp_matrix(sys%gr, xx, yy, ddot)
1255 do ist = 1, nst
1256 do jst = 1, nst
1257 write(message(jst), '(a,2i3,3x,e23.16)') 'Dotp_matrix_ab states', ist, jst, ddot(ist,jst)
1258 end do
1259 call messages_info(nst, namespace=sys%namespace)
1260 end do
1261 safe_deallocate_a(ddot)
1262 else
1263 safe_allocate(zff_mf(1:sys%gr%np))
1264 do ip = 1, sys%gr%np
1265 zff_mf(ip) = cmplx(0.5_real64 + 0.1_real64*real(sys%gr%x(1, ip), real64), &
1266 0.3_real64*real(sys%gr%x(1, ip), real64), real64)
1267 end do
1268 call batch_mul_mf(sys%gr%np, zff_mf, xx, yy)
1269 safe_deallocate_a(zff_mf)
1270
1271 safe_allocate(zdot(1:nst, 1:nst))
1272 call zmesh_batch_dotp_matrix(sys%gr, xx, yy, zdot)
1273 do ist = 1, nst
1274 do jst = 1, nst
1275 write(message(jst), '(a,2i3,3x,2(e23.16,1x))') 'Dotp_matrix_ab states', ist, jst, zdot(ist,jst)
1276 end do
1277 call messages_info(nst, namespace=sys%namespace)
1278 end do
1279 safe_deallocate_a(zdot)
1280 end if
1281
1282 call xx%end()
1283 call yy%end()
1284 end if
1285
1286 if (bitand(ops, option__testbatchops__ops_dotp_vector) /= 0) then
1287
1288 message(1) = 'Info: Testing dotp_vector'
1289 call messages_info(1, namespace=sys%namespace)
1290
1291 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1292 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .true.)
1293
1294 nst = sys%st%group%psib(1, 1)%nst
1295
1296 if (states_are_real(sys%st)) then
1297 safe_allocate(ddotv(1:nst))
1298 do itime = 1, param%repetitions
1299 call dmesh_batch_dotp_vector(sys%gr, xx, yy, ddotv)
1300 end do
1301
1302 do ist = 1, nst
1303 write(message(ist), '(a,i3,3x,e23.16)') 'Dotp_vector state', ist, ddotv(ist)
1304 end do
1305 call messages_info(nst, namespace=sys%namespace)
1306 safe_deallocate_a(ddotv)
1307 else
1308 safe_allocate(zdotv(1:nst))
1309 do itime = 1, param%repetitions
1310 call zmesh_batch_dotp_vector(sys%gr, xx, yy, zdotv)
1311 end do
1312
1313 do ist = 1, nst
1314 write(message(ist), '(a,i3,3x,2(e23.16,1x))') 'Dotp_vector state', ist, zdotv(ist)
1315 end do
1316 call messages_info(nst, namespace=sys%namespace)
1317 safe_deallocate_a(zdotv)
1318 end if
1319
1320 call xx%end()
1321 call yy%end()
1322 end if
1323
1324 if (bitand(ops, option__testbatchops__ops_dotp_self) /= 0) then
1325
1326 message(1) = 'Info: Testing dotp_self'
1327 call messages_info(1, namespace=sys%namespace)
1328
1329 call sys%st%group%psib(1, 1)%copy_to(xx, copy_data = .true.)
1330
1331 nst = sys%st%group%psib(1, 1)%nst
1332
1333 if (states_are_real(sys%st)) then
1334 safe_allocate(ddot(1:nst, 1:nst))
1335 do itime = 1, param%repetitions
1336 call dmesh_batch_dotp_self(sys%gr, xx, ddot)
1337 end do
1338
1339 do ist = 1, nst
1340 do jst = 1, nst
1341 write(message(jst), '(a,2i3,3x,e23.16)') 'Dotp_self states', ist, jst, ddot(ist,jst)
1342 end do
1343 call messages_info(nst, namespace=sys%namespace)
1344 end do
1345 safe_deallocate_a(ddot)
1346 else
1347 safe_allocate(zdot(1:nst, 1:nst))
1348 do itime = 1, param%repetitions
1349 call zmesh_batch_dotp_self(sys%gr, xx, zdot)
1350 end do
1351
1352 do ist = 1, nst
1353 do jst = 1, nst
1354 write(message(jst), '(a,2i3,3x,2(e23.16,1x))') 'Dotp_self states', ist, jst, zdot(ist,jst)
1355 end do
1356 call messages_info(nst, namespace=sys%namespace)
1357 end do
1358 safe_deallocate_a(zdot)
1359 end if
1360
1361 call xx%end()
1362 end if
1363
1364 if (bitand(ops, option__testbatchops__ops_ax_function_py) /= 0) then
1365 message(1) = 'Info: Testing ax_function_py'
1366 call messages_info(1, namespace=sys%namespace)
1367
1368 call sys%st%group%psib(1, 1)%copy_to(yy, copy_data = .true.)
1369
1370 if (states_are_real(sys%st)) then
1371 safe_allocate(df(sys%gr%np, sys%st%d%dim))
1372 safe_allocate(dweight(1:sys%st%group%psib(1, 1)%nst_linear))
1373 dweight = 0.1_real64
1374 call batch_get_state(sys%st%group%psib(1, 1), 1, sys%gr%np, df)
1375
1376 do itime = 1, param%repetitions
1377 call dbatch_ax_function_py(sys%gr%np, dweight, df, yy)
1378 end do
1379 safe_deallocate_a(df)
1380 safe_deallocate_a(dweight)
1381 else
1382 safe_allocate(zf(sys%gr%np, sys%st%d%dim))
1383 safe_allocate(zweight(1:sys%st%group%psib(1, 1)%nst_linear))
1384 zweight = cmplx(0.1_real64, m_zero, real64)
1385 call batch_get_state(sys%st%group%psib(1, 1), 1, sys%gr%np, zf)
1386
1387 do itime = 1, param%repetitions
1388 call zbatch_ax_function_py(sys%gr%np, zweight, zf, yy)
1389 end do
1390 safe_deallocate_a(zf)
1391 safe_deallocate_a(zweight)
1392 end if
1393 call test_prints_info_batch(sys%st, sys%gr, yy, string = "ax_function_py")
1394
1395 call yy%end()
1396 end if
1397
1398 if (bitand(ops, option__testbatchops__ops_get_points) /= 0) then
1399 ! The blocking over grid points mimics the actual use in the code
1400
1401 if (.not. accel_is_enabled()) then
1402
1403 call sys%st%group%psib(1, 1)%copy_to(yy)
1404
1405 if (states_are_real(sys%st)) then
1406
1407 block_size = max(40, cpu_hardware%l2%size/(2 * sizeof_real64 * sys%st%nst))
1408
1409 safe_allocate(dpoints(1:sys%st%nst, 1:sys%st%d%dim, 1:block_size))
1410
1411 do itime = 1, param%repetitions
1412 do sp = 1, sys%gr%np, block_size
1413 size = min(block_size, sys%gr%np - sp + 1)
1414 call batch_get_points(sys%st%group%psib(1, 1), sp, sp + size - 1, dpoints)
1415 call batch_set_points(yy, sp, sp + size - 1, dpoints)
1416 end do
1417 end do
1418 safe_deallocate_a(dpoints)
1419
1420 else
1421
1422 block_size = max(20, cpu_hardware%l2%size/(2 * sizeof_complex64 * sys%st%nst))
1423
1424 safe_allocate(zpoints(1:sys%st%nst, 1:sys%st%d%dim, 1:block_size))
1425
1426 do itime = 1, param%repetitions
1427 do sp = 1, sys%gr%np, block_size
1428 size = min(block_size, sys%gr%np - sp + 1)
1429 call batch_get_points(sys%st%group%psib(1, 1), sp, sp + size - 1, zpoints)
1430 call batch_set_points(yy, sp, sp + size - 1, zpoints)
1431 end do
1432 end do
1433 safe_deallocate_a(zpoints)
1434
1435 end if
1436
1437 call test_prints_info_batch(sys%st, sys%gr, yy, string = "get_points")
1438
1439 call yy%end()
1440 end if
1441 end if
1442
1443 call states_elec_deallocate_wfns(sys%st)
1444 safe_deallocate_p(sys)
1445
1446 pop_sub(test_batch_ops)
1447 end subroutine test_batch_ops
1448
1449! ---------------------------------------------------------
1450 subroutine test_derivatives(param, namespace, grp)
1451 type(test_parameters_t), intent(in) :: param
1452 type(namespace_t), intent(in) :: namespace
1453 type(mpi_grp_t), intent(in) :: grp
1454
1455 type(electrons_t), pointer :: sys
1456
1457 push_sub(test_derivatives)
1458
1459 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
1460
1461 message(1) = 'Info: Testing the finite-differences derivatives.'
1462 message(2) = ''
1463 call messages_info(2, namespace=namespace)
1464
1465 if (param%type == option__testtype__all .or. param%type == option__testtype__real) then
1466 call dderivatives_test(sys%gr%der, sys%namespace, param%repetitions, param%min_blocksize, param%max_blocksize)
1467 end if
1468
1469 if (param%type == option__testtype__all .or. param%type == option__testtype__complex) then
1470 call zderivatives_test(sys%gr%der, sys%namespace, param%repetitions, param%min_blocksize, param%max_blocksize)
1471 end if
1472
1473 if (sys%space%dim > 1) then
1474 call derivatives_advanced_benchmark(sys%gr%der, sys%namespace)
1475 end if
1476
1477 safe_deallocate_p(sys)
1478
1479 pop_sub(test_derivatives)
1480 end subroutine test_derivatives
1481
1482 ! ---------------------------------------------------------
1483
1484 subroutine test_orthogonalization(param, namespace, grp)
1485 type(test_parameters_t), intent(in) :: param
1486 type(namespace_t), intent(in) :: namespace
1487 type(mpi_grp_t), intent(in) :: grp
1488
1489 type(electrons_t), pointer :: sys
1490 integer :: itime
1491
1492 push_sub(test_orthogonalization)
1493
1494 call calc_mode_par%set_parallelization(p_strategy_states, default = .false.)
1495 call calc_mode_par%set_scalapack_compat()
1496
1497 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
1498
1499 message(1) = 'Info: Testing orthogonalization.'
1500 message(2) = ''
1501 call messages_info(2, namespace=namespace)
1502
1503 if (param%type == option__testtype__all .or. param%type == option__testtype__real) then
1504 message(1) = 'Info: Real wave-functions.'
1505 call messages_info(1, namespace=namespace)
1506 do itime = 1, param%repetitions
1507 call dstates_elec_calc_orth_test(sys%st, sys%namespace, sys%gr, sys%kpoints)
1508 end do
1509 end if
1510
1511 if (param%type == option__testtype__all .or. param%type == option__testtype__complex) then
1512 message(1) = 'Info: Complex wave-functions.'
1513 call messages_info(1, namespace=namespace)
1514 do itime = 1, param%repetitions
1515 call zstates_elec_calc_orth_test(sys%st, sys%namespace, sys%gr, sys%kpoints)
1516 end do
1517 end if
1518
1519 safe_deallocate_p(sys)
1520
1521 pop_sub(test_orthogonalization)
1522 end subroutine test_orthogonalization
1523
1524 ! ---------------------------------------------------------
1525
1526 subroutine test_interpolation(param, namespace, grp)
1527 type(test_parameters_t), intent(in) :: param
1528 type(namespace_t), intent(in) :: namespace
1529 type(mpi_grp_t), intent(in) :: grp
1530
1531 type(electrons_t), pointer :: sys
1532
1533 push_sub(test_interpolation)
1534
1535 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
1536
1537 if (param%type == option__testtype__all .or. param%type == option__testtype__real) then
1538 call messages_write('Info: Testing real interpolation routines')
1539 call messages_new_line()
1540 call messages_new_line()
1541 call messages_info(namespace=namespace)
1542
1543 call dmesh_interpolation_test(sys%gr)
1544 end if
1545
1546 if (param%type == option__testtype__all .or. param%type == option__testtype__complex) then
1547 call messages_new_line()
1548 call messages_write('Info: Testing complex interpolation routines')
1549 call messages_new_line()
1550 call messages_new_line()
1551 call messages_info(namespace=namespace)
1552
1553 call zmesh_interpolation_test(sys%gr)
1554 end if
1555
1556 safe_deallocate_p(sys)
1557
1558 pop_sub(test_interpolation)
1559 end subroutine test_interpolation
1560
1561
1562 ! ---------------------------------------------------------
1563
1564 subroutine test_ion_interaction(namespace, grp)
1565 type(namespace_t), intent(in) :: namespace
1566 type(mpi_grp_t), intent(in) :: grp
1567
1568 type(electrons_t), pointer :: sys
1569
1570 push_sub(test_ion_interaction)
1571
1572 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
1573
1574 call ion_interaction_test(sys%space, sys%ions%latt, sys%ions%atom, sys%ions%natoms, sys%ions%pos, &
1575 sys%gr%box%bounding_box_l, namespace, sys%mc)
1576
1577 safe_deallocate_p(sys)
1578
1579 pop_sub(test_ion_interaction)
1580 end subroutine test_ion_interaction
1581
1582 ! ---------------------------------------------------------
1583
1584 subroutine test_prints_info_batch(st, gr, psib, string)
1585 type(states_elec_t), intent(in) :: st
1586 type(grid_t), intent(in) :: gr
1587 class(batch_t), intent(inout) :: psib
1588 character(*), optional, intent(in) :: string
1589
1590 integer :: itime
1591 complex(real64), allocatable :: zpsi(:, :)
1592 real(real64), allocatable :: dpsi(:, :)
1593
1594 character(80) :: string_
1595
1596 string_ = optional_default(string, "")
1597
1598 push_sub(test_prints_info_batch)
1599
1600 if (states_are_real(st)) then
1601 safe_allocate(dpsi(1:gr%np, 1:st%d%dim))
1602 else
1603 safe_allocate(zpsi(1:gr%np, 1:st%d%dim))
1604 end if
1605
1606 do itime = 1, psib%nst
1607 if (states_are_real(st)) then
1608 call batch_get_state(psib, itime, gr%np, dpsi)
1609 write(message(1),'(a,i2,3x,e23.16)') "Norm state "//trim(string_)//" ", itime, dmf_nrm2(gr, st%d%dim, dpsi)
1610 else
1611 call batch_get_state(psib, itime, gr%np, zpsi)
1612 write(message(1),'(a,i2,3x,e23.16)') "Norm state "//trim(string_)//" ", itime, zmf_nrm2(gr, st%d%dim, zpsi)
1613 end if
1614 call messages_info(1)
1615 end do
1616
1617 if (states_are_real(st)) then
1618 safe_deallocate_a(dpsi)
1619 else
1620 safe_deallocate_a(zpsi)
1621 end if
1622
1623 pop_sub(test_prints_info_batch)
1624
1625 end subroutine test_prints_info_batch
1626
1627
1628 ! ---------------------------------------------------------
1629 subroutine test_clock()
1630
1631 type(iteration_counter_t) :: clock, other_clock
1632
1633 push_sub(test_clock)
1634
1635 ! Clock operations
1636 write(message(1),'(a)') " Operation Counter Time Global step"
1637 call messages_info(1)
1638
1639 clock = clock_t(time_step=m_two, initial_iteration=2)
1640 call write_clock("initialization")
1641
1642 clock = clock + 5
1643 call write_clock("addition (+5)")
1644
1645 clock = clock - 3
1646 call write_clock("subtraction (-3)")
1647
1648 call clock%reset()
1649 call write_clock("reset (=0)")
1650
1651 call clock%set(clock_t(time_step=m_four, initial_iteration=1))
1652 call write_clock("set (=2)")
1653
1654 ! Clock comparisons:
1655 write(message(1),'(a)') " Clock comparisons:"
1656 call messages_info(1)
1657
1658 other_clock = clock_t(time_step=m_one, initial_iteration=5)
1659
1660 call write_condition_result("4 < 5", clock < other_clock)
1661 call write_condition_result("4 <= 5", clock <= other_clock)
1662 call write_condition_result("4 > 5", clock > other_clock)
1663 call write_condition_result("4 >= 5", clock >= other_clock)
1664 call write_condition_result("4 == 5", clock == other_clock)
1665 call write_condition_result("4 /= 5", clock /= other_clock)
1666
1667 pop_sub(test_clock)
1668 contains
1669 subroutine write_clock(operation)
1670 character(len=*), intent(in) :: operation
1671
1672 write(message(1),'(a17,1x,i6,1x,f10.1,1x,i16)') operation, clock%counter(), clock%value(), clock%global_step()
1673 call messages_info(1)
1674 end subroutine write_clock
1675
1676 subroutine write_condition_result(condition, result)
1677 character(len=*), intent(in) :: condition
1678 logical, intent(in) :: result
1679
1680 write(message(1),'(a10," = ",i1," (",l1,")")') condition, abs(transfer(result, 0)), result
1681 call messages_info(1)
1682 end subroutine write_condition_result
1683
1684 end subroutine test_clock
1685
1686
1687 ! ---------------------------------------------------------
1688 subroutine test_cgal()
1689
1690 type(cgal_polyhedra_t) :: cgal_poly
1691
1692 push_sub(test_cgal)
1693
1694 call cgal_polyhedron_init(cgal_poly, "28-cgal.02-X.off", verbose = .true.)
1695
1696 if (cgal_polyhedron_point_inside(cgal_poly, 30._real64, 10._real64, 30._real64)) then
1697 message(1) = "cgal_polyhedron_point_inside"
1698 call messages_info(1)
1699 end if
1700
1701 call cgal_polyhedron_end(cgal_poly)
1702
1703 pop_sub(test_cgal)
1704 end subroutine test_cgal
1705
1706
1707 ! ---------------------------------------------------------
1708 subroutine test_dense_eigensolver()
1709 integer :: N, ii, jj, N_list(4), i_N
1710 real(real64), allocatable :: matrix(:, :), eigenvectors(:, :), eigenvalues(:), test(:)
1711 real(real64), allocatable :: differences(:)
1712
1713 push_sub(test_dense_eigensolver)
1714
1715 n_list = [15, 32, 100, 500]
1716
1717 do i_n = 1, 4
1718 n = n_list(i_n)
1719 safe_allocate(matrix(1:n, 1:n))
1720 safe_allocate(eigenvectors(1:n, 1:n))
1721 safe_allocate(eigenvalues(1:n))
1722 safe_allocate(test(1:n))
1723 safe_allocate(differences(1:n))
1724
1725
1726 ! set up symmetrix matrix
1727 do jj = 1, n
1728 do ii = 1, n
1729 matrix(ii, jj) = ii * jj
1730 end do
1731 end do
1732
1733 ! parallel solver
1734 eigenvectors(1:n, 1:n) = matrix(1:n, 1:n)
1735 call lalg_eigensolve_parallel(n, eigenvectors, eigenvalues)
1736
1737 do ii = 1, n
1738 test(:) = matmul(matrix, eigenvectors(:, ii)) - eigenvalues(ii) * eigenvectors(:, ii)
1739 differences(ii) = sum(abs(test)) / sum(abs(eigenvectors(:, ii)))
1740 end do
1741 write(message(1), "(A, I3, A, E13.6)") "Parallel solver - N: ", n, &
1742 ", average difference: ", sum(differences)/n
1743 call messages_info(1)
1744
1745 ! serial solver
1746 eigenvectors(1:n, 1:n) = matrix(1:n, 1:n)
1747 call lalg_eigensolve(n, eigenvectors, eigenvalues)
1748
1749 do ii = 1, n
1750 test(:) = matmul(matrix, eigenvectors(:, ii)) - eigenvalues(ii) * eigenvectors(:, ii)
1751 differences(ii) = sum(abs(test)) / sum(abs(eigenvectors(:, ii)))
1752 end do
1753 write(message(1), "(A, I3, A, E13.6)") "Serial solver - N: ", n, &
1754 ", average difference: ", sum(differences)/n
1755 call messages_info(1)
1756
1757 safe_deallocate_a(matrix)
1758 safe_deallocate_a(eigenvectors)
1759 safe_deallocate_a(eigenvalues)
1760 safe_deallocate_a(test)
1761 safe_deallocate_a(differences)
1762 end do
1763
1764 pop_sub(test_dense_eigensolver)
1765 end subroutine test_dense_eigensolver
1766
1767 subroutine test_batch_set_gaussian(psib, mesh)
1768 class(batch_t), intent(inout) :: psib
1769 class(mesh_t), intent(in) :: mesh
1770
1771 real(real64), allocatable :: dff(:)
1772 complex(real64), allocatable :: zff(:)
1773 integer :: ist, ip
1774 real(real64) :: da, db, dc
1775 complex(real64) :: za, zb, zc
1776
1777 push_sub(test_batch_set_gaussian)
1778
1779 ! use a similar function as in the derivatives test
1780 da = m_one/mesh%box%bounding_box_l(1)
1781 db = 10.0_real64
1782 dc = 100.0_real64
1783
1784 if (type_is_complex(psib%type())) then
1785 ! we make things more "complex"
1786 za = da + m_zi*0.01_real64
1787 zb = db*exp(m_zi*0.345_real64)
1788 zc = dc - m_zi*50.0_real64
1789
1790 safe_allocate(zff(1:mesh%np))
1791 do ist = 1, psib%nst_linear
1792 za = za * ist
1793 zb = zb / ist
1794 do ip = 1, mesh%np
1795 zff(ip) = zb*exp(-za*sum(mesh%x(:, ip)**2)) + zc
1796 end do
1797 call batch_set_state(psib, ist, mesh%np, zff)
1798 end do
1799 call zmesh_batch_normalize(mesh, psib)
1800 safe_deallocate_a(zff)
1801 else
1802 safe_allocate(dff(1:mesh%np))
1803 do ist = 1, psib%nst_linear
1804 da = da * ist
1805 db = db / ist
1806 do ip = 1, mesh%np
1807 dff(ip) = db*exp(-da*sum(mesh%x(:, ip)**2)) + dc
1808 end do
1809 call batch_set_state(psib, ist, mesh%np, dff)
1810 end do
1811 call dmesh_batch_normalize(mesh, psib)
1812 safe_deallocate_a(dff)
1813 end if
1814
1815 pop_sub(test_batch_set_gaussian)
1816 end subroutine test_batch_set_gaussian
1817
1818 ! ---------------------------------------------------------
1819 subroutine test_grid_interpolation(grp)
1820 type(mpi_grp_t), intent(in) :: grp
1821 type(electrons_t), pointer :: sys
1822 type(multigrid_t) :: mgrid
1823
1824 push_sub(test_grid_interpolation)
1825
1826 sys => electrons_t(global_namespace, grp, int(option__calculationmode__dummy, int32))
1827
1828 call multigrid_init(mgrid, global_namespace, sys%space, sys%gr, sys%gr%der, &
1829 sys%gr%stencil, sys%mc, nlevels=3)
1830
1831 call multigrid_test_interpolation(mgrid, sys%space)
1832
1833 call multigrid_end(mgrid)
1834
1835 safe_deallocate_p(sys)
1836
1837 pop_sub(test_grid_interpolation)
1838 end subroutine test_grid_interpolation
1839
1840 subroutine test_iihash()
1841
1842 type(iihash_t) :: h
1843 integer :: value
1844 logical :: found
1845
1846 call iihash_init(h)
1847 call iihash_insert(h, 1,5)
1848 call iihash_insert(h, 3,7)
1849
1850 value = iihash_lookup(h, 1, found)
1851 write(message(1),*) 'hash[1] :', found, value
1852 call messages_info()
1853
1854 value = iihash_lookup(h, 2, found)
1855 write(message(1),*) 'hash[2] :', found, value
1856 call messages_info()
1857
1858 value = iihash_lookup(h, 3, found)
1859 write(message(1),*) 'hash[3] :', found, value
1860 call messages_info()
1861
1862 call iihash_end(h)
1863
1864 end subroutine test_iihash
1865
1866 subroutine test_sihash()
1867 type(sihash_t) :: h
1868 type(sihash_iterator_t) :: it
1869
1870 integer :: counter
1871 integer :: value, sum
1872 logical :: found
1873
1874 call sihash_init(h)
1875 call sihash_insert(h, "one",5)
1876 call sihash_insert(h, "three",7)
1877 call sihash_insert(h, "the answer", 42)
1878
1879 value = sihash_lookup(h, "one", found)
1880 write(message(1),*) 'hash["one"]: ', found, value
1881 call messages_info()
1882
1883 value = sihash_lookup(h, "two", found)
1884 write(message(1),*) 'hash["two"]: ', found, value
1885 call messages_info()
1886
1887 value = sihash_lookup(h, "three", found)
1888 write(message(1),*) 'hash["three"]: ', found, value
1889 call messages_info()
1890
1891 sum = 0
1892 counter = 1
1893 call it%start(h)
1894
1895 do while (it%has_next())
1896 value = it%get_next()
1897 sum = sum + value
1898 write(message(1),'(I3,A,I5)') counter,': hash[...] = ',value
1899 call messages_info()
1900 counter = counter + 1
1901 end do
1902 write(message(1),*) 'counter = ', counter
1903 write(message(2),*) 'sum = ', sum
1904 call messages_info(2)
1905
1906
1907 call sihash_end(h)
1908
1909 end subroutine test_sihash
1910
1911
1912 subroutine test_sphash(namespace)
1913 type(namespace_t), intent(in) :: namespace
1914
1915 type(sphash_t) :: h
1916 type(sphash_iterator_t) :: it
1917
1918 logical :: found
1919
1920
1921 type(iteration_counter_t) :: clock_1
1922 type(iteration_counter_t), allocatable :: clock_2
1923 type(space_t) :: space_1
1924
1925 class(*), pointer :: value
1926
1927 integer :: count_clock, count_space
1928
1929 safe_allocate(clock_2)
1930
1931 clock_1 = clock_t(1.d-5, 1)
1932 clock_2 = clock_t(2.d-5, 2)
1933 space_1 = space_t(namespace)
1934
1935 call sphash_init(h)
1936 call sphash_insert(h, "one", clock_1)
1937 call sphash_insert(h, "two", space_1)
1938 call sphash_insert(h, "three", clock_2, clone=.true.)
1939
1940 value => sphash_lookup(h, "one", found)
1941 select type(value)
1942 type is (iteration_counter_t)
1943 write(message(1),*) 'hash["one"]: ', found, value%counter()
1944 call messages_info()
1945 type is (space_t)
1946 write(message(1),*) 'hash["one"]: ', found, value%short_info()
1947 call messages_info()
1948 class default
1949 write(message(1),*) 'wrong type. found = ', found
1950 call messages_info()
1951 end select
1952
1953 value => sphash_lookup(h, "two", found)
1954 select type(value)
1955 type is (iteration_counter_t)
1956 write(message(1),*) 'hash["two"]: ', found, value%counter()
1957 call messages_info()
1958 type is (space_t)
1959 write(message(1),*) 'hash["two"]: ', found, value%short_info()
1960 call messages_info()
1961 class default
1962 write(message(1),*) 'wrong type. found = ',found
1963 call messages_info()
1964 end select
1965
1966 safe_deallocate_a(clock_2)
1967
1968 value => sphash_lookup(h, "three", found)
1969 select type(value)
1970 type is (iteration_counter_t)
1971 write(message(1),*) 'hash["three"]: ', found, value%counter()
1972 call messages_info()
1973 type is (space_t)
1974 write(message(1),*) 'hash["three"]: ', found, value%short_info()
1975 call messages_info()
1976 class default
1977 write(message(1),*) 'wrong type. found = ',found
1978 call messages_info()
1979 end select
1980
1981 count_clock = 0
1982 count_space = 0
1983
1984 call it%start(h)
1985
1986 do while (it%has_next())
1987 value => it%get_next()
1988 select type(value)
1989 type is (iteration_counter_t)
1990 count_clock = count_clock + 1
1991 type is (space_t)
1992 count_space = count_space + 1
1993 end select
1994 end do
1995
1996 write(message(1), *) 'Count_clock = ', count_clock
1997 write(message(2), *) 'Count_space = ', count_space
1998 call messages_info(2)
1999
2000 call sphash_end(h)
2001
2002 end subroutine test_sphash
2003
2004 ! ---------------------------------------------------------
2005 subroutine test_regridding(namespace, grp)
2006 type(namespace_t), intent(in) :: namespace
2007 type(mpi_grp_t), intent(in) :: grp
2008
2009 type(electrons_t), pointer :: sysA, sysB
2010 type(regridding_t), pointer :: regridding
2011 real(real64), allocatable :: ff_A(:), ff_A_reference(:), ff_B(:), ff_B_reference(:), diff_A(:), diff_B(:)
2012 real(real64) :: norm_ff, norm_diff
2013 integer :: ip, ierr
2014
2015 push_sub(test_regridding)
2016
2017 call calc_mode_par%set_parallelization(p_strategy_states, default=.true.)
2018
2019 sysa => electrons_t(namespace_t("A", namespace), grp, int(option__calculationmode__dummy, int32))
2020 sysb => electrons_t(namespace_t("B", namespace), grp, int(option__calculationmode__dummy, int32))
2021
2022 safe_allocate(ff_a(1:sysa%gr%np))
2023 safe_allocate(ff_a_reference(1:sysa%gr%np))
2024 safe_allocate(diff_a(1:sysa%gr%np))
2025 safe_allocate(ff_b(1:sysb%gr%np))
2026 safe_allocate(ff_b_reference(1:sysb%gr%np))
2027 safe_allocate(diff_b(1:sysb%gr%np))
2028
2029 do ip = 1, sysa%gr%np
2030 ff_a_reference(ip) = values(sysa%gr%x(:, ip))
2031 end do
2032 do ip = 1, sysb%gr%np
2033 ff_b_reference(ip) = values(sysb%gr%x(:, ip))
2034 end do
2035
2036 ! forward mapping A->B
2037 regridding => regridding_t(sysb%gr, sysa%gr, sysa%space, sysa%namespace)
2038 call regridding%do_transfer(ff_b, ff_a_reference)
2039 safe_deallocate_p(regridding)
2040
2041 ! check that mapped function is exactly the same for the points that are available on both grids
2042 do ip = 1, sysb%gr%np
2043 if (abs(ff_b(ip)) <= m_epsilon) then
2044 ff_b_reference(ip) = m_zero
2045 diff_b(ip) = m_zero
2046 else
2047 diff_b(ip) = abs(ff_b_reference(ip) - ff_b(ip))
2048 end if
2049 end do
2050 norm_ff = dmf_nrm2(sysb%gr, ff_b_reference)
2051 norm_diff = dmf_nrm2(sysb%gr, diff_b)
2052
2053 write(message(1),'(a, E14.6)') "Forward: difference of reference to mapped function (rel.): ", &
2054 norm_diff/norm_ff
2055 call messages_info(1, namespace=namespace)
2056
2057 call dio_function_output(io_function_fill_how('AxisX'), ".", "forward_reference", namespace, sysb%space, &
2058 sysb%gr, ff_b_reference, unit_one, ierr)
2059 call dio_function_output(io_function_fill_how('AxisX'), ".", "forward_mapped", namespace, sysb%space, &
2060 sysb%gr, ff_b, unit_one, ierr)
2061 call dio_function_output(io_function_fill_how('AxisX'), ".", "forward_original", namespace, sysa%space, &
2062 sysa%gr, ff_a_reference, unit_one, ierr)
2063
2064 ! backward mapping B->A
2065 regridding => regridding_t(sysa%gr, sysb%gr, sysb%space, sysb%namespace)
2066 call regridding%do_transfer(ff_a, ff_b_reference)
2067 safe_deallocate_p(regridding)
2068
2069 do ip = 1, sysa%gr%np
2070 if (abs(ff_a(ip)) <= m_epsilon) then
2071 ff_a_reference(ip) = m_zero
2072 diff_a(ip) = m_zero
2073 else
2074 diff_a(ip) = abs(ff_a_reference(ip) - ff_a(ip))
2075 end if
2076 end do
2077 norm_ff = dmf_nrm2(sysa%gr, ff_a_reference)
2078 norm_diff = dmf_nrm2(sysa%gr, diff_a)
2079
2080 write(message(1),'(a, E14.6)') "Backward: difference of reference to mapped function (rel.): ", &
2081 norm_diff/norm_ff
2082 call messages_info(1, namespace=namespace)
2083
2084 call dio_function_output(io_function_fill_how('AxisX'), ".", "backward_reference", namespace, sysa%space, &
2085 sysa%gr, ff_a_reference, unit_one, ierr)
2086 call dio_function_output(io_function_fill_how('AxisX'), ".", "backward_mapped", namespace, sysa%space, &
2087 sysa%gr, ff_a, unit_one, ierr)
2088 call dio_function_output(io_function_fill_how('AxisX'), ".", "backward_original", namespace, sysb%space, &
2089 sysb%gr, ff_b_reference, unit_one, ierr)
2090
2091 safe_deallocate_a(ff_a)
2092 safe_deallocate_a(ff_a_reference)
2093 safe_deallocate_a(ff_b)
2094 safe_deallocate_a(ff_b_reference)
2095 safe_deallocate_a(diff_a)
2096 safe_deallocate_a(diff_b)
2097 safe_deallocate_p(sysa)
2098 safe_deallocate_p(sysb)
2099
2100 call messages_info(1, namespace=namespace)
2101
2102 pop_sub(test_regridding)
2103 contains
2104 real(real64) function values(xx)
2105 real(real64), intent(in) :: xx(:)
2106 real(real64) :: xx0(1:size(xx, dim=1))
2107 real(real64), parameter :: aa = m_half
2108 real(real64), parameter :: bb = m_four
2109
2110 ! no push_sub/pop_sub because this function is called often
2111 xx0(:) = m_one
2112 values = bb * exp(-aa*sum((xx-xx0)**2))
2113 end function values
2114 end subroutine test_regridding
2115
2116 ! ---------------------------------------------------------
2117 subroutine test_mesh_generation(namespace, grp)
2118 type(namespace_t), intent(in) :: namespace
2119 type(mpi_grp_t), intent(in) :: grp
2120
2121 type(namespace_t) :: namespace_a, namespace_b
2122 type(electrons_t), pointer :: sysA, sysB
2123 type(regridding_t), pointer :: regridding
2124 real(real64), allocatable :: ff_a(:), ff_b(:), lapl_a(:), lapl_b(:), lapl_b_on_a(:), diff_lapl(:)
2125 real(real64) :: abs_diff, rel_diff, norm_lapl
2126 integer :: ip, dim
2127 real(real64), parameter :: alpha = 0.5_real64
2128
2129 push_sub(test_mesh_generation)
2130
2131 call calc_mode_par%set_parallelization(p_strategy_states, default=.false.)
2132
2133 namespace_a = namespace_t("A", namespace)
2134 namespace_b = namespace_t("B", namespace)
2135
2136 sysa => electrons_t(namespace_a, grp, int(option__calculationmode__dummy, int32))
2137 sysb => electrons_t(namespace_b, grp, int(option__calculationmode__dummy, int32))
2138
2139 assert(sysa%gr%np_global == sysb%gr%np_global)
2140 assert(sysa%st%nst == sysb%st%nst)
2141 assert(sysa%ions%natoms == sysb%ions%natoms)
2142
2143 dim = sysa%gr%box%dim
2144 safe_allocate(ff_a(1:sysa%gr%np_part))
2145 safe_allocate(ff_b(1:sysb%gr%np_part))
2146 safe_allocate(lapl_a(1:sysa%gr%np))
2147 safe_allocate(lapl_b(1:sysb%gr%np))
2148 safe_allocate(lapl_b_on_a(1:sysa%gr%np))
2149 safe_allocate(diff_lapl(1:sysa%gr%np))
2150
2151 do ip = 1, sysa%gr%np_part
2152 ff_a(ip) = gaussian(sysa%gr%x(1:dim, ip), alpha)
2153 end do
2154 do ip = 1, sysb%gr%np_part
2155 ff_b(ip) = gaussian(sysb%gr%x(1:dim, ip), alpha)
2156 end do
2157
2158 call dderivatives_lapl(sysa%gr%der, ff_a, lapl_a, set_bc=.false.)
2159 call dderivatives_lapl(sysb%gr%der, ff_b, lapl_b, set_bc=.false.)
2160
2161 regridding => regridding_t(sysa%gr, sysb%gr, sysb%space, sysb%namespace)
2162 call regridding%do_transfer(lapl_b_on_a, lapl_b)
2163 safe_deallocate_p(regridding)
2164
2165 diff_lapl = lapl_a - lapl_b_on_a
2166 abs_diff = dmf_nrm2(sysa%gr, diff_lapl)
2167 norm_lapl = dmf_nrm2(sysa%gr, lapl_a)
2168 rel_diff = abs_diff/max(norm_lapl, m_epsilon)
2169
2170 assert(rel_diff <= 1.0e-12_real64)
2171
2172 write(message(1), '(a)') "System A initialized"
2173 call messages_info(1, namespace=namespace_a)
2174 write(message(1), '(a)') "System B initialized"
2175 call messages_info(1, namespace=namespace_b)
2176 write(message(1), '(a,es17.10)') "Gaussian Laplacian absolute difference = ", abs_diff
2177 write(message(2), '(a,es17.10)') "Gaussian Laplacian relative difference = ", rel_diff
2178 write(message(3), '(a,es17.10)') "Gaussian Laplacian reference norm = ", norm_lapl
2179 call messages_info(3, namespace=namespace)
2180
2181 safe_deallocate_a(ff_a)
2182 safe_deallocate_a(ff_b)
2183 safe_deallocate_a(lapl_a)
2184 safe_deallocate_a(lapl_b)
2185 safe_deallocate_a(lapl_b_on_a)
2186 safe_deallocate_a(diff_lapl)
2187 safe_deallocate_p(sysa)
2188 safe_deallocate_p(sysb)
2189
2190 pop_sub(test_mesh_generation)
2191 contains
2192 real(real64) function gaussian(x, a) result(val)
2193 real(real64), intent(in) :: x(:)
2194 real(real64), intent(in) :: a
2195
2196 val = exp(-a*sum(x**2))
2197 end function gaussian
2198
2199 end subroutine test_mesh_generation
2200
2206 subroutine test_vecpot_analytical(namespace, grp)
2207 type(namespace_t), intent(in) :: namespace
2208 type(mpi_grp_t), intent(in) :: grp
2209
2210 class(maxwell_t), pointer :: maxwell_system
2211
2212 real(real64), allocatable :: magnetic_field(:,:)
2213 real(real64), allocatable :: vector_potential_mag(:,:)
2214 real(real64), allocatable :: vector_potential_analytical(:,:)
2215 real(real64), allocatable :: delta(:,:)
2216 real(real64) :: exp_factor
2217 real(real64) :: xx
2218 real(real64) :: yy
2219 real(real64) :: zz
2220 real(real64) :: sigma
2221 integer :: ip, j, ierr, nn
2222 integer(int64) :: out_how
2223 character(len=MAX_PATH_LEN) :: fname, fname2, fname3
2224
2225 out_how = 32
2226 maxwell_system => maxwell_t(namespace, grp)
2227 sigma = maxwell_system%gr%box%bounding_box_l(1)/10_real64 ! this is exponential width
2228
2229 safe_allocate(magnetic_field(1:maxwell_system%gr%np_part, 1:3))
2230 safe_allocate(vector_potential_mag(1:maxwell_system%gr%np_part, 1:3))
2231 safe_allocate(vector_potential_analytical(1:maxwell_system%gr%np_part, 1:3))
2232 safe_allocate(delta(1:maxwell_system%gr%np, 1:3))
2233
2234 !%Variable TestVectorPotentialType
2235 !%Type integer
2236 !%Default bounded
2237 !%Section Calculation Modes::Test
2238 !%Description
2239 !% Select whether bounded or unbounded type will be used for vector potential tests
2240 !%Option bounded 1
2241 !% Analytical Vector Potential formulation is bounded by spatial gaussian
2242 !%Option unbounded 2
2243 !% Analytical Vector Potential is not bounded
2244 !%End
2245 call parse_variable(namespace, 'TestVectorPotentialType', option__testvectorpotentialtype__bounded, nn)
2246
2247 select case (nn)
2248 case (option__testvectorpotentialtype__bounded) ! bounded input
2249 do ip = 1, maxwell_system%gr%np_part
2250 xx = maxwell_system%gr%x(1, ip)
2251 yy = maxwell_system%gr%x(2, ip)
2252 zz = maxwell_system%gr%x(3, ip)
2253 exp_factor = exp((-xx**2 - yy**2 - zz**2)*1/(2*sigma**2))
2254 magnetic_field(ip, 1) = exp_factor*yy*(1 - (-xx**2 + yy**2)/(3*sigma**2) - zz**2/(3*sigma**2))
2255 magnetic_field(ip, 2) = exp_factor * xx * (1 + (-xx**2 + yy**2)/(3*sigma**2) - zz**2/(3*sigma**2))
2256 magnetic_field(ip, 3) = exp_factor * 2 * xx * yy * zz * 1/(3*sigma**2)
2257
2258 vector_potential_analytical(ip, 1) = m_third * xx * zz * exp_factor
2259 vector_potential_analytical(ip, 2) = - m_third * yy * zz * exp_factor
2260 vector_potential_analytical(ip, 3) = m_third * (-xx**2 + yy**2) * exp_factor
2261 end do
2262 case (option__testvectorpotentialtype__unbounded) ! unbounded input,
2263 ! TODO this unit test requires implementation of BCs for Helmholtz decomposition
2264 do ip = 1, maxwell_system%gr%np_part
2265 magnetic_field(ip, 1) = maxwell_system%gr%x(2, ip)
2266 magnetic_field(ip, 2) = maxwell_system%gr%x(1, ip)
2267 magnetic_field(ip, 3) = m_zero
2268
2269 vector_potential_analytical(ip, 1) = m_third * maxwell_system%gr%x(1, ip) * maxwell_system%gr%x(3, ip)
2270 vector_potential_analytical(ip, 2) = - m_third * maxwell_system%gr%x(2, ip) * maxwell_system%gr%x(3, ip)
2271 vector_potential_analytical(ip, 3) = - m_third * (maxwell_system%gr%x(1, ip)**2 - maxwell_system%gr%x(2, ip)**2)
2272 end do
2273 end select
2274 call maxwell_system%helmholtz%get_vector_potential(namespace, vector_potential_mag, magnetic_field)
2275
2276 do j = 1, 3
2277 delta(:,:) = m_zero
2278 do ip = 1, maxwell_system%gr%np
2279 delta(ip,j) = vector_potential_analytical(ip, j) - vector_potential_mag(ip, j)
2280 end do
2281 end do
2282
2283 do j = 1, 3
2284 write(message(j),*) 'j, norm2(delta)', j, norm2(delta(:,j))
2285 end do
2286 call messages_info(3)
2287
2288 write(fname, '(a)') 'deviation_from_analytical_formulation' ! use messages later
2289 call io_function_output_vector(out_how , './' , trim(fname), namespace, maxwell_system%space, maxwell_system%gr, &
2290 delta, unit_one, ierr)
2291 write(fname2, '(a)') 'vector_potential_analytical'
2292 call io_function_output_vector(out_how , './' , trim(fname2), namespace, maxwell_system%space, maxwell_system%gr, &
2293 vector_potential_analytical, unit_one, ierr)
2294 write(fname3, '(a)') 'vector_potential_mag'
2295 call io_function_output_vector(out_how , './' , trim(fname3), namespace, maxwell_system%space, maxwell_system%gr, &
2296 vector_potential_mag, unit_one, ierr)
2297
2298 safe_deallocate_a(magnetic_field)
2299 safe_deallocate_a(vector_potential_mag)
2300 safe_deallocate_a(vector_potential_analytical)
2301 safe_deallocate_a(delta)
2302 safe_deallocate_p(maxwell_system)
2303
2304 end subroutine test_vecpot_analytical
2305
2306 ! ---------------------------------------------------------
2307 subroutine multigrid_test_interpolation(mgrid, space)
2308 type(multigrid_t), intent(in) :: mgrid
2309 class(space_t), intent(in) :: space
2310
2311 real(real64), allocatable :: guess0(:), res0(:), guess1(:)
2312 type(mesh_t), pointer :: mesh0, mesh1
2313 real(real64) :: delta, xx(3,2), alpha, beta, rr
2314 integer :: nn, ip, ierr
2315
2316 push_sub(multigrid_test_interpolation)
2317
2318 message(1) = 'Info: Testing the grid interpolation.'
2319 message(2) = ''
2320 call messages_info(2)
2321
2322 mesh0 => mgrid%level(0)%mesh
2323 mesh1 => mgrid%level(1)%mesh
2324
2325 safe_allocate(guess0(1:mesh0%np_part))
2326 safe_allocate(res0(1:mesh0%np_part))
2327 safe_allocate(guess1(1:mesh1%np_part))
2328
2329 alpha = m_four*mesh0%spacing(1)
2330 beta = m_one / (alpha**space%dim * sqrt(m_pi)**space%dim)
2331
2332 ! Set the centers of the Gaussians by hand
2333 xx(1, 1) = m_one
2334 xx(2, 1) = -m_half
2335 xx(3, 1) = m_two
2336 xx(1, 2) = -m_two
2337 xx(2, 2) = m_zero
2338 xx(3, 2) = -m_one
2339 xx = xx * alpha
2340
2341 ! Density as sum of Gaussians
2342 guess0 = m_zero
2343 do nn = 1, 2
2344 do ip = 1, mesh0%np
2345 call mesh_r(mesh0, ip, rr, origin = xx(:, nn))
2346 guess0(ip) = guess0(ip) + (-1)**nn * beta*exp(-(rr/alpha)**2)
2347 end do
2348 end do
2349
2350 call dio_function_output (io_function_fill_how('AxisX'), ".", "interpolation_target", global_namespace, &
2351 space, mesh0, guess0, unit_one, ierr)
2352 call dio_function_output (io_function_fill_how('AxisZ'), ".", "interpolation_target", global_namespace, &
2353 space, mesh0, guess0, unit_one, ierr)
2354 call dio_function_output (io_function_fill_how('PlaneZ'), ".", "interpolation_target", global_namespace, &
2355 space, mesh0, guess0, unit_one, ierr)
2356
2357 ! We start by testing the interpolation scheme. For this, we generate a function on the fine grid
2358 ! and we inject it on the coarse grid. Then we interpolate it back to the fine grid and we compare
2359 ! This allows for testing the quality of the interpolation scheme
2360
2361 ! move to level 1
2362 call dmultigrid_fine2coarse(mgrid%level(1)%tt, mgrid%level(0)%der, mesh1, guess0, guess1, injection)
2363 ! back to level 0
2364 call dmultigrid_coarse2fine(mgrid%level(1)%tt, mgrid%level(1)%der, mesh0, guess1, res0)
2365
2366 call dio_function_output (io_function_fill_how('AxisX'), ".", "interpolation_result", global_namespace, &
2367 space, mesh0, res0, unit_one, ierr)
2368 call dio_function_output (io_function_fill_how('AxisZ'), ".", "interpolation_result", global_namespace, &
2369 space, mesh0, res0, unit_one, ierr)
2370 call dio_function_output (io_function_fill_how('PlaneZ'), ".", "interpolation_result", global_namespace, &
2371 space, mesh0, res0, unit_one, ierr)
2372
2373 delta = dmf_nrm2(mesh0, guess0(1:mesh0%np)-res0(1:mesh0%np))
2374 write(message(1),'(a,e13.6)') 'Interpolation test (abs.) = ', delta
2375
2376 ! Now we test if the restriction+interpolation combination returns the original result or not
2377 ! This allows to test if restriction and interpolation are adjoint operators or not.
2378
2379 ! move to level 1
2380 call dmultigrid_fine2coarse(mgrid%level(1)%tt, mgrid%level(0)%der, mesh1, guess0, guess1, fullweight)
2381 ! back to level 0
2382 call dmultigrid_coarse2fine(mgrid%level(1)%tt, mgrid%level(1)%der, mesh0, guess1, res0)
2383
2384 call dio_function_output (io_function_fill_how('AxisX'), ".", "restriction_result", global_namespace, &
2385 space, mesh0, res0, unit_one, ierr)
2386 call dio_function_output (io_function_fill_how('AxisZ'), ".", "restriction_result", global_namespace, &
2387 space, mesh0, res0, unit_one, ierr)
2388 call dio_function_output (io_function_fill_how('PlaneZ'), ".", "restriction_result", global_namespace, &
2389 space, mesh0, res0, unit_one, ierr)
2390
2391 delta = dmf_nrm2(mesh0, guess0(1:mesh0%np)-res0(1:mesh0%np))
2392 write(message(2),'(a,e13.6)') 'Restriction test (abs.) = ', delta
2393 call messages_info(2)
2394
2395 safe_deallocate_a(guess0)
2396 safe_deallocate_a(res0)
2397 safe_deallocate_a(guess1)
2398
2399 pop_sub(multigrid_test_interpolation)
2400 end subroutine multigrid_test_interpolation
2401
2402
2404 subroutine test_current_density(namespace, grp)
2405 type(namespace_t), intent(in) :: namespace
2406 type(mpi_grp_t), intent(in) :: grp
2407 type(electrons_t), pointer :: sys
2408
2409 type(current_t) :: current
2410 character(len=MAX_PATH_LEN) :: fname
2411 integer :: ierr, ip, idir
2412 integer(int64) :: out_how
2413 real(real64), allocatable :: current_para_ref(:,:), current_dia_ref(:,:), current_mag_ref(:,:), delta(:)
2414 real(real64) :: xx(3), rr, a0, mag_curr, sin_thet, sin_phi, cos_phi, vec_pot_slope
2415 complex(real64) :: alpha
2416
2417 sys => electrons_t(namespace, grp, int(option__calculationmode__dummy, int32))
2418
2419 alpha = (0.0_real64, 0.5_real64)
2420 a0 = m_one
2421 vec_pot_slope = 0.4_real64 ! units of B field
2422
2423 call states_elec_allocate_wfns(sys%st, sys%gr, wfs_type = type_cmplx)
2424 call set_hydrogen_states(sys%st%group%psib(1, 1), sys%gr, namespace, alpha, a0)
2425
2426 call current_init(current, namespace)
2427
2428 call hamiltonian_elec_epot_generate(sys%hm, sys%namespace, sys%space, sys%gr, sys%ions, sys%ext_partners, sys%st)
2429
2430 safe_allocate(sys%hm%hm_base%vector_potential(1:3, 1:sys%gr%np))
2431
2432 sys%hm%hm_base%vector_potential = m_zero
2433 do ip = 1, sys%gr%np
2434 xx = sys%gr%x(1:3, ip)
2435 sys%hm%hm_base%vector_potential(2, ip) = vec_pot_slope * xx(1) / p_c ! vector potential is here devided by c_0
2436 end do
2437
2438 call states_elec_allocate_current(sys%st, sys%space, sys%gr)
2439 call density_calc(sys%st, sys%gr, sys%st%rho)
2440 ! paramagnetic + diamagnetic current densities
2441 call current_calculate(current, namespace, sys%gr, sys%hm, sys%space, sys%st)
2442
2443 safe_allocate(current_para_ref(1:sys%gr%np,1:3))
2444 safe_allocate(current_dia_ref(1:sys%gr%np,1:3))
2445 safe_allocate(current_mag_ref(1:sys%gr%np,1:3))
2446 safe_allocate(delta(1:sys%gr%np))
2447
2448 ! analytic paramagnetic
2449 current_para_ref(:,:) = m_zero
2450 do ip = 1, sys%gr%np
2451 call mesh_r(sys%gr, ip, rr)
2452 xx = sys%gr%x(1:3, ip)
2453 if (rr > r_small) then
2454 current_para_ref(ip,1:3) = - (psi_2s(rr, a0) * dr_psi_1s(rr, a0) - &
2455 psi_1s(rr, a0) * dr_psi_2s(rr, a0) ) * aimag(alpha) / (1 + abs(alpha)**2) * xx(1:3) / rr
2456 end if
2457 end do
2458
2459 write(fname, '(a)') 'current_para'
2460 out_how = io_function_fill_how("PlaneZ")
2461 call io_function_output_vector(out_how , './' , trim(fname), namespace, sys%space, sys%gr, &
2462 sys%st%current_para(:,:,1), unit_one, ierr)
2463
2464 write(fname, '(a)') 'current_para-ref'
2465 out_how = io_function_fill_how("PlaneZ")
2466 call io_function_output_vector(out_how , './' , trim(fname), namespace, sys%space, sys%gr, &
2467 current_para_ref(:,:), unit_one, ierr)
2468
2469 do idir = 1, 3
2470 delta = m_zero
2471 delta(:) = current_para_ref(1:sys%gr%np, idir) - sys%st%current_para(1:sys%gr%np, idir, 1)
2472 write(message(idir),*) 'idir =',idir,', norm2(delta paramagnetic)',norm2(delta)
2473 end do
2474 call messages_info(3)
2475
2476 ! analytic diamagnetic
2477 current_dia_ref(:,:) = m_zero
2478 do ip = 1, sys%gr%np
2479 call mesh_r(sys%gr, ip, rr)
2480 current_dia_ref(ip,1:3) = - sys%hm%hm_base%vector_potential(1:3,ip) *&
2481 m_one / (1 + abs(alpha)**2) * abs(lc_hydrogen_state(rr, alpha, a0))**2
2482 end do
2483
2484 write(fname, '(a)') 'current_dia'
2485 out_how = io_function_fill_how("PlaneZ")
2486 call io_function_output_vector(out_how , './' , trim(fname), namespace, sys%space, sys%gr, &
2487 sys%st%current_dia(:,:,1), unit_one, ierr)
2488
2489 write(fname, '(a)') 'current_dia-ref'
2490 out_how = io_function_fill_how("PlaneZ")
2491 call io_function_output_vector(out_how , './' , trim(fname), namespace, sys%space, sys%gr, &
2492 current_dia_ref(:,:), unit_one, ierr)
2493
2494 do idir = 1, 3
2495 delta = m_zero
2496 delta(:) = current_dia_ref(1:sys%gr%np, idir) - sys%st%current_dia(1:sys%gr%np, idir, 1)
2497 write(message(idir),*) 'idir =',idir,', norm2(delta diamagnetic)',norm2(delta)
2498 end do
2499 call messages_info(3)
2500
2501 ! magnetization current
2502 call current_calculate_mag(sys%gr%der, sys%st)
2503
2504 ! analytic magnetization
2505 current_mag_ref(:,:) = m_zero
2506 do ip = 1, sys%gr%np
2507 call mesh_r(sys%gr, ip, rr)
2508 xx = sys%gr%x(1:3, ip)
2509 if (norm2(xx(1:2)) > r_small .and. rr > r_small) then
2510 sin_thet = norm2(xx(1:2)) / rr
2511 sin_phi = xx(2) / norm2(xx(1:2))
2512 cos_phi = xx(1) / norm2(xx(1:2))
2513 mag_curr = m_two * (psi_1s(rr,a0)*dr_psi_1s(rr,a0) + real(alpha)*(psi_1s(rr,a0)*dr_psi_2s(rr,a0) &
2514 + dr_psi_1s(rr,a0)*psi_2s(rr,a0)) + abs(alpha)**2 * psi_2s(rr,a0)*dr_psi_2s(rr,a0))
2515 ! minus signs are reversed because of the charge of the electron
2516 current_mag_ref(ip,1) = m_half * mag_curr * sin_thet * sin_phi / (1+abs(alpha)**2)
2517 current_mag_ref(ip,2) = -m_half * mag_curr * sin_thet * cos_phi / (1+abs(alpha)**2)
2518 endif
2519 end do
2520
2521 write(fname, '(a)') 'current_mag'
2522 out_how = io_function_fill_how("PlaneZ")
2523 call io_function_output_vector(out_how , './' , trim(fname), namespace, sys%space, sys%gr, &
2524 sys%st%current_mag(:,:,1), unit_one, ierr)
2525
2526 write(fname, '(a)') 'current_mag-ref'
2527 out_how = io_function_fill_how("PlaneZ")
2528 call io_function_output_vector(out_how , './' , trim(fname), namespace, sys%space, sys%gr, &
2529 current_mag_ref(:,:), unit_one, ierr)
2530
2531 do idir = 1, 3
2532 delta = m_zero
2533 delta(:) = current_mag_ref(1:sys%gr%np, idir) - sys%st%current_mag(1:sys%gr%np, idir, 1)
2534 write(message(idir),*) 'idir =',idir,', norm2(delta magnetization)',norm2(delta)
2535 end do
2536 call messages_info(3)
2537
2538 safe_deallocate_a(current_para_ref)
2539 safe_deallocate_a(current_dia_ref)
2540 safe_deallocate_a(current_mag_ref)
2541 safe_deallocate_a(delta)
2542 safe_deallocate_p(sys)
2543
2544 end subroutine test_current_density
2545
2546
2547 subroutine set_hydrogen_states(psib, mesh, namespace, alpha, a0)
2548 class(batch_t), intent(inout) :: psib
2549 class(mesh_t), intent(in) :: mesh
2550 type(namespace_t), intent(in) :: namespace
2551 complex(real64), intent(in) :: alpha
2552 real(real64), intent(in) :: a0
2553
2554 complex(real64), allocatable :: zff(:)
2555 integer :: ip
2556 real(real64) :: rr
2557
2558 push_sub(set_hydrogen_states)
2559
2560 ! psi = phi_{1s} + \alpha \phi_{2s}
2561 safe_allocate(zff(1:mesh%np))
2562 if (type_is_complex(psib%type())) then
2563 do ip = 1, mesh%np
2564 call mesh_r(mesh, ip, rr)
2565 zff(ip) = m_one / sqrt(1 + abs(alpha)**2) * lc_hydrogen_state(rr, alpha, a0)
2566 end do
2567 call batch_set_state(psib, 1, mesh%np, zff)
2568 else
2569 write(message(1),*) "States should be complex for the linear combination of hydrogenic states to work"
2570 call messages_info(1, namespace=namespace)
2571 end if
2572 safe_deallocate_a(zff)
2573
2574 pop_sub(set_hydrogen_states)
2575 end subroutine set_hydrogen_states
2576
2577 complex(real64) function lc_hydrogen_state(rr, alpha, a0)
2578 real(real64), intent(in) :: a0, rr
2579 complex(real64), intent(in) :: alpha
2580
2581 lc_hydrogen_state = psi_1s(rr, a0) + alpha * psi_2s(rr, a0)
2582 end function lc_hydrogen_state
2583
2584 ! phi_{1s} = R_{10}(r)Y_{00}(\theta, \phi) = 1/(2\sqrt(\pi)) * 2 a0^{-3/2} * \exp(-r/a0)
2585 real(real64) function psi_1s(rr, a0)
2586 real(real64), intent(in) :: a0, rr
2587
2588 psi_1s = a0**(-m_three/m_two) * exp(-rr / a0) / sqrt(m_pi)
2589 end function psi_1s
2590
2591 ! phi_{2s} = R_{20}(r)Y_{00}(\theta, \phi) = 1/(2\sqrt(\pi)) * \sqrt(2) a0^{-3/2} (1 - r/(2a0)) \exp(-r/(2a0))
2592 real(real64) function psi_2s(rr, a0)
2593 real(real64), intent(in) :: a0, rr
2594
2595 psi_2s = sqrt(m_two) * a0**(-m_three/m_two) * (m_one - rr/(m_two * a0)) * exp(-rr/(m_two * a0)) &
2596 / (m_two * sqrt(m_pi))
2597 end function psi_2s
2598
2599 real(real64) function dr_psi_1s(rr, a0)
2600 real(real64), intent(in) :: a0, rr
2601
2602 dr_psi_1s = -(m_one / a0) * psi_1s(rr, a0)
2603 end function dr_psi_1s
2604
2605 real(real64) function dr_psi_2s(rr, a0)
2606 real(real64), intent(in) :: a0, rr
2607
2608 dr_psi_2s = -(m_half / a0) * psi_2s(rr, a0) - &
2609 a0**(-m_five/m_two) * exp(-rr/(m_two * a0)) / (sqrt(m_two) * m_two * sqrt(m_pi))
2610 end function dr_psi_2s
2611
2612
2613 subroutine test_csv_input(namespace)
2614 type(namespace_t), intent(in) :: namespace
2615
2616 integer :: ierr
2617 integer(int64) :: i, j, k
2618 integer(int64) :: dims(3)
2619 character(len=MAX_PATH_LEN) :: fname
2620
2621 real(real64), allocatable :: ff(:)
2622
2623 push_sub(test_csv_input)
2624
2625 !%Variable TestCSVFileName
2626 !%Type string
2627 !%Section Calculation Modes::Test
2628 !%Description
2629 !% The file name of the csv file to be read
2630 !%End
2631 call parse_variable(namespace, 'TestCSVFileName', "", fname)
2632
2633 message(1) = "Attempting to probe "//trim(fname)
2634 call messages_info(1, namespace=namespace)
2635
2636 call io_csv_get_info(fname, dims, ierr)
2637 if (ierr == 0) then
2638 message(1) = "Probing successful."
2639 write(message(2), '("found dimensions: ",3I20)') dims
2640 call messages_info(2, namespace=namespace)
2641 else
2642 write(message(1),'("Probing failed. ierr = ",I5)') ierr
2643 call messages_fatal(1, namespace=namespace)
2644 end if
2645
2646 safe_allocate(ff(1:dims(1)*dims(2)*dims(3)))
2647
2648 message(1) = "Attempting to read "//trim(fname)
2649 call messages_info(1, namespace=namespace)
2650
2651 call dread_csv(fname, dims(1)*dims(2)*dims(3), ff, ierr)
2652 if (ierr == 0) then
2653 message(1) = "Reading successful."
2654 call messages_info(1, namespace=namespace)
2655
2656 do k=1, min(4_int64, dims(3))
2657 do j=1, min(4_int64, dims(2))
2658 write(message(j), '("data ",2I5, 1X, 4F8.2)') k, j, &
2659 (ff(i + dims(1)*(j-1) + dims(1)* dims(2)*(k-1)), i=1, min(4_int64, dims(1)))
2660 end do
2661 write(message(int(j, int32)), '("")')
2662 call messages_info(int(j, int32), namespace=namespace)
2663 end do
2664
2665 else
2666 message(1) = "Reading failed."
2667 call messages_fatal(1, namespace=namespace)
2668 end if
2669
2670 safe_deallocate_a(ff)
2671
2672 pop_sub(test_csv_input)
2673 end subroutine test_csv_input
2674
2675end module test_oct_m
2676
2677!! Local Variables:
2678!! mode: f90
2679!! coding: utf-8
2680!! End:
batch_ops_oct_m::batch_axpby
batchified axpby:
Definition: batch_ops.F90:175
batch_ops_oct_m::batch_axpy
batchified version of the BLAS axpy routine:
Definition: batch_ops.F90:159
batch_ops_oct_m::batch_get_points
Definition: batch_ops.F90:236
batch_ops_oct_m::batch_get_state
Definition: batch_ops.F90:227
batch_ops_oct_m::batch_mul_mf
batchified multiplication by mesh function with optional conjugation:
Definition: batch_ops.F90:254
batch_ops_oct_m::batch_pointwise_mul
Definition: batch_ops.F90:248
batch_ops_oct_m::batch_scal2v
batchified scale with optional conjugation:
Definition: batch_ops.F90:181
batch_ops_oct_m::batch_scal
scale a batch by a constant or vector
Definition: batch_ops.F90:167
batch_ops_oct_m::batch_set_points
Definition: batch_ops.F90:242
batch_ops_oct_m::batch_set_state
There are several ways how to call batch_set_state and batch_get_state:
Definition: batch_ops.F90:218
boundaries_oct_m::boundaries_set
Definition: boundaries.F90:217
clock_oct_m::clock_t
Definition: clock.F90:149
global_oct_m::bitand
Definition: global.F90:318
global_oct_m::optional_default
Definition: global.F90:299
lalg_adv_oct_m::lalg_eigensolve_parallel
Definition: lalg_adv.F90:196
lalg_adv_oct_m::lalg_eigensolve
Definition: lalg_adv.F90:188
mesh_function_oct_m::dmf_nrm2
Definition: mesh_function.F90:199
mesh_function_oct_m::zmf_nrm2
Definition: mesh_function.F90:203
messages_oct_m::messages_print_var_option
Definition: messages.F90:199
messages_oct_m::messages_print_var_value
Prints out to iunit a message in the form: ["InputVariable" = value] where "InputVariable" is given b...
Definition: messages.F90:182
messages_oct_m::messages_write
Definition: messages.F90:190
parser_oct_m::parse_variable
Definition: parser.F90:253
solvers_oct_m::dconjugate_gradients
Definition: solvers.F90:212
exp
double exp(double __x) __attribute__((__nothrow__
sqrt
double sqrt(double __x) __attribute__((__nothrow__
accel_oct_m
Definition: accel.F90:120
accel_oct_m::accel_is_enabled
pure logical function, public accel_is_enabled()
Definition: accel.F90:403
batch_oct_m
This module implements batches of mesh functions.
Definition: batch.F90:135
batch_ops_oct_m
This module implements common operations on batches of mesh functions.
Definition: batch_ops.F90:118
batch_ops_oct_m::dbatch_ax_function_py
subroutine, public dbatch_ax_function_py(np, aa, psi, yy)
This routine performs a set of axpy operations adding the same function psi to all functions of a bat...
Definition: batch_ops.F90:1198
batch_ops_oct_m::batch_set_zero
subroutine, public batch_set_zero(this, np, async)
fill all mesh functions of the batch with zero
Definition: batch_ops.F90:265
batch_ops_oct_m::zbatch_ax_function_py
subroutine, public zbatch_ax_function_py(np, aa, psi, yy)
This routine performs a set of axpy operations adding the same function psi to all functions of a bat...
Definition: batch_ops.F90:3099
boundaries_oct_m
Module implementing boundary conditions in Octopus.
Definition: boundaries.F90:124
box_factory_oct_m
Module, implementing a factory for boxes.
Definition: box_factory.F90:122
box_oct_m
Definition: box.F90:116
calc_mode_par_oct_m
This module handles the calculation mode.
Definition: calc_mode_par.F90:124
calc_mode_par_oct_m::p_strategy_kpoints
integer, parameter, public p_strategy_kpoints
parallelization in k-points
Definition: calc_mode_par.F90:136
calc_mode_par_oct_m::p_strategy_domains
integer, parameter, public p_strategy_domains
parallelization in domains
Definition: calc_mode_par.F90:136
calc_mode_par_oct_m::calc_mode_par
type(calc_mode_par_t), public calc_mode_par
Singleton instance of parallel calculation mode.
Definition: calc_mode_par.F90:181
calc_mode_par_oct_m::p_strategy_states
integer, parameter, public p_strategy_states
parallelization in states
Definition: calc_mode_par.F90:136
cgal_polyhedra_oct_m
Definition: cgal_polyhedra.F90:115
cgal_polyhedra_oct_m::cgal_polyhedron_init
subroutine, public cgal_polyhedron_init(cgal_poly, fname, verbose)
Definition: cgal_polyhedra.F90:147
cgal_polyhedra_oct_m::cgal_polyhedron_point_inside
logical function, public cgal_polyhedron_point_inside(cgal_poly, xq, yq, zq)
Definition: cgal_polyhedra.F90:220
cgal_polyhedra_oct_m::cgal_polyhedron_end
subroutine, public cgal_polyhedron_end(cgal_poly)
Definition: cgal_polyhedra.F90:237
chebyshev_filter_bounds_oct_m
Definition: chebyshev_filter_bounds.F90:104
chebyshev_filter_bounds_oct_m::center
pure real(real64) function center(this)
Center of the filter interval.
Definition: chebyshev_filter_bounds.F90:205
clock_oct_m
Definition: clock.F90:116
current_oct_m
Definition: current.F90:116
debug_oct_m
Definition: debug.F90:116
density_oct_m
This module implements a calculator for the density and defines related functions.
Definition: density.F90:122
density_oct_m::density_calc
subroutine, public density_calc(st, gr, density, istin)
Computes the density from the orbitals in st.
Definition: density.F90:653
derivatives_oct_m
This module calculates the derivatives (gradients, Laplacians, etc.) of a function.
Definition: derivatives.F90:123
derivatives_test_oct_m
This module provides unit tests for the derivatives.
Definition: derivatives_test.F90:118
derivatives_test_oct_m::zderivatives_test
subroutine, public zderivatives_test(this, namespace, repetitions, min_blocksize, max_blocksize)
unit test for derivatives
Definition: derivatives_test.F90:727
derivatives_test_oct_m::dderivatives_test
subroutine, public dderivatives_test(this, namespace, repetitions, min_blocksize, max_blocksize)
unit test for derivatives
Definition: derivatives_test.F90:373
derivatives_test_oct_m::derivatives_advanced_benchmark
subroutine, public derivatives_advanced_benchmark(this, namespace)
Further unit tests design to challenge numerical stability of the finite differences.
Definition: derivatives_test.F90:150
eigen_chebyshev_oct_m
Definition: eigen_chebyshev.F90:104
eigen_chebyshev_oct_m::zchebyshev_filter
subroutine, public zchebyshev_filter(namespace, mesh, st, hm, degree, bounds, ik, normalize)
Chebyshev Filter.
Definition: eigen_chebyshev.F90:743
eigen_chebyshev_oct_m::dchebyshev_filter
subroutine, public dchebyshev_filter(namespace, mesh, st, hm, degree, bounds, ik, normalize)
Chebyshev Filter.
Definition: eigen_chebyshev.F90:398
electron_space_oct_m
Definition: electron_space.F90:116
electron_space_oct_m::spin_polarized
integer, parameter, public spin_polarized
Definition: electron_space.F90:141
electrons_oct_m
Definition: electrons.F90:119
exponential_oct_m
Definition: exponential.F90:117
exponential_oct_m::exponential_init
subroutine, public exponential_init(te, namespace, full_batch)
Definition: exponential.F90:212
global_oct_m
Definition: global.F90:116
global_oct_m::m_two
real(real64), parameter, public m_two
Definition: global.F90:202
global_oct_m::m_zero
real(real64), parameter, public m_zero
Definition: global.F90:200
global_oct_m::m_four
real(real64), parameter, public m_four
Definition: global.F90:204
global_oct_m::m_pi
real(real64), parameter, public m_pi
some mathematical constants
Definition: global.F90:198
global_oct_m::m_zi
complex(real64), parameter, public m_zi
Definition: global.F90:214
global_oct_m::m_epsilon
real(real64), parameter, public m_epsilon
Definition: global.F90:216
global_oct_m::m_half
real(real64), parameter, public m_half
Definition: global.F90:206
global_oct_m::m_one
real(real64), parameter, public m_one
Definition: global.F90:201
grid_oct_m
This module implements the underlying real-space grid.
Definition: grid.F90:119
hamiltonian_elec_base_oct_m
Definition: hamiltonian_elec_base.F90:116
hamiltonian_elec_oct_m
Definition: hamiltonian_elec.F90:117
hamiltonian_elec_oct_m::zhamiltonian_elec_apply_batch
subroutine, public zhamiltonian_elec_apply_batch(hm, namespace, mesh, psib, hpsib, terms, set_bc)
Definition: hamiltonian_elec.F90:2594
hamiltonian_elec_oct_m::hamiltonian_elec_end
subroutine, public hamiltonian_elec_end(hm)
Definition: hamiltonian_elec.F90:768
hamiltonian_elec_oct_m::hamiltonian_elec_epot_generate
subroutine, public hamiltonian_elec_epot_generate(this, namespace, space, gr, ions, ext_partners, st, time)
Definition: hamiltonian_elec.F90:1522
hamiltonian_elec_oct_m::hamiltonian_elec_copy
subroutine, public hamiltonian_elec_copy(hm_out, hm_in)
Deep-copy a hamiltonian_elec_t snapshot.
Definition: hamiltonian_elec.F90:902
hamiltonian_elec_oct_m::dhamiltonian_elec_apply_batch
subroutine, public dhamiltonian_elec_apply_batch(hm, namespace, mesh, psib, hpsib, terms, set_bc)
Definition: hamiltonian_elec.F90:2098
hardware_oct_m
Definition: hardware.F90:24
hardware_oct_m::cpu_hardware
type(hardware_t), public cpu_hardware
Global instance of CPU hardware specification.
Definition: hardware.F90:61
hardware_oct_m::sizeof_real64
integer, parameter, public sizeof_real64
Number of bytes to store a variable of type real(real64)
Definition: hardware.F90:55
hardware_oct_m::sizeof_complex64
integer, parameter, public sizeof_complex64
Number of bytes to store a variable of type complex(real64)
Definition: hardware.F90:58
helmholtz_decomposition_m
The Helmholtz decomposition is intended to contain "only mathematical" functions and procedures to co...
Definition: helmholtz_decomposition.F90:158
helmholtz_decomposition_test_m
Test suit for the Helmholtz decomposition module.
Definition: helmholtz_decomposition_test.F90:119
helmholtz_decomposition_test_m::gaussian_test
subroutine, public gaussian_test(helmholtz, sys_grid, namespace, space)
Definition: helmholtz_decomposition_test.F90:349
helmholtz_decomposition_test_m::hertzian_dipole_test
subroutine, public hertzian_dipole_test(helmholtz, sys_grid, namespace, space)
Definition: helmholtz_decomposition_test.F90:211
iihash_oct_m
This module implements a simple hash table for non-negative integer keys and integer values.
Definition: iihash.F90:127
iihash_oct_m::iihash_end
subroutine, public iihash_end(h)
Free a hash table.
Definition: iihash.F90:186
iihash_oct_m::iihash_insert
subroutine, public iihash_insert(h, key, val)
Insert a (key, val) pair into the hash table h.
Definition: iihash.F90:208
iihash_oct_m::iihash_lookup
integer function, public iihash_lookup(h, key, found)
Look up a value in the hash table h. If found is present, it indicates if key could be found in the t...
Definition: iihash.F90:233
iihash_oct_m::iihash_init
subroutine, public iihash_init(h)
Initialize a hash table h.
Definition: iihash.F90:163
io_csv_oct_m
Definition: io_csv_f.F90:118
io_function_oct_m
Definition: io_function.F90:116
io_function_oct_m::dio_function_output
subroutine, public dio_function_output(how, dir, fname, namespace, space, mesh, ff, unit, ierr, pos, atoms, grp, root)
Top-level IO routine for functions defined on the mesh.
Definition: io_function.F90:1677
io_function_oct_m::io_function_fill_how
integer(int64) function, public io_function_fill_how(where)
Use this function to quickly plot functions for debugging purposes: call dio_function_output(io_funct...
Definition: io_function.F90:817
io_oct_m
Definition: io.F90:116
ion_interaction_oct_m
Definition: ion_interaction.F90:117
ion_interaction_oct_m::ion_interaction_test
subroutine, public ion_interaction_test(space, latt, atom, natoms, pos, lsize, namespace, mc)
Definition: ion_interaction.F90:1392
iteration_counter_oct_m
Definition: iteration_counter.F90:27
lalg_adv_oct_m
Definition: lalg_adv.F90:116
lalg_adv_tests_oct_m
This modules takes care of testing some linear algebra routines.
Definition: lalg_adv_tests.F90:117
lalg_adv_tests_oct_m::test_exponential_matrix
subroutine, public test_exponential_matrix(namespace)
Unit tests for the exponential of a matrix.
Definition: lalg_adv_tests.F90:139
lalg_basic_oct_m
Definition: lalg_basic.F90:116
lapl_operator_oct_m
Computes and , suitable as an operator callback for iterative solvers (CG, QMR, etc....
Definition: lapl_operator.F90:119
lapl_operator_oct_m::dshifted_laplacian_op
subroutine, public dshifted_laplacian_op(x, lx, userdata)
Computes the shifted Laplacian operator: .
Definition: lapl_operator.F90:427
maxwell_oct_m
Definition: maxwell.F90:116
mesh_batch_oct_m
This module defines functions over batches of mesh functions.
Definition: mesh_batch.F90:118
mesh_batch_oct_m::dmesh_batch_dotp_matrix
subroutine, public dmesh_batch_dotp_matrix(mesh, aa, bb, dot, reduce)
Calculate the overlap matrix of two batches.
Definition: mesh_batch.F90:275
mesh_batch_oct_m::dmesh_batch_dotp_self
subroutine, public dmesh_batch_dotp_self(mesh, aa, dot, reduce)
calculate the overlap matrix of a batch with itself
Definition: mesh_batch.F90:539
mesh_batch_oct_m::mesh_batch_nrm2
subroutine, public mesh_batch_nrm2(mesh, aa, nrm2, reduce)
Calculate the norms (norm2, not the square!) of a batch of mesh functions.
Definition: mesh_batch.F90:178
mesh_batch_oct_m::zmesh_batch_dotp_vector
subroutine, public zmesh_batch_dotp_vector(mesh, aa, bb, dot, reduce, cproduct)
A simple switch between specialized kernels and generic kernels.
Definition: mesh_batch.F90:2087
mesh_batch_oct_m::zmesh_batch_dotp_matrix
subroutine, public zmesh_batch_dotp_matrix(mesh, aa, bb, dot, reduce)
Calculate the overlap matrix of two batches.
Definition: mesh_batch.F90:1707
mesh_batch_oct_m::zmesh_batch_normalize
subroutine, public zmesh_batch_normalize(mesh, psib, norm)
Normalize a batch.
Definition: mesh_batch.F90:3071
mesh_batch_oct_m::zmesh_batch_dotp_self
subroutine, public zmesh_batch_dotp_self(mesh, aa, dot, reduce)
calculate the overlap matrix of a batch with itself
Definition: mesh_batch.F90:1971
mesh_batch_oct_m::dmesh_batch_normalize
subroutine, public dmesh_batch_normalize(mesh, psib, norm)
Normalize a batch.
Definition: mesh_batch.F90:1601
mesh_batch_oct_m::dmesh_batch_dotp_vector
subroutine, public dmesh_batch_dotp_vector(mesh, aa, bb, dot, reduce, cproduct)
A simple switch between specialized kernels and generic kernels.
Definition: mesh_batch.F90:655
mesh_function_oct_m
This module defines various routines, operating on mesh functions.
Definition: mesh_function.F90:118
mesh_function_oct_m::dmf_dotp_aux
real(real64) function, public dmf_dotp_aux(f1, f2)
dot product between two vectors (mesh functions)
Definition: mesh_function.F90:1201
mesh_function_oct_m::mesh_init_mesh_aux
subroutine, public mesh_init_mesh_aux(mesh)
Initialise a pointer to the grid/mesh, that is globally exposed, such that low level mesh operations ...
Definition: mesh_function.F90:225
mesh_interpolation_oct_m
Definition: mesh_interpolation.F90:116
mesh_interpolation_oct_m::zmesh_interpolation_test
subroutine, public zmesh_interpolation_test(mesh)
Definition: mesh_interpolation.F90:1033
mesh_interpolation_oct_m::dmesh_interpolation_test
subroutine, public dmesh_interpolation_test(mesh)
Definition: mesh_interpolation.F90:549
mesh_oct_m
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:120
mesh_oct_m::mesh_r
pure subroutine, public mesh_r(mesh, ip, rr, origin, coords)
return the distance to the origin for a given grid point
Definition: mesh.F90:342
messages_oct_m
Definition: messages.F90:117
messages_oct_m::messages_print_with_emphasis
subroutine, public messages_print_with_emphasis(msg, iunit, namespace)
Definition: messages.F90:898
messages_oct_m::msg
character(len=512), private msg
Definition: messages.F90:167
messages_oct_m::messages_obsolete_variable
subroutine, public messages_obsolete_variable(namespace, name, rep)
Definition: messages.F90:1000
messages_oct_m::messages_new_line
subroutine, public messages_new_line()
Definition: messages.F90:1089
messages_oct_m::message
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:162
messages_oct_m::messages_fatal
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
Definition: messages.F90:410
messages_oct_m::messages_info
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:594
minimizer_tests_oct_m
This modules takes care of testing optimizers using standard test functions.
Definition: minimizer_tests.F90:117
minimizer_tests_oct_m::test_optimizers
subroutine, public test_optimizers(namespace)
Unit tests for different optimizers.
Definition: minimizer_tests.F90:149
mix_tests_oct_m
This module implements unit tests for the mixing methods.
Definition: mix_tests.F90:117
mix_tests_oct_m::mix_tests_run
subroutine, public mix_tests_run()
Definition: mix_tests.F90:138
mpi_oct_m
Definition: mpi.F90:116
mpi_test_oct_m
Definition: mpi_test.F90:116
mpi_test_oct_m::test_mpiwrappers
subroutine, public test_mpiwrappers
Definition: mpi_test.F90:130
multicomm_oct_m
This module handles the communicators for the various parallelization strategies.
Definition: multicomm.F90:147
multigrid_oct_m
Definition: multigrid.F90:117
multigrid_oct_m::multigrid_end
subroutine, public multigrid_end(mgrid)
Definition: multigrid.F90:518
multigrid_oct_m::multigrid_init
subroutine, public multigrid_init(mgrid, namespace, space, mesh, der, stencil, mc, nlevels)
Definition: multigrid.F90:187
namespace_oct_m
Definition: namespace.F90:105
namespace_oct_m::global_namespace
type(namespace_t), public global_namespace
Definition: namespace.F90:135
orbitalbasis_oct_m
Definition: orbitalbasis.F90:115
orbitalbasis_oct_m::orbitalbasis_end
subroutine, public orbitalbasis_end(this)
Definition: orbitalbasis.F90:294
orbitalbasis_oct_m::zorbitalbasis_build
subroutine, public zorbitalbasis_build(this, namespace, ions, mesh, kpt, ndim, skip_s_orb, use_all_orb, verbose)
This routine is an interface for constructing the orbital basis.
Definition: orbitalbasis.F90:916
orbitalbasis_oct_m::dorbitalbasis_build
subroutine, public dorbitalbasis_build(this, namespace, ions, mesh, kpt, ndim, skip_s_orb, use_all_orb, verbose)
This routine is an interface for constructing the orbital basis.
Definition: orbitalbasis.F90:417
orbitalbasis_oct_m::orbitalbasis_init
subroutine, public orbitalbasis_init(this, namespace, periodic_dim)
Definition: orbitalbasis.F90:179
orbitalset_oct_m
Definition: orbitalset.F90:116
orbitalset_oct_m::dorbitalset_add_to_batch
subroutine, public dorbitalset_add_to_batch(os, ndim, psib, weight)
Definition: orbitalset.F90:864
orbitalset_oct_m::zorbitalset_add_to_batch
subroutine, public zorbitalset_add_to_batch(os, ndim, psib, weight)
Definition: orbitalset.F90:1845
orbitalset_oct_m::orbitalset_update_phase
subroutine, public orbitalset_update_phase(os, dim, kpt, kpoints, spin_polarized, vec_pot, vec_pot_var, kpt_max)
Build the phase correction to the global phase in case the orbital crosses the border of the simulato...
Definition: orbitalset.F90:285
orbitalset_oct_m::dorbitalset_get_coeff_batch
subroutine, public dorbitalset_get_coeff_batch(os, ndim, psib, dot, reduce)
Definition: orbitalset.F90:589
orbitalset_oct_m::zorbitalset_get_coeff_batch
subroutine, public zorbitalset_get_coeff_batch(os, ndim, psib, dot, reduce)
Definition: orbitalset.F90:1517
parser_oct_m
Definition: parser.F90:116
poisson_test_oct_m
Definition: poisson_test.F90:117
poisson_test_oct_m::poisson_test
subroutine, public poisson_test(this, space, mesh, latt, namespace, repetitions)
This routine checks the Hartree solver selected in the input file by calculating numerically and anal...
Definition: poisson_test.F90:152
preconditioners_oct_m
Definition: preconditioners.F90:116
preconditioners_oct_m::preconditioner_end
subroutine, public preconditioner_end(this)
Definition: preconditioners.F90:356
preconditioners_oct_m::preconditioner_init
subroutine, public preconditioner_init(this, namespace, gr, mc, space)
Definition: preconditioners.F90:181
profiling_oct_m
Definition: profiling.F90:118
projector_oct_m
Definition: projector.F90:116
projector_oct_m::zproject_psi_batch
subroutine, public zproject_psi_batch(mesh, bnd, pj, npj, dim, psib, ppsib)
To optimize the application of the non-local operator in parallel, the projectors are applied in step...
Definition: projector.F90:1256
regridding_oct_m
Implementation details for regridding.
Definition: regridding.F90:172
sihash_oct_m
This module implements a simple hash table for string valued keys and integer values using the C++ ST...
Definition: sihash.F90:120
sihash_oct_m::sihash_insert
subroutine, public sihash_insert(h, key, val)
Insert a (key, val) pair into the hash table h.
Definition: sihash.F90:205
sihash_oct_m::sihash_init
subroutine, public sihash_init(h)
Initialize a hash table h with size entries. Since we use separate chaining, the number of entries in...
Definition: sihash.F90:164
sihash_oct_m::sihash_lookup
integer function, public sihash_lookup(h, key, found)
Look up a value in the hash table h. If found is present, it indicates if key could be found in the t...
Definition: sihash.F90:231
sihash_oct_m::sihash_end
subroutine, public sihash_end(h)
Free a hash table.
Definition: sihash.F90:185
solvers_oct_m
This module is intended to contain "only mathematical" functions and procedures.
Definition: solvers.F90:117
space_oct_m
Definition: space.F90:116
sphash_oct_m
This module implements a simple hash table for string valued keys and integer values using the C++ ST...
Definition: sphash.F90:120
sphash_oct_m::sphash_init
subroutine, public sphash_init(h)
Initialize a hash table h with size entries. Since we use separate chaining, the number of entries in...
Definition: sphash.F90:223
sphash_oct_m::sphash_insert
subroutine, public sphash_insert(h, key, val, clone)
Insert a (key, val) pair into the hash table h. If clone=.true., the object will be copied.
Definition: sphash.F90:293
sphash_oct_m::sphash_end
subroutine, public sphash_end(h)
Free a hash table.
Definition: sphash.F90:249
sphash_oct_m::sphash_lookup
class(*) function, pointer, public sphash_lookup(h, key, found)
Look up a value in the hash table h. If found is present, it indicates if key could be found in the t...
Definition: sphash.F90:325
states_abst_oct_m
Definition: states_abst.F90:115
states_abst_oct_m::states_are_real
pure logical function, public states_are_real(st)
Definition: states_abst.F90:210
states_elec_calc_oct_m
Definition: states_elec_calc.F90:116
states_elec_calc_oct_m::zstates_elec_calc_orth_test
subroutine, public zstates_elec_calc_orth_test(st, namespace, mesh, kpoints)
Definition: states_elec_calc.F90:3124
states_elec_calc_oct_m::dstates_elec_calc_orth_test
subroutine, public dstates_elec_calc_orth_test(st, namespace, mesh, kpoints)
Definition: states_elec_calc.F90:1286
states_elec_dim_oct_m
This module handles spin dimensions of the states and the k-point distribution.
Definition: states_elec_dim.F90:122
states_elec_oct_m
Definition: states_elec.F90:115
states_elec_oct_m::states_elec_deallocate_wfns
subroutine, public states_elec_deallocate_wfns(st)
Deallocates the KS wavefunctions defined within a states_elec_t structure.
Definition: states_elec.F90:1332
states_elec_oct_m::states_elec_allocate_wfns
subroutine, public states_elec_allocate_wfns(st, mesh, wfs_type, skip, packed)
Allocates the KS wavefunctions defined within a states_elec_t structure.
Definition: states_elec.F90:1124
states_mxll_oct_m
Definition: states_mxll.F90:115
submesh_oct_m
Definition: submesh.F90:116
subspace_oct_m
Definition: subspace.F90:116
subspace_oct_m::dsubspace_diag
subroutine, public dsubspace_diag(this, namespace, mesh, st, hm, ik, eigenval, diff, nonortho)
Diagonalises the Hamiltonian in the subspace defined by the states.
Definition: subspace.F90:298
subspace_oct_m::zsubspace_diag
subroutine, public zsubspace_diag(this, namespace, mesh, st, hm, ik, eigenval, diff, nonortho)
Diagonalises the Hamiltonian in the subspace defined by the states.
Definition: subspace.F90:769
test_isdf_oct_m
Integration tests for ISDF.
Definition: test_isdf.F90:116
test_isdf_oct_m::test_isdf
subroutine, public test_isdf(namespace, serial)
Set up an electron system, compute some optimal centroid positions, and use these to build a set of I...
Definition: test_isdf.F90:156
test_kmeans_clustering_oct_m
Definition: test_kmeans_clustering.F90:116
test_kmeans_clustering_oct_m::test_weighted_kmeans
subroutine, public test_weighted_kmeans(namespace)
Test weighted kmeans algorithm for a finite system.
Definition: test_kmeans_clustering.F90:145
test_oct_m
This module implements a unit-test like runmode for Octopus.
Definition: test.F90:117
test_oct_m::psi_2s
real(real64) function psi_2s(rr, a0)
Definition: test.F90:2688
test_oct_m::test_sihash
subroutine test_sihash()
Definition: test.F90:1962
test_oct_m::test_cgal
subroutine test_cgal()
Definition: test.F90:1784
test_oct_m::test_ion_interaction
subroutine test_ion_interaction(namespace, grp)
Definition: test.F90:1660
test_oct_m::test_hartree
subroutine test_hartree(param, namespace, grp)
Definition: test.F90:402
test_oct_m::test_interpolation
subroutine test_interpolation(param, namespace, grp)
Definition: test.F90:1622
test_oct_m::test_helmholtz_decomposition
subroutine test_helmholtz_decomposition(namespace, grp)
Definition: test.F90:421
test_oct_m::test_iihash
subroutine test_iihash()
Definition: test.F90:1936
test_oct_m::test_current_density
subroutine test_current_density(namespace, grp)
Here we test the different contributions to the total electronic current density.
Definition: test.F90:2500
test_oct_m::test_derivatives
subroutine test_derivatives(param, namespace, grp)
Definition: test.F90:1546
test_oct_m::psi_1s
real(real64) function psi_1s(rr, a0)
Definition: test.F90:2681
test_oct_m::test_batch_set_gaussian
subroutine test_batch_set_gaussian(psib, mesh)
Definition: test.F90:1863
test_oct_m::test_boundaries
subroutine test_boundaries(param, namespace, grp)
Definition: test.F90:875
test_oct_m::test_sphash
subroutine test_sphash(namespace)
Definition: test.F90:2008
test_oct_m::test_hamiltonian
subroutine test_hamiltonian(param, namespace, grp)
Definition: test.F90:681
test_oct_m::test_dense_eigensolver
subroutine test_dense_eigensolver()
Definition: test.F90:1804
test_oct_m::multigrid_test_interpolation
subroutine multigrid_test_interpolation(mgrid, space)
Definition: test.F90:2403
test_oct_m::test_projector
subroutine test_projector(param, namespace, grp)
Definition: test.F90:515
test_oct_m::dr_psi_2s
real(real64) function dr_psi_2s(rr, a0)
Definition: test.F90:2701
test_oct_m::test_composition_chebyshev
subroutine test_composition_chebyshev(namespace, grp)
Test the composition rule for Chebyshev polynomials.
Definition: test.F90:914
test_oct_m::lc_hydrogen_state
complex(real64) function lc_hydrogen_state(rr, alpha, a0)
Definition: test.F90:2673
test_oct_m::test_subspace_diagonalization
subroutine test_subspace_diagonalization(param, namespace, grp)
Definition: test.F90:1057
test_oct_m::test_run
subroutine, public test_run(namespace, grp)
Components and integration test runner.
Definition: test.F90:215
test_oct_m::test_density_calc
subroutine test_density_calc(param, namespace, grp)
Definition: test.F90:837
test_oct_m::test_linear_solver
subroutine test_linear_solver(namespace, grp)
Definition: test.F90:452
test_oct_m::test_prints_info_batch
subroutine test_prints_info_batch(st, gr, psib, string)
Definition: test.F90:1680
test_oct_m::test_orthogonalization
subroutine test_orthogonalization(param, namespace, grp)
Definition: test.F90:1580
test_oct_m::dr_psi_1s
real(real64) function dr_psi_1s(rr, a0)
Definition: test.F90:2695
test_oct_m::test_regridding
subroutine test_regridding(namespace, grp)
Definition: test.F90:2101
test_oct_m::test_batch_ops
subroutine test_batch_ops(param, namespace, grp)
Definition: test.F90:1110
test_oct_m::test_vecpot_analytical
subroutine test_vecpot_analytical(namespace, grp)
Here, analytical formulation for vector potential and B field are used. Ref: Sangita Sen and Erik I....
Definition: test.F90:2302
test_oct_m::test_exponential
subroutine test_exponential(param, namespace, grp)
Definition: test.F90:1008
test_oct_m::test_mesh_generation
subroutine test_mesh_generation(namespace, grp)
Definition: test.F90:2213
test_oct_m::set_hydrogen_states
subroutine set_hydrogen_states(psib, mesh, namespace, alpha, a0)
Definition: test.F90:2643
test_oct_m::test_dft_u
subroutine test_dft_u(param, namespace, grp)
Definition: test.F90:567
test_oct_m::test_grid_interpolation
subroutine test_grid_interpolation(grp)
Definition: test.F90:1915
test_oct_m::test_csv_input
subroutine test_csv_input(namespace)
Definition: test.F90:2709
test_oct_m::test_clock
subroutine test_clock()
Definition: test.F90:1725
test_parameters_oct_m
Definition: test_parameters.F90:118
test_poisson_fft_batch_oct_m
Definition: test_poisson_fft_batch.F90:115
test_poisson_fft_batch_oct_m::test_poisson_fft_batch
subroutine, public test_poisson_fft_batch(namespace)
Definition: test_poisson_fft_batch.F90:167
types_oct_m
Definition: types.F90:116
types_oct_m::type_cmplx
type(type_t), parameter, public type_cmplx
Definition: types.F90:136
types_oct_m::type_is_complex
logical pure function, public type_is_complex(this)
Definition: types.F90:182
unit_system_oct_m
This module defines the unit system, used for input and output.
Definition: unit_system.F90:128
unit_system_oct_m::unit_one
type(unit_t), public unit_one
some special units required for particular quantities
Definition: unit_system.F90:166
v_ks_oct_m
Definition: v_ks.F90:116
v_ks_oct_m::v_ks_calc
subroutine, public v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_eigenval, time, calc_energy, calc_current, force_semilocal)
Definition: v_ks.F90:708
wfs_elec_oct_m
Definition: wfs_elec.F90:116
batch_oct_m::batch_t
Class defining batches of mesh functions.
Definition: batch.F90:161
cgal_polyhedra_oct_m::cgal_polyhedra_t
Definition: cgal_polyhedra.F90:137
chebyshev_filter_bounds_oct_m::chebyshev_filter_bounds_t
Overload default constructor.
Definition: chebyshev_filter_bounds.F90:137
electrons_oct_m::electrons_t
Class describing the electron system.
Definition: electrons.F90:222
exponential_oct_m::exponential_t
Definition: exponential.F90:164
grid_oct_m::grid_t
Description of the grid, containing information on derivatives, stencil, and symmetries.
Definition: grid.F90:171
hamiltonian_elec_oct_m::hamiltonian_elec_t
Definition: hamiltonian_elec.F90:226
helmholtz_decomposition_m::helmholtz_decomposition_t
Definition: helmholtz_decomposition.F90:189
iihash_oct_m::iihash_t
Definition: iihash.F90:148
iteration_counter_oct_m::iteration_counter_t
This class implements the iteration counter used by the multisystem algorithms. As any iteration coun...
Definition: iteration_counter.F90:66
mesh_oct_m::mesh_t
Describes mesh distribution to nodes.
Definition: mesh.F90:187
mpi_oct_m::mpi_grp_t
This is defined even when running serial.
Definition: mpi.F90:144
multigrid_oct_m::multigrid_t
Definition: multigrid.F90:171
namespace_oct_m::namespace_t
Definition: namespace.F90:117
orbitalbasis_oct_m::orbitalbasis_t
Definition: orbitalbasis.F90:153
regridding_oct_m::regridding_t
contains the information of the meshes and provides the transfer functions
Definition: regridding.F90:198
sihash_oct_m::sihash_iterator_t
Definition: sihash.F90:143
sihash_oct_m::sihash_t
Definition: sihash.F90:137
space_oct_m::space_t
Definition: space.F90:132
sphash_oct_m::sphash_iterator_t
Definition: sphash.F90:163
sphash_oct_m::sphash_t
Definition: sphash.F90:154
states_elec_oct_m::states_elec_t
The states_elec_t class contains all electronic wave functions.
Definition: states_elec.F90:228
subspace_oct_m::subspace_t
Definition: subspace.F90:162
test_parameters_oct_m::test_parameters_t
Definition: test_parameters.F90:130
wfs_elec_oct_m::wfs_elec_t
batches of electronic states
Definition: wfs_elec.F90:141
true
int true(void)
Definition: symmetries_finite.c:3150
write_clock
subroutine write_clock(operation)
Definition: test.F90:1765
gaussian
real(real64) function gaussian(x, a)
Definition: test.F90:2288
write_condition_result
subroutine write_condition_result(condition, result)
Definition: test.F90:1772
values
real(real64) function values(xx)
Definition: test.F90:2200