Octopus
states_mxll_restart.F90
Go to the documentation of this file.
1!! Copyright (C) 2019 R. Jestaedt, F. Bonafe, H. Appel, A. Rubio
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
19#include "global.h"
20
22 use batch_oct_m
25 use debug_oct_m
27 use global_oct_m
28 use grid_oct_m
31 use loct_oct_m
32 use mesh_oct_m
35 use mpi_oct_m
39 use parser_oct_m
43 use space_oct_m
47 use string_oct_m
48 use types_oct_m
50 use unit_oct_m
51
52 implicit none
53
54
55 private
56 public :: &
60
61contains
62
63 !----------------------------------------------------------
64 subroutine states_mxll_read_user_def(namespace, space, mesh, der, st, bc, user_def_rs_state)
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(:,:)
72
73 type(block_t) :: blk
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
80 character(1) :: cdim
81
82 integer, parameter :: &
83 STATE_FROM_FORMULA = 1, &
84 state_from_incident_waves = 2, &
85 state_from_file = -10010
86
88
89 call profiling_in('STATES_MXLL_READ_USER_DEF')
90
91 !%Variable UserDefinedInitialMaxwellStates
92 !%Type block
93 !%Section Maxwell
94 !%Description
95 !% The initial electromagnetic fields can be set by the user
96 !% with the <tt>UserDefinedMaxwellStates</tt> block variable.
97 !% The electromagnetic fields have to fulfill the
98 !% Maxwells equations in vacuum.
99 !%
100 !% Example:
101 !%
102 !% <tt>%UserDefinedMaxwellStates
103 !% <br>&nbsp;&nbsp; formula | 2 | "magnetic_field" | "-1/P_c * sin(x)"
104 !% <br>&nbsp;&nbsp; formula | 3 | "electric_field" | " sin(x) "
105 !% <br>%</tt>
106 !%
107 !% The second column specifies the component of the dimension of
108 !% the electric field and magnetic field. The first column
109 !% indicates that column four should be interpreted
110 !% as a formula for the corresponding state. P_c is the
111 !% speed of light constant.
112 !%
113 !% Alternatively, if column one states <tt>file</tt> the
114 !% electric field and magnetic field will be read from
115 !% the files given in column four.
116 !%
117 !% <tt>%UserDefinedMaxwellStates
118 !% <br>&nbsp;&nbsp; file | 3 | electric_field | "/path/to/file_electric_field_of_dimension_3"
119 !% <br>&nbsp;&nbsp; file | 2 | magnetic_field | "/path/to/file_magnetic_field_of_dimension_2"
120 !% <br>%</tt>
121 !%
122 !% The third option to define the initial state inside the box is to extend
123 !% the plane waves used as incident waves in the <tt>MaxwellIncidentWaves</tt> block,
124 !% as follows:
125 !%
126 !% <tt>%UserDefinedMaxwellStates
127 !% <br>&nbsp;&nbsp; use_incident_waves
128 !% <br>%</tt>
129 !%
130 !%Option file -10010
131 !% Read initial orbital from file.
132 !% Accepted file formats: obf, ncdf and csv.
133 !%Option formula 1
134 !% Calculate initial orbital by given analytic expression.
135 !%Option use_incident_waves 2
136 !% Extend the plane waves given in the <tt>MaxwellIncidentWaves</tt> block inside the box.
137 !%Option electric_field 1
138 !% This row defines the electric field component of the corresponding dimension
139 !%Option magnetic_field 2
140 !% This row defines the magnetic field component of the corresponding dimension
141 !%End
142
143 if (parse_block(namespace, 'UserDefinedInitialMaxwellStates', blk) == 0) then
144
145 safe_allocate(rs_state_add(1:mesh%np_part))
146 safe_allocate(rs_state(1:mesh%np, 1:3))
147
148 ! Set electromagnetic field equal to zero in the whole simulation box.
149 user_def_rs_state(:,:) = m_zero
150
151 ! find out how many lines (i.e. states) the block has
152 nlines = parse_block_n(blk)
153
154 write(message(1), '(a,i5)') 'Maxwell electromagnetic fields are added.'
155 write(message(2), '(a,i5)') ''
156 call messages_info(2, namespace=namespace)
158
159 safe_allocate(total_efield(1:mesh%np, 1:3))
160 safe_allocate(total_bfield(1:mesh%np, 1:3))
161 total_efield = m_zero
162 total_bfield = m_zero
163
164 ! read all lines
165 do il = 1, nlines
166 ! Check that number of columns is one or four.
167 ncols = parse_block_cols(blk, il - 1)
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.'
171 call messages_fatal(2, namespace=namespace)
172 end if
173
174 ! Calculate from expression, use or read from file?
175 call parse_block_integer(blk, il - 1, 0, state_from)
176
177 if (ncols > 1) then
178 call parse_block_integer(blk, il - 1, 1, idim)
179 write(cdim,'(I1)')idim
180 end if
181
182 rs_state_add(:) = m_zero
183
184 select case (state_from)
185
186 case (state_from_formula)
187
188 assert(ncols == 4)
189
190 ! parse formula string
191 call parse_block_integer(blk, il - 1, 2, maxwell_field)
192 if (maxwell_field == option__userdefinedinitialmaxwellstates__electric_field) then
193 call parse_block_string( blk, il - 1, 3, st%user_def_e_field(idim))
194 call messages_write(" E-field in dimension "//trim(cdim)//" : "//trim(st%user_def_e_field(idim)), fmt='(a,i1,2a)')
195 call conv_to_c_string(st%user_def_e_field(idim))
196 else if (maxwell_field == option__userdefinedinitialmaxwellstates__magnetic_field) then
197 call parse_block_string( blk, il - 1, 3, st%user_def_b_field(idim))
198 call messages_write(" B-field in dimension "//trim(cdim)//" : "//trim(st%user_def_b_field(idim)), fmt='(a,i1,2a)')
199 call conv_to_c_string(st%user_def_b_field(idim))
200 end if
201
202 if (maxwell_field == option__userdefinedinitialmaxwellstates__electric_field) then
203
204 ! fill Maxwell states with user-defined formulas
205 do ip = 1, mesh%np
206 xx = mesh%x(ip, :)
207 rr = norm2(xx)
208 ! parse user-defined expressions
209 call parse_expression(e_value, dummy, st%dim, xx, rr, m_zero, &
210 st%user_def_e_field(idim))
211 total_efield(ip, idim) = total_efield(ip, idim) + e_value
212 end do
213
214 else if (maxwell_field == option__userdefinedinitialmaxwellstates__magnetic_field) then
215 ! fill Maxwell states with user-defined formulas
216 do ip = 1, mesh%np
217 xx = mesh%x(ip, :)
218 rr = norm2(xx)
219
220 call parse_expression(b_value, dummy, st%dim, xx, rr, m_zero, &
221 st%user_def_b_field(idim))
222 total_bfield(ip, idim) = total_bfield(ip, idim) + b_value
223 end do
224
225 end if
226
227 case (state_from_file)
228
229 assert(ncols == 4)
230 ! The input format can be coded in column four now. As it is
231 ! not used now, we just say "file".
232 ! Read the filename.
233 call parse_block_integer(blk, il - 1, 2, maxwell_field)
234 safe_allocate(e_field(1:mesh%np))
235 safe_allocate(b_field(1:mesh%np))
236 e_field = m_zero
237 b_field = m_zero
238 if (maxwell_field == option__userdefinedinitialmaxwellstates__electric_field) then
239 call parse_block_string(blk, il - 1, 3, filename_e_field)
240 call messages_write(" E-field in dimension "//trim(cdim)//" : "//trim(filename_e_field), fmt='(a,i1,2a)')
241 call dio_function_input(filename_e_field, namespace, space, mesh, e_field(:), ierr)
242 if (ierr > 0) then
243 message(1) = 'Could not read the file!'
244 write(message(2),'(a,i1)') 'Error code: ', ierr
245 call messages_fatal(2, namespace=namespace)
246 end if
247 else if (maxwell_field == option__userdefinedinitialmaxwellstates__magnetic_field) then
248 call parse_block_string(blk, il - 1, 3, filename_b_field)
249 call messages_write(" B-field in dimension "//trim(cdim)//" : "//trim(filename_b_field), fmt='(a,i1,2a)')
250 call dio_function_input(filename_b_field, namespace, space, mesh, b_field(:), ierr)
251 if (ierr > 0) then
252 message(1) = 'Could not read the file!'
253 write(message(2),'(a,i1)') 'Error code: ', ierr
254 call messages_fatal(2, namespace=namespace)
255 end if
256 end if
257 ! fill state
258 call build_rs_vector(e_field(:), b_field(:), st%rs_sign, rs_state_add(:), &
259 st%ep(ip), st%mu(ip))
260
261 safe_deallocate_a(e_field)
262 safe_deallocate_a(b_field)
263
264 call lalg_axpy(mesh%np, m_one, rs_state_add, user_def_rs_state(:,idim))
265
266 case (state_from_incident_waves)
267
268 assert(ncols == 1)
269 call plane_waves_in_box_calculation(bc, m_zero, space, mesh, der, st, user_def_rs_state)
270
271 case default
272 message(1) = 'Wrong entry in UserDefinedMaxwellStates, column 2.'
273 message(2) = 'You may state "formula", "file" or "use_incident_waves" here.'
274 call messages_fatal(2, namespace=namespace)
275 end select
276
277 end do
278
279 if (state_from == state_from_formula) then
280 ! fill state
281 call build_rs_state(total_efield, total_bfield, st%rs_sign, rs_state, mesh, st%ep, st%mu)
282 call lalg_axpy(mesh%np, 3, m_one, rs_state, user_def_rs_state)
283 end if
284
285
286 safe_deallocate_a(total_efield)
287 safe_deallocate_a(total_bfield)
288
289 safe_deallocate_a(rs_state)
290 safe_deallocate_a(rs_state_add)
291 call parse_block_end(blk)
292 !call messages_print_with_emphasis(namespace=namespace)
293
294 else
295 call messages_variable_is_block(namespace, 'UserDefineInitialdStates')
296 end if
297
298 call profiling_out('STATES_MXLL_READ_USER_DEF')
299
301 end subroutine states_mxll_read_user_def
302
303
304 !----------------------------------------------------------
305 subroutine states_mxll_dump(restart, st, space, mesh, zff, zff_dim, ierr, iter, st_start_writing, verbose)
306 type(restart_t), intent(in) :: restart
307 type(states_mxll_t), intent(in) :: st
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
316
317 integer :: iunit_wfns, iunit_states
318 integer :: err, err2(2), ist, idim, itot
319 integer :: root(1:p_strategy_max)
320 character(len=MAX_PATH_LEN) :: filename
321 character(len=300) :: lines(3)
322 logical :: should_write, verbose_
323
324 push_sub(states_mxll_dump)
325
326 verbose_ = optional_default(verbose, .true.)
327
328 ierr = 0
329
330 if (restart_skip(restart)) then
331 pop_sub(states_mxll_dump)
332 return
333 end if
334
335 if (verbose_) then
336 message(1) = "Info: Writing Maxwell states."
337 call print_date(trim(message(1))//' ')
338 end if
339
340 call profiling_in("MAXWELL_RESTART_WRITE")
341
343
344 iunit_states = restart_open(restart, 'maxwell_states')
345 write(lines(1),*) zff_dim
346 call restart_write(restart, iunit_states, lines, 1, err)
347 if (err /= 0) ierr = ierr + 1
348 call restart_close(restart, iunit_states)
349
350 iunit_wfns = restart_open(restart, 'wfns')
351 lines(1) = '# #dim filename'
352 lines(2) = '%RS States'
353 call restart_write(restart, iunit_wfns, lines, 2, err)
354 if (err /= 0) ierr = ierr + 2
355
356
357 itot = 1
358 root = 0
359 err2 = 0
360 ist = 1
361
362 do idim = 1, zff_dim
363 itot = itot + 1
364
365 root(p_strategy_domains) = mod(itot - 1, mesh%mpi_grp%size)
366 write(filename,'(i10.10)') itot
367
368 write(lines(1), '(i8,3a)') idim, ' | "', trim(filename), '"'
369 call restart_write(restart, iunit_wfns, lines, 1, err)
370 if (err /= 0) err2(1) = err2(1) + 1
371
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.
375 end if
376
377 if (should_write) then
378 call restart_write_mesh_function(restart, space, filename, mesh, zff(:,idim), err, root = root)
379 if (err /= 0) err2(2) = err2(2) + 1
380 end if
381 end do ! zff_dim
382
383 if (err2(1) /= 0) ierr = ierr + 8
384 if (err2(2) /= 0) ierr = ierr + 16
385
386 lines(1) = '%'
387 call restart_write(restart, iunit_wfns, lines, 1, err)
388 if (err /= 0) ierr = ierr + 64
389 if (present(iter)) then
390 write(lines(1),'(a,i7)') 'Iter = ', iter
391 call restart_write(restart, iunit_wfns, lines, 1, err)
392 if (err /= 0) ierr = ierr + 128
393 end if
394
395 call restart_close(restart, iunit_wfns)
396
397 if (verbose_) then
398 message(1) = "Info: Finished writing Maxwell states."
399 call print_date(trim(message(1))//' ')
400 end if
401
403
404 call profiling_out("MAXWELL_RESTART_WRITE")
405 pop_sub(states_mxll_dump)
406 return
407
408 end subroutine states_mxll_dump
409
410 !----------------------------------------------------------
411 subroutine states_mxll_load(restart, st, mesh, namespace, space, zff, zff_dim, ierr, iter, lowest_missing, label, verbose)
412 type(restart_t), intent(in) :: restart
413 type(states_mxll_t), intent(inout) :: st
414 class(mesh_t), intent(in) :: mesh
415 type(namespace_t), intent(in) :: namespace
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
424
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_
431
432 logical :: verbose_
433 character(len=256), allocatable :: restart_file(:, :)
434 logical, allocatable :: restart_file_present(:, :)
435
436 push_sub(states_mxll_load)
437
438 ierr = 0
439 dim = zff_dim
440 ist = 1
441
442 ! make sure these intent(out)`s are initialized no matter what
443 if (present(lowest_missing)) lowest_missing = 1
444 if (present(iter)) iter = 0
445
446 if (restart_skip(restart)) then
447 ierr = -1
448 pop_sub(states_mxll_load)
449 return
450 end if
451
452 call profiling_in("RESTART_READ")
453
454 verbose_ = optional_default(verbose, .true.)
455
456 if (present(label)) then
457 label_ = trim(label)
458 end if
459
460 message(1) = 'Info: Reading Maxwell states'
461 if (len(trim(label_)) > 0) then
462 message(1) = trim(message(1)) // trim(label_)
463 end if
464 message(1) = trim(message(1)) // "."
465 if (verbose_) call print_date(trim(message(1))//' ')
466
467 states_file = restart_open(restart, 'maxwell_states')
468 call restart_read(restart, states_file, lines, 1, err)
469 if (err /= 0) then
470 ierr = ierr - 2
471 else
472 read(lines(1), *) idim
473 end if
474 call restart_close(restart, states_file)
475
476 ! open files to read
477 wfns_file = restart_open(restart, 'wfns')
478 call restart_read(restart, wfns_file, lines, 2, err)
479 if (err /= 0) then
480 ierr = ierr - 2**5
481 end if
482
483 ! If any error occured up to this point then it is not worth continuing,
484 ! as there something fundamentally wrong with the restart files
485 if (ierr /= 0) then
486 call restart_close(restart, wfns_file)
487 call profiling_out("RESTART_READ")
488 pop_sub(states_mxll_load)
489 return
490 end if
491
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.
495
496 ! Next we read the list of states from the files.
497 ! Errors in reading the information of a specific state from the files are ignored
498 ! at this point, because later we will skip reading the wavefunction of that state.
499 do
500 call restart_read(restart, wfns_file, lines, 1, err)
501 if (err == 0) then
502 read(lines(1), '(a)') char
503 if (char == '%') then
504 !We reached the end of the file
505 exit
506 else
507 read(lines(1), *) idim, char, filename
508 end if
509 end if
510
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.
514 end if
515
516 end do
517 if (present(iter)) then
518 call restart_read(restart, wfns_file, lines, 1, err)
519 if (err /= 0) then
520 ierr = ierr - 2**8
521 else
522 read(lines(1), *) filename, filename, iter
523 end if
524 end if
525
526 call restart_close(restart, wfns_file)
527
528 ! Now we read the wavefunctions. At this point we need to have all the information from the
529 ! states, and wfns files in order to avoid serialisation of reads, as restart_read
530 ! works as a barrier.
531
532 mx_st_start=st%st_start
533 mx_st_end=st%st_end
534 safe_allocate(filled(1:zff_dim,mx_st_start:mx_st_end))
535 filled = .false.
536
537 if (present(lowest_missing)) lowest_missing = st%nst + 1
538
539 iread = 0
540 if (mpi_grp_is_root(mpi_world) .and. verbose_) then
541 idone = 0
542 ntodo = st%lnst*zff_dim
543 call loct_progress_bar(-1, ntodo)
544 end if
545
546 ist = 1
547 do idim = 1, 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)
551 end if
552 cycle
553 end if
554
555 call restart_read_mesh_function(restart, space, restart_file(idim, ist), mesh, &
556 zff(:,idim), err)
557
558 if (err == 0) then
559 filled(idim, ist) = .true.
560 iread = iread + 1
561 else if (present(lowest_missing)) then
562 lowest_missing(idim) = min(lowest_missing(idim), ist)
563 end if
564
565 if (mpi_grp_is_root(mpi_world) .and. verbose_) then
566 idone = idone + 1
567 call loct_progress_bar(idone, ntodo)
568 end if
569
570 end do
571
572 safe_deallocate_a(restart_file)
573 safe_deallocate_a(restart_file_present)
574 safe_deallocate_a(filled)
575
576 if (mpi_grp_is_root(mpi_world) .and. verbose_) then
577 call messages_new_line()
578 end if
579
580 if (ierr == 0 .and. iread /= st%nst * zff_dim) then
581 if (iread > 0) then
582 ierr = iread
583 else
584 ierr = -1
585 end if
586 ! otherwise ierr = 0 would mean either all was read correctly, or nothing at all was read!
587
588 call messages_print_with_emphasis(msg='Reading Maxwell states.', namespace=namespace)
589 write(message(1),'(a,i6,a,i6,a)') 'Only ', iread,' files out of ', &
590 st%nst * zff_dim, ' could be read.'
591 call messages_info(1, namespace=namespace)
592 call messages_print_with_emphasis(namespace=namespace)
593 end if
594
595 message(1) = 'Info: Maxwell states reading done.'
596 if (verbose_) call print_date(trim(message(1))//' ')
597
598 call profiling_out("RESTART_READ")
599 pop_sub(states_mxll_load)
600
601 end subroutine states_mxll_load
602
604
605
606!! Local Variables:
607!! mode: f90
608!! coding: utf-8
609!! End:
constant times a vector plus a vector
Definition: lalg_basic.F90:171
This module implements batches of mesh functions.
Definition: batch.F90:133
This module implements common operations on batches of mesh functions.
Definition: batch_ops.F90:116
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
Definition: global.F90:188
real(real64), parameter, public m_one
Definition: global.F90:189
This module implements the underlying real-space grid.
Definition: grid.F90:117
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.
Definition: mesh.F90:118
subroutine, public messages_print_with_emphasis(msg, iunit, namespace)
Definition: messages.F90:903
character(len=512), private msg
Definition: messages.F90:165
subroutine, public messages_variable_is_block(namespace, name)
Definition: messages.F90:1052
subroutine, public print_date(str)
Definition: messages.F90:988
subroutine, public messages_new_line()
Definition: messages.F90:1117
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_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:599
logical function mpi_grp_is_root(grp)
Is the current MPI process of grpcomm, root.
Definition: mpi.F90:370
type(mpi_grp_t), public mpi_world
Definition: mpi.F90:270
This module handles the communicators for the various parallelization strategies.
Definition: multicomm.F90:145
integer function, public parse_block(namespace, name, blk, check_varinfo_)
Definition: parser.F90:618
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
Definition: profiling.F90:623
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
Definition: profiling.F90:552
subroutine, public plane_waves_in_box_calculation(bc, time, space, mesh, der, st, rs_state)
subroutine, public restart_read(restart, iunit, lines, nlines, ierr)
Definition: restart.F90:936
subroutine, public restart_close(restart, iunit)
Close a file previously opened with restart_open.
Definition: restart.F90:953
subroutine, public restart_write(restart, iunit, lines, nlines, ierr)
Definition: restart.F90:909
logical pure function, public restart_skip(restart)
Returns true if the restart information should neither be read nor written. This might happen because...
Definition: restart.F90:970
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,...
Definition: restart.F90:862
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
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.
Describes mesh distribution to nodes.
Definition: mesh.F90:186
int true(void)