65 type(namespace_t),
intent(in) :: namespace
66 class(space_t),
intent(in) :: space
67 class(mesh_t),
intent(inout) :: mesh
68 type(derivatives_t),
intent(in) :: der
69 type(states_mxll_t),
intent(inout) :: st
70 type(bc_mxll_t),
intent(inout) :: bc
71 complex(real64),
contiguous,
intent(inout) :: user_def_rs_state(:,:)
74 integer :: il, nlines, idim, ncols, ip, state_from, ierr, maxwell_field
75 real(real64) :: xx(space%dim), rr, e_value, dummy, b_value
76 real(real64),
allocatable :: e_field(:), b_field(:)
77 real(real64),
allocatable :: total_efield(:,:), total_bfield(:,:)
78 complex(real64),
allocatable :: rs_state_add(:), rs_state(:,:)
79 character(len=150),
pointer :: filename_e_field, filename_b_field
82 integer,
parameter :: &
83 STATE_FROM_FORMULA = 1, &
84 state_from_incident_waves = 2, &
85 state_from_file = -10010
143 if (
parse_block(namespace,
'UserDefinedInitialMaxwellStates', blk) == 0)
then
145 safe_allocate(rs_state_add(1:mesh%np_part))
146 safe_allocate(rs_state(1:mesh%np, 1:3))
149 user_def_rs_state(:,:) =
m_zero
154 write(
message(1),
'(a,i5)')
'Maxwell electromagnetic fields are added.'
159 safe_allocate(total_efield(1:mesh%np, 1:3))
160 safe_allocate(total_bfield(1:mesh%np, 1:3))
168 if (ncols /= 1 .and. ncols /= 4)
then
169 message(1) =
'Each line in the UserDefinedMaxwellStates block must have'
170 message(2) =
'one or four columns.'
179 write(cdim,
'(I1)')idim
184 select case (state_from)
186 case (state_from_formula)
192 if (maxwell_field == option__userdefinedinitialmaxwellstates__electric_field)
then
194 call messages_write(
" E-field in dimension "//trim(cdim)//
" : "//trim(st%user_def_e_field(idim)), fmt=
'(a,i1,2a)')
196 else if (maxwell_field == option__userdefinedinitialmaxwellstates__magnetic_field)
then
198 call messages_write(
" B-field in dimension "//trim(cdim)//
" : "//trim(st%user_def_b_field(idim)), fmt=
'(a,i1,2a)')
202 if (maxwell_field == option__userdefinedinitialmaxwellstates__electric_field)
then
210 st%user_def_e_field(idim))
211 total_efield(ip, idim) = total_efield(ip, idim) + e_value
214 else if (maxwell_field == option__userdefinedinitialmaxwellstates__magnetic_field)
then
221 st%user_def_b_field(idim))
222 total_bfield(ip, idim) = total_bfield(ip, idim) + b_value
227 case (state_from_file)
234 safe_allocate(e_field(1:mesh%np))
235 safe_allocate(b_field(1:mesh%np))
238 if (maxwell_field == option__userdefinedinitialmaxwellstates__electric_field)
then
240 call messages_write(
" E-field in dimension "//trim(cdim)//
" : "//trim(filename_e_field), fmt=
'(a,i1,2a)')
243 message(1) =
'Could not read the file!'
244 write(
message(2),
'(a,i1)')
'Error code: ', ierr
247 else if (maxwell_field == option__userdefinedinitialmaxwellstates__magnetic_field)
then
249 call messages_write(
" B-field in dimension "//trim(cdim)//
" : "//trim(filename_b_field), fmt=
'(a,i1,2a)')
252 message(1) =
'Could not read the file!'
253 write(
message(2),
'(a,i1)')
'Error code: ', ierr
258 call build_rs_vector(e_field(:), b_field(:), st%rs_sign, rs_state_add(:), &
259 st%ep(ip), st%mu(ip))
261 safe_deallocate_a(e_field)
262 safe_deallocate_a(b_field)
264 call lalg_axpy(mesh%np,
m_one, rs_state_add, user_def_rs_state(:,idim))
266 case (state_from_incident_waves)
272 message(1) =
'Wrong entry in UserDefinedMaxwellStates, column 2.'
273 message(2) =
'You may state "formula", "file" or "use_incident_waves" here.'
279 if (state_from == state_from_formula)
then
281 call build_rs_state(total_efield, total_bfield, st%rs_sign, rs_state, mesh, st%ep, st%mu)
286 safe_deallocate_a(total_efield)
287 safe_deallocate_a(total_bfield)
289 safe_deallocate_a(rs_state)
290 safe_deallocate_a(rs_state_add)
305 subroutine states_mxll_dump(restart, st, space, mesh, zff, zff_dim, ierr, iter, st_start_writing, verbose)
308 class(
space_t),
intent(in) :: space
309 class(
mesh_t),
intent(in) :: mesh
310 complex(real64),
intent(in) :: zff(:,:)
311 integer,
intent(in) :: zff_dim
312 integer,
intent(out) :: ierr
313 integer,
optional,
intent(in) :: iter
314 integer,
optional,
intent(in) :: st_start_writing
315 logical,
optional,
intent(in) :: verbose
317 integer :: iunit_wfns, iunit_states
318 integer :: err, err2(2), ist, idim, itot
320 character(len=MAX_PATH_LEN) :: filename
321 character(len=300) :: lines(3)
322 logical :: should_write, verbose_
336 message(1) =
"Info: Writing Maxwell states."
345 write(lines(1),*) zff_dim
347 if (err /= 0) ierr = ierr + 1
351 lines(1) =
'# #dim filename'
352 lines(2) =
'%RS States'
354 if (err /= 0) ierr = ierr + 2
366 write(filename,
'(i10.10)') itot
368 write(lines(1),
'(i8,3a)') idim,
' | "', trim(filename),
'"'
370 if (err /= 0) err2(1) = err2(1) + 1
372 should_write = st%st_start <= ist .and. ist <= st%st_end
373 if (should_write .and.
present(st_start_writing))
then
374 if (ist < st_start_writing) should_write = .false.
377 if (should_write)
then
379 if (err /= 0) err2(2) = err2(2) + 1
383 if (err2(1) /= 0) ierr = ierr + 8
384 if (err2(2) /= 0) ierr = ierr + 16
388 if (err /= 0) ierr = ierr + 64
389 if (
present(iter))
then
390 write(lines(1),
'(a,i7)')
'Iter = ', iter
392 if (err /= 0) ierr = ierr + 128
398 message(1) =
"Info: Finished writing Maxwell states."
411 subroutine states_mxll_load(restart, st, mesh, namespace, space, zff, zff_dim, ierr, iter, lowest_missing, label, verbose)
414 class(
mesh_t),
intent(in) :: mesh
416 class(
space_t),
intent(in) :: space
417 complex(real64),
contiguous,
intent(inout) :: zff(:,:)
418 integer,
intent(in) :: zff_dim
419 integer,
intent(out) :: ierr
420 integer,
optional,
intent(out) :: iter
421 integer,
optional,
intent(out) :: lowest_missing(:)
422 character(len=*),
optional,
intent(in) :: label
423 logical,
optional,
intent(in) :: verbose
425 integer :: states_file, wfns_file, err, ist, idim, dim, mx_st_start, mx_st_end
426 integer :: idone, iread, ntodo
427 character(len=12) :: filename
428 character(len=1) :: char
429 logical,
allocatable :: filled(:, :)
430 character(len=256) :: lines(3), label_
433 character(len=256),
allocatable :: restart_file(:, :)
434 logical,
allocatable :: restart_file_present(:, :)
443 if (
present(lowest_missing)) lowest_missing = 1
444 if (
present(iter)) iter = 0
456 if (
present(label))
then
460 message(1) =
'Info: Reading Maxwell states'
461 if (len(trim(label_)) > 0)
then
472 read(lines(1), *) idim
492 safe_allocate(restart_file(1:zff_dim, st%st_start:st%st_end))
493 safe_allocate(restart_file_present(1:zff_dim,st%st_start:st%st_end))
494 restart_file_present = .false.
502 read(lines(1),
'(a)') char
503 if (char ==
'%')
then
507 read(lines(1), *) idim, char, filename
511 if (ist >= st%st_start .and. ist <= st%st_end)
then
512 restart_file(idim, ist) = trim(filename)
513 restart_file_present(idim, ist) = .
true.
517 if (
present(iter))
then
522 read(lines(1), *) filename, filename, iter
532 mx_st_start=st%st_start
534 safe_allocate(filled(1:zff_dim,mx_st_start:mx_st_end))
537 if (
present(lowest_missing)) lowest_missing = st%nst + 1
542 ntodo = st%lnst*zff_dim
548 if (.not. restart_file_present(idim, ist))
then
549 if (
present(lowest_missing))
then
550 lowest_missing(idim) = min(lowest_missing(idim), ist)
559 filled(idim, ist) = .
true.
561 else if (
present(lowest_missing))
then
562 lowest_missing(idim) = min(lowest_missing(idim), ist)
572 safe_deallocate_a(restart_file)
573 safe_deallocate_a(restart_file_present)
574 safe_deallocate_a(filled)
580 if (ierr == 0 .and. iread /= st%nst * zff_dim)
then
589 write(
message(1),
'(a,i6,a,i6,a)')
'Only ', iread,
' files out of ', &
590 st%nst * zff_dim,
' could be read.'
595 message(1) =
'Info: Maxwell states reading done.'
constant times a vector plus a vector
This module implements batches of mesh functions.
This module implements common operations on batches of mesh functions.
This module handles the calculation mode.
integer, parameter, public p_strategy_max
integer, parameter, public p_strategy_domains
parallelization in domains
This module calculates the derivatives (gradients, Laplacians, etc.) of a function.
real(real64), parameter, public m_zero
real(real64), parameter, public m_one
This module implements the underlying real-space grid.
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 various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
subroutine, public messages_print_with_emphasis(msg, iunit, namespace)
character(len=512), private msg
subroutine, public messages_variable_is_block(namespace, name)
subroutine, public print_date(str)
subroutine, public messages_new_line()
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
logical function mpi_grp_is_root(grp)
Is the current MPI process of grpcomm, root.
type(mpi_grp_t), public mpi_world
This module handles the communicators for the various parallelization strategies.
integer function, public parse_block(namespace, name, blk, check_varinfo_)
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
subroutine, public plane_waves_in_box_calculation(bc, time, space, mesh, der, st, rs_state)
subroutine, public restart_read(restart, iunit, lines, nlines, ierr)
subroutine, public restart_close(restart, iunit)
Close a file previously opened with restart_open.
subroutine, public restart_write(restart, iunit, lines, nlines, ierr)
logical pure function, public restart_skip(restart)
Returns true if the restart information should neither be read nor written. This might happen because...
integer function, public restart_open(restart, filename, status, position, silent)
Open file "filename" found inside the current restart directory. Depending on the type of restart,...
This module handles spin dimensions of the states and the k-point distribution.
This module handles reading and writing restart information for the states_elec_t.
subroutine, public build_rs_vector(e_vector, b_vector, rs_sign, rs_vector, ep_element, mu_element)
subroutine, public build_rs_state(e_field, b_field, rs_sign, rs_state, mesh, ep_field, mu_field, np)
subroutine, public states_mxll_load(restart, st, mesh, namespace, space, zff, zff_dim, ierr, iter, lowest_missing, label, verbose)
subroutine, public states_mxll_dump(restart, st, space, mesh, zff, zff_dim, ierr, iter, st_start_writing, verbose)
subroutine, public states_mxll_read_user_def(namespace, space, mesh, der, st, bc, user_def_rs_state)
subroutine, public conv_to_c_string(str)
converts to c string
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
This module defines the unit system, used for input and output.
Describes mesh distribution to nodes.