Octopus
external_potential.F90
Go to the documentation of this file.
1!! Copyright (C) 2020 N. Tancogne-Dejean
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#include "global.h"
19
21 use debug_oct_m
22 use global_oct_m
23 use iihash_oct_m
31 use mesh_oct_m
34 use parser_oct_m
38 use space_oct_m
39 use string_oct_m
40 use unit_oct_m
43 implicit none
44
45 private
46 public :: &
49
51 private
52
53 type(space_t) :: space
54
55 integer, public :: type
56
57 character(len=1024) :: potential_formula
58 character(len=200) :: density_formula
59 character(len=MAX_PATH_LEN) :: filename
60 real(real64) :: omega
61
62 real(real64), allocatable, public :: pot(:)
63
64 real(real64), allocatable, public :: b_field(:)
65 integer :: gauge_2D
66 real(real64), allocatable, public :: a_static(:,:)
67 real(real64), allocatable, public :: e_field(:)
68 !Auxiliary arrays for the electrons only
69 !TODO: Suppress once electrons fully use the new framework
70 real(real64), allocatable, public :: v_ext(:)
71
72 contains
73 procedure :: calculate => external_potential_calculate
74 procedure :: allocate_memory => external_potential_allocate
75 procedure :: deallocate_memory => external_potential_deallocate
76 procedure :: init_interaction_as_partner => external_potential_init_interaction_as_partner
77 procedure :: copy_quantities_to_interaction => external_potential_copy_quantities_to_interaction
80
81 integer, public, parameter :: &
82 EXTERNAL_POT_USDEF = 201, & !< user-defined function for local potential
87
88
89 interface external_potential_t
90 module procedure external_potential_init
91 end interface external_potential_t
92
93contains
94
95 function external_potential_init(namespace) result(this)
96 class(external_potential_t), pointer :: this
97 type(namespace_t), intent(in) :: namespace
98
100
101 allocate(this)
102
103 this%namespace = namespace_t("ExternalPotential", parent=namespace)
104 this%space = space_t(namespace)
105
106 this%type = -1
107
108 allocate(this%supported_interactions_as_partner(0))
109
110 call this%quantities%add(quantity_t("E field", always_available = .true., updated_on_demand = .false., &
111 iteration = iteration_counter_t()))
112 call this%quantities%add(quantity_t("B field", always_available = .true., updated_on_demand = .false., &
113 iteration = iteration_counter_t()))
114
116 end function external_potential_init
117
118 ! ---------------------------------------------------------
119 subroutine external_potential_finalize(this)
120 type(external_potential_t), intent(inout) :: this
121
123
124 call this%deallocate_memory()
125
127 end subroutine external_potential_finalize
128
129 ! ---------------------------------------------------------
130 subroutine external_potential_allocate(this, mesh)
131 class(external_potential_t), intent(inout) :: this
132 class(mesh_t), intent(in) :: mesh
133
135
136 select case (this%type)
137 case (external_pot_usdef, external_pot_from_file, external_pot_charge_density)
138 safe_allocate(this%pot(1:mesh%np))
140 safe_allocate(this%a_static(1:mesh%np, 1:this%space%dim))
142 if (this%space%periodic_dim < this%space%dim) then
143 safe_allocate(this%pot(1:mesh%np))
144 safe_allocate(this%v_ext(1:mesh%np_part))
145 end if
146 end select
147
149 end subroutine external_potential_allocate
151 ! ---------------------------------------------------------
153 class(external_potential_t), intent(inout) :: this
154
156
157 safe_deallocate_a(this%pot)
158 safe_deallocate_a(this%b_field)
159 safe_deallocate_a(this%a_static)
160 safe_deallocate_a(this%e_field)
161 safe_deallocate_a(this%v_ext)
162
164 end subroutine external_potential_deallocate
165
166 ! ---------------------------------------------------------
167 subroutine external_potential_init_interaction_as_partner(partner, interaction)
168 class(external_potential_t), intent(in) :: partner
169 class(interaction_surrogate_t), intent(inout) :: interaction
172
173 select type (interaction)
175 ! Nothing to be initialized for the Lorentz force.
176 class default
177 message(1) = "Unsupported interaction."
178 call messages_fatal(1, namespace=partner%namespace)
179 end select
180
183
184 ! ---------------------------------------------------------
185 subroutine external_potential_copy_quantities_to_interaction(partner, interaction)
186 class(external_potential_t), intent(inout) :: partner
187 class(interaction_surrogate_t), intent(inout) :: interaction
189 integer :: ip
190
192
193 select type (interaction)
194 type is (lorentz_force_t)
195 if (partner%type == external_pot_static_efield) then
196 do ip = 1, interaction%system_np
197 interaction%partner_e_field(:, ip) = partner%e_field
198 interaction%partner_b_field(:, ip) = m_zero
199 end do
200 else if (partner%type == external_pot_static_bfield) then
201 do ip = 1, interaction%system_np
202 interaction%partner_e_field(:, ip) = m_zero
203 interaction%partner_b_field(:, ip) = partner%b_field
204 end do
205 else
206 assert(.false.) !This should never occur.
207 end if
208
209 class default
210 message(1) = "Unsupported interaction."
211 call messages_fatal(1, namespace=partner%namespace)
212 end select
213
216
217
218 ! ---------------------------------------------------------
219 subroutine external_potential_calculate(this, namespace, mesh, poisson)
220 class(external_potential_t), intent(inout) :: this
221 type(namespace_t), intent(in) :: namespace
222 class(mesh_t), intent(in) :: mesh
223 type(poisson_t), intent(in) :: poisson
224
225 real(real64) :: pot_re, pot_im, r, xx(this%space%dim)
226 real(real64), allocatable :: den(:), grx(:)
227 integer :: ip, err
228
230
231 select case (this%type)
232
233 case (external_pot_usdef)
234 assert(allocated(this%pot))
235
236 do ip = 1, mesh%np
237 call mesh_r(mesh, ip, r, coords = xx)
238 call parse_expression(pot_re, pot_im, this%space%dim, xx, r, m_zero, this%potential_formula)
239 this%pot(ip) = pot_re
240 end do
241
243 assert(allocated(this%pot))
244
245 call dio_function_input(trim(this%filename), namespace, this%space, mesh, this%pot, err)
246 if (err /= 0) then
247 write(message(1), '(a)') 'Error loading file '//trim(this%filename)//'.'
248 write(message(2), '(a,i4)') 'Error code returned = ', err
249 call messages_fatal(2, namespace=namespace)
250 end if
251
253 assert(allocated(this%pot))
254
255 safe_allocate(den(1:mesh%np))
256
257 do ip = 1, mesh%np
258 call mesh_r(mesh, ip, r, coords = xx)
259 call parse_expression(pot_re, pot_im, this%space%dim, xx, r, m_zero, this%potential_formula)
260 den(ip) = pot_re
261 end do
262
263 call dpoisson_solve(poisson, namespace, this%pot, den, all_nodes = .false.)
264
265 safe_deallocate_a(den)
266
268 assert(allocated(this%b_field))
269
270 ! Compute the vector potential from a uniform B field.
271 ! The sign is determined by the relation $\vec{B} = \nabla \times \vec{A}$.
272 ! This leads to $\vec{A} = -\frac{1}{2}\vec{r}\times\vec{B}$.
273 ! A factor 1/c is already added, to avoid adding it everytime we update the Hamiltonian
274 safe_allocate(grx(1:this%space%dim))
275
276 select case (this%space%dim)
277 case (2)
278 select case (this%gauge_2d)
279 case (0) ! linear_xy
280 if (this%space%periodic_dim == 1) then
281 message(1) = "For 2D system, 1D-periodic, StaticMagneticField can only be "
282 message(2) = "applied for StaticMagneticField2DGauge = linear_y."
283 call messages_fatal(2, namespace=namespace)
284 end if
285 !$omp parallel do private(grx)
286 do ip = 1, mesh%np
287 grx(1:this%space%dim) = mesh%x(ip, 1:this%space%dim)
288 this%a_static(ip, 1) = -m_half/p_c * grx(2) * this%b_field(3)
289 this%a_static(ip, 2) = m_half/p_c * grx(1) * this%b_field(3)
290 end do
291 case (1) ! linear y
292 !$omp parallel do private(grx)
293 do ip = 1, mesh%np
294 grx(1:this%space%dim) = mesh%x(ip, 1:this%space%dim)
295 this%a_static(ip, 1) = -m_one/p_c * grx(2) * this%b_field(3)
296 this%a_static(ip, 2) = m_zero
297 end do
298 end select
299 case (3)
300 !$omp parallel do private(grx)
301 do ip = 1, mesh%np
302 grx(1:this%space%dim) = mesh%x(ip, 1:this%space%dim)
303 this%a_static(ip, 1) = -m_half/p_c*(grx(2) * this%b_field(3) - grx(3) * this%b_field(2))
304 this%a_static(ip, 2) = -m_half/p_c*(grx(3) * this%b_field(1) - grx(1) * this%b_field(3))
305 this%a_static(ip, 3) = -m_half/p_c*(grx(1) * this%b_field(2) - grx(2) * this%b_field(1))
306 end do
307 case default
308 assert(.false.)
309 end select
310
311 safe_deallocate_a(grx)
314 assert(allocated(this%e_field))
315
316 if (this%space%periodic_dim < this%space%dim) then
317 ! Compute the scalar potential
318 !
319 ! Note that the -1 sign is missing. This is because we
320 ! consider the electrons with +1 charge. The electric field
321 ! however retains the sign because we also consider protons to
322 ! have +1 charge when calculating the force.
323 !
324 ! NTD: This comment is very confusing and prone to error
325 ! TODO: Fix this to have physically sound quantities and interactions
326 do ip = 1, mesh%np
327 this%pot(ip) = sum(mesh%x(ip, this%space%periodic_dim + 1:this%space%dim) &
328 * this%e_field(this%space%periodic_dim + 1:this%space%dim))
329 end do
330 ! The following is needed to make interpolations.
331 ! It is used by PCM.
332 this%v_ext(1:mesh%np) = this%pot(1:mesh%np)
333 do ip = mesh%np+1, mesh%np_part
334 this%v_ext(ip) = sum(mesh%x(ip, this%space%periodic_dim + 1:this%space%dim) &
335 * this%e_field(this%space%periodic_dim + 1:this%space%dim))
336 end do
337 end if
338
339 end select
340
342 end subroutine external_potential_calculate
343
344 subroutine load_external_potentials(external_potentials, namespace)
345 class(partner_list_t), intent(inout) :: external_potentials
346 type(namespace_t), intent(in) :: namespace
347
348 integer :: n_pot_block, row, read_data
349 type(block_t) :: blk
350 class(external_potential_t), pointer :: pot
351
352 integer :: dim, periodic_dim, idir
353
355
356 !%Variable StaticExternalPotentials
357 !%Type block
358 !%Section System
359 !%Description
360 !% An static external potential is a model potential added to the local potential of the Hamiltonian
361 !%
362 !% The format of this block is the following:
363 !% The first field defines the type of species (the valid options are detailed
364 !% below).
365 !%
366 !% Then a list of parameters follows. The parameters are specified
367 !% by a first field with the parameter name and the field that
368 !% follows with the value of the parameter. Some parameters are
369 !% specific to a certain species while others are accepted by all
370 !% species. These are <tt>mass</tt>, <tt>max_spacing</tt>, and <tt>min_radius</tt>.
371 !%
372 !% These are examples of possible species:
373 !%
374 !% <tt>%ExternalPotential
375 !% <br>&nbsp;&nbsp; potential_user_defined | potential_formula | "1/2*r^2"
376 !% <br>%</tt>
377 !%Option potential_from_file -202
378 !% The potential is read from a file. Accepted file formats, detected by extension: obf, ncdf and csv.
379 !%Option potential_user_defined -201
380 !% Species with user-defined potential. The potential for the
381 !% species is defined by the formula given by the <tt>potential_formula</tt>
382 !% parameter.
383 !%Option potential_charge_density -203
384 !% The potential for this species is created from the distribution
385 !% of charge given by the <tt>density_formula</tt> parameter.
386 !%Option file -10010
387 !% The path for the file that describes the species.
388 !%Option potential_formula -10012
389 !% Mathematical expression that defines the potential for <tt>species_user_defined</tt>. You can use
390 !% any of the <i>x</i>, <i>y</i>, <i>z</i> or <i>r</i> variables.
391 !%Option density_formula -10013
392 !% Mathematical expression that defines the charge density for <tt>species_charge_density</tt>. You can use
393 !% any of the <i>x</i>, <i>y</i>, <i>z</i> or <i>r</i> variables.
394 !%End
395
396 ! First, find out if there is a Species block.
397 n_pot_block = 0
398 if (parse_block(namespace, 'StaticExternalPotentials', blk) == 0) then
399 n_pot_block = parse_block_n(blk)
400
401 do row = 0, n_pot_block-1
402 !Create a potential
403 pot => external_potential_t(namespace)
404 !Parse the information from the block
405 call read_from_block(pot, namespace, blk, row, read_data)
406 assert(read_data > 0)
407 !Add this to the list
408 call external_potentials%add(pot)
409 end do
410 call parse_block_end(blk)
411 end if
412
413
414 !Here I am parsing the variables Dimensions et PeriodicDimensions because we do not have access
415 !to this information here.
416 !TODO: This needs to be removed and replaced by something better
417 call parse_variable(namespace, 'Dimensions', 3, dim)
418 call parse_variable(namespace, 'PeriodicDimensions', 0, periodic_dim)
419
420
421 !%Variable StaticMagneticField
422 !%Type block
423 !%Section Hamiltonian
424 !%Description
425 !% A static constant magnetic field may be added to the usual Hamiltonian,
426 !% by setting the block <tt>StaticMagneticField</tt>.
427 !% The three possible components of the block (which should only have one
428 !% line) are the three components of the magnetic field vector. Note that
429 !% if you are running the code in 1D mode, this will not work, and if you
430 !% are running the code in 2D mode the magnetic field will have to be in
431 !% the <i>z</i>-direction, so that the first two columns should be zero.
432 !% Possible in periodic system only in these cases: 2D system, 1D periodic,
433 !% with <tt>StaticMagneticField2DGauge = linear_y</tt>;
434 !% 3D system, 1D periodic, field is zero in <i>y</i>- and <i>z</i>-directions (given
435 !% currently implemented gauges).
436 !%
437 !% The magnetic field should always be entered in atomic units, regardless
438 !% of the <tt>Units</tt> variable. Note that we use the "Gaussian" system
439 !% meaning 1 au[B] = <math> 2.350517568\times 10^9</math> Gauss, which corresponds to
440 !% <math>2.3505175678\times 10^5</math> Tesla.
441 !%End
442 if (parse_block(namespace, 'StaticMagneticField', blk) == 0) then
443 !Create a potential
444 pot => external_potential_t(namespace)
446 pot%supported_interactions_as_partner = [lorentz_force]
447
448 !%Variable StaticMagneticField2DGauge
449 !%Type integer
450 !%Default linear_xy
451 !%Section Hamiltonian
452 !%Description
453 !% The gauge of the static vector potential <math>A</math> when a magnetic field
454 !% <math>B = \left( 0, 0, B_z \right)</math> is applied to a 2D-system.
455 !%Option linear_xy 0
456 !% Linear gauge with <math>A = \frac{1}{2c} \left( -y, x \right) B_z</math>. (Cannot be used for periodic systems.)
457 !%Option linear_y 1
458 !% Linear gauge with <math>A = \frac{1}{c} \left( -y, 0 \right) B_z</math>. Can be used for <tt>PeriodicDimensions = 1</tt>
459 !% but not <tt>PeriodicDimensions = 2</tt>.
460 !%End
461 call parse_variable(namespace, 'StaticMagneticField2DGauge', 0, pot%gauge_2d)
462 if (.not. varinfo_valid_option('StaticMagneticField2DGauge', pot%gauge_2d)) then
463 call messages_input_error(namespace, 'StaticMagneticField2DGauge')
464 end if
465
466 safe_allocate(pot%b_field(1:3))
467 do idir = 1, 3
468 call parse_block_float(blk, 0, idir - 1, pot%b_field(idir))
469 end do
470 select case (dim)
471 case (1)
472 call messages_input_error(namespace, 'StaticMagneticField')
473 case (2)
474 if (periodic_dim == 2) then
475 message(1) = "StaticMagneticField cannot be applied in a 2D, 2D-periodic system."
476 call messages_fatal(1, namespace=namespace)
477 end if
478 if (pot%b_field(1)**2 + pot%b_field(2)**2 > m_zero) then
479 call messages_input_error(namespace, 'StaticMagneticField')
480 end if
481 case (3)
482 ! Consider cross-product below: if grx(1:this%space%periodic_dim) is used, it is not ok.
483 ! Therefore, if idir is periodic, b_field for all other directions must be zero.
484 ! 1D-periodic: only Bx. 2D-periodic or 3D-periodic: not allowed. Other gauges could allow 2D-periodic case.
485 if (periodic_dim >= 2) then
486 message(1) = "In 3D, StaticMagneticField cannot be applied when the system is 2D- or 3D-periodic."
487 call messages_fatal(1, namespace=namespace)
488 else if (periodic_dim == 1 .and. any(abs(pot%b_field(2:3)) > m_zero)) then
489 message(1) = "In 3D, 1D-periodic, StaticMagneticField must be zero in the y- and z-directions."
490 call messages_fatal(1, namespace=namespace)
491 end if
492 end select
493 call parse_block_end(blk)
494
495 if (dim > 3) call messages_not_implemented('Magnetic field for dim > 3', namespace=namespace)
496
497 !Add this to the list
498 call external_potentials%add(pot)
499
500 !The corresponding A field on the mesh is computed in the routine external_potential_calculate
501
502 end if
503
504 !%Variable StaticElectricField
505 !%Type block
506 !%Default 0
507 !%Section Hamiltonian
508 !%Description
509 !% A static constant electric field may be added to the usual Hamiltonian,
510 !% by setting the block <tt>StaticElectricField</tt>.
511 !% The three possible components of the block (which should only have one
512 !% line) are the three components of the electric field vector.
513 !% It can be applied in a periodic direction of a large supercell via
514 !% the single-point Berry phase.
515 !%End
516 if (parse_block(namespace, 'StaticElectricField', blk) == 0) then
517 !Create a potential
518 pot => external_potential_t(namespace)
520 pot%supported_interactions_as_partner = [lorentz_force]
521
522 safe_allocate(pot%e_field(1:dim))
523 do idir = 1, dim
524 call parse_block_float(blk, 0, idir - 1, pot%e_field(idir), units_inp%energy / units_inp%length)
525
526 !Electron-specific checks (k-points) are done in the hamiltonian_elec.F90 file
527 if (idir <= periodic_dim .and. abs(pot%e_field(idir)) > m_epsilon) then
528 message(1) = "Applying StaticElectricField in a periodic direction is only accurate for large supercells."
529 call messages_warning(1, namespace=namespace)
530 end if
531 end do
532 call parse_block_end(blk)
533
534 !Add this to the list
535 call external_potentials%add(pot)
536
537 !The corresponding A field on the mesh is computed in the routine external_potential_calculate
538 end if
539
540
542 end subroutine load_external_potentials
543
544 ! ---------------------------------------------------------
545 subroutine read_from_block(pot, namespace, blk, row, read_data)
546 type(external_potential_t), intent(inout) :: pot
547 type(namespace_t), intent(in) :: namespace
548 type(block_t), intent(in) :: blk
549 integer, intent(in) :: row
550 integer, intent(out) :: read_data
551
552 integer :: ncols, icol, flag
553 type(iihash_t) :: read_parameters
554
555
556 push_sub(read_from_block)
557
558 ncols = parse_block_cols(blk, row)
559 read_data = 0
560
561 call parse_block_integer(blk, row, 0, pot%type)
562
563 ! To detect the old species block format, options are represented
564 ! as negative values. If we get a non-negative value we know we
565 ! are reading a mass.
566 if (pot%type >= 0) then
567 message(1) = 'Error in reading the ExternalPotentials block'
568 call messages_fatal(1, namespace=namespace)
569 end if
570
571 ! now we convert back to positive
572 pot%type = -pot%type
573
574 read_data = 1
575
576 if (pot%type /= external_pot_charge_density .and. pot%type /= external_pot_usdef .and. pot%type /= external_pot_from_file) then
577 call messages_input_error(namespace, 'ExternalPotentials', "Unknown type of external potential")
578 end if
579
580 call iihash_init(read_parameters)
581
582 icol = read_data
583 do
584 if (icol >= ncols) exit
585
586 call parse_block_integer(blk, row, icol, flag)
587
588 select case (flag)
589
590 case (option__staticexternalpotentials__file)
591 call check_duplication(option__staticexternalpotentials__file)
592 call parse_block_string(blk, row, icol + 1, pot%filename)
593
594 case (option__staticexternalpotentials__potential_formula)
595 call check_duplication(option__staticexternalpotentials__potential_formula)
596 call parse_block_string(blk, row, icol + 1, pot%potential_formula)
597 call conv_to_c_string(pot%potential_formula)
598
599 if (pot%type /= external_pot_usdef) then
600 call messages_input_error(namespace, 'ExternalPotentials', 'potential_formula can only be used with user_defined')
601 end if
602
603 case (option__staticexternalpotentials__density_formula)
604 call check_duplication(option__staticexternalpotentials__density_formula)
605 call parse_block_string(blk, row, icol + 1, pot%density_formula)
606 call conv_to_c_string(pot%density_formula)
607
608 if (pot%type /= external_pot_charge_density) then
609 call messages_input_error(namespace, 'ExternalPotentials', 'density_formula can only be used with charge_density')
610 end if
611
612 case default
613 call messages_input_error(namespace, 'ExternalPotentials', "Unknown parameter ")
614
615 end select
616
617 icol = icol + 2
618 end do
619 ! CHECK THAT WHAT WE PARSED MAKES SENSE
620
621
622 if (pot%type == external_pot_usdef .and. &
623 .not. parameter_defined(option__staticexternalpotentials__potential_formula)) then
624 call messages_input_error(namespace, 'ExternalPotentials', "The 'potential_formula' parameter is missing.")
625 end if
626
627 if (pot%type == external_pot_charge_density .and. &
628 .not. parameter_defined(option__staticexternalpotentials__density_formula)) then
629 call messages_input_error(namespace, 'ExternalPotentials', "The 'density_formula' parameter is missing.")
630 end if
631
632 if (pot%type == external_pot_from_file .and. .not. (parameter_defined(option__staticexternalpotentials__file))) then
633 call messages_input_error(namespace, 'ExternalPotentials', "The 'file' parameter is missing.")
634 end if
635
636 call iihash_end(read_parameters)
637
639
640 contains
641
642 logical function parameter_defined(param) result(defined)
643 integer(int64), intent(in) :: param
644
645 integer :: tmp
646
648
649 tmp = iihash_lookup(read_parameters, int(-param), defined)
650
652 end function parameter_defined
653
654 !------------------------------------------------------
655
656 subroutine check_duplication(param)
657 integer(int64), intent(in) :: param
658
660
661 if (parameter_defined(param)) then
662 call messages_input_error(namespace, 'ExternalPotentials', "Duplicated parameter in external potential.")
663 end if
664
665 call iihash_insert(read_parameters, int(-param), 1)
666
668 end subroutine check_duplication
669
670 end subroutine read_from_block
671 ! ---------------------------------------------------------
672
673
675
676!! Local Variables:
677!! mode: f90
678!! coding: utf-8
679!! End:
subroutine check_duplication(param)
logical function parameter_defined(param)
integer, parameter, public external_pot_from_file
potential, defined in a file
subroutine, public load_external_potentials(external_potentials, namespace)
class(external_potential_t) function, pointer external_potential_init(namespace)
subroutine read_from_block(pot, namespace, blk, row, read_data)
integer, parameter, public external_pot_charge_density
user-defined function for charge density
subroutine external_potential_finalize(this)
subroutine external_potential_calculate(this, namespace, mesh, poisson)
subroutine external_potential_init_interaction_as_partner(partner, interaction)
subroutine external_potential_deallocate(this)
subroutine external_potential_copy_quantities_to_interaction(partner, interaction)
integer, parameter, public external_pot_static_efield
Static electric field.
integer, parameter, public external_pot_static_bfield
Static magnetic field.
subroutine external_potential_allocate(this, mesh)
real(real64), parameter, public m_zero
Definition: global.F90:188
real(real64), parameter, public m_epsilon
Definition: global.F90:204
real(real64), parameter, public m_half
Definition: global.F90:194
real(real64), parameter, public p_c
Electron gyromagnetic ratio, see Phys. Rev. Lett. 130, 071801 (2023)
Definition: global.F90:224
real(real64), parameter, public m_one
Definition: global.F90:189
This module implements a simple hash table for non-negative integer keys and integer values.
Definition: iihash.F90:125
subroutine, public iihash_end(h)
Free a hash table.
Definition: iihash.F90:184
subroutine, public iihash_insert(h, key, val)
Insert a (key, val) pair into the hash table h.
Definition: iihash.F90:206
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:231
subroutine, public iihash_init(h)
Initialize a hash table h.
Definition: iihash.F90:161
integer, parameter, public lorentz_force
This module defines classes and functions for interaction partners.
subroutine, public dio_function_input(filename, namespace, space, mesh, ff, ierr, map)
Reads a mesh function from file filename, and puts it into ff. If the map argument is passed,...
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:118
pure subroutine, public mesh_r(mesh, ip, rr, origin, coords)
return the distance to the origin for a given grid point
Definition: mesh.F90:336
subroutine, public messages_not_implemented(feature, namespace)
Definition: messages.F90:1113
subroutine, public messages_warning(no_lines, all_nodes, namespace)
Definition: messages.F90:537
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:160
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
Definition: messages.F90:414
subroutine, public messages_input_error(namespace, var, details, row, column)
Definition: messages.F90:713
integer function, public parse_block(namespace, name, blk, check_varinfo_)
Definition: parser.F90:618
subroutine, public dpoisson_solve(this, namespace, pot, rho, all_nodes, kernel, reset)
Calculates the Poisson equation. Given the density returns the corresponding potential.
Definition: poisson.F90:892
This module defines the quantity_t class and the IDs for quantities, which can be exposed by a system...
Definition: quantity.F90:138
subroutine, public conv_to_c_string(str)
converts to c string
Definition: string.F90:252
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
Definition: unit.F90:132
This module defines the unit system, used for input and output.
type(unit_system_t), public units_inp
the units systems for reading and writing
abstract class for general interaction partners
surrogate interaction class to avoid circular dependencies between modules.
This class implements the iteration counter used by the multisystem algorithms. As any iteration coun...
Lorenz force between a systems of particles and an electromagnetic field.
Describes mesh distribution to nodes.
Definition: mesh.F90:186
Systems (system_t) can expose quantities that can be used to calculate interactions with other system...
Definition: quantity.F90:171
int true(void)