45#if defined(HAVE_NETCDF)
79#if defined(HAVE_NETCDF)
86 integer,
parameter,
private :: &
91 character(len=3),
parameter :: &
92 index2label(3) = (/
're ',
'im ',
'abs' /)
117 what_tag_in, how_tag_in, output_interval_tag_in, ignore_error)
118 type(namespace_t),
intent(in) :: namespace
119 class(space_t),
intent(in) :: space
120 logical,
intent(inout) :: what(MAX_OUTPUT_TYPES)
121 integer(int64),
intent(out) :: how(0:MAX_OUTPUT_TYPES)
122 integer,
intent(out) :: output_interval(0:MAX_OUTPUT_TYPES)
123 character(len=*),
optional,
intent(in) :: what_tag_in
124 character(len=*),
optional,
intent(in) :: how_tag_in
125 character(len=*),
optional,
intent(in) :: output_interval_tag_in
126 logical,
optional,
intent(in) :: ignore_error
129 integer :: ncols, nrows, iout, column_index, what_i
130 integer(int64) :: what_no_how(12)
131 character(len=80) :: what_tag
132 character(len=80) :: how_tag
133 character(len=80) :: output_interval_tag
134 character(len=80) :: output_column_marker
138 output_interval = 0_8
142 output_interval_tag =
optional_default(output_interval_tag_in,
'OutputInterval')
471 what_no_how = (/ option__output__matrix_elements, option__output__berkeleygw, option__output__dos, &
472 option__output__tpa, option__output__mmb_den, option__output__j_flow, option__output__occ_matrices, &
473 option__output__effectiveu, option__output__magnetization, option__output__kanamoriu, option__output__stress, option__output__jdos /)
475 if (
parse_block(namespace, what_tag, blk) == 0)
then
477 do iout = 0, nrows - 1
493 what(what_i) = .
true.
494 call parse_variable(namespace, output_interval_tag, 50, output_interval(what_i))
495 if (((what_tag ==
'Output') .and. (.not. any(what_no_how == what_i)))&
496 .or. (what_tag /=
'Output'))
then
504 else if (ncols == 2)
then
517 what(what_i) = .
true.
518 call parse_variable(namespace, output_interval_tag, 50, output_interval(what_i))
519 if (((what_tag ==
'Output') .and. (.not. any(what_no_how == what_i)))&
520 .or. (what_tag /=
'Output'))
then
522 if (how(what_i) == 0)
call parse_variable(namespace, how_tag, 0, how(what_i))
549 what(what_i) = .
true.
550 do column_index = 0, 1
552 if (output_column_marker ==
'output_interval')
then
554 else if (output_column_marker ==
'output_format')
then
555 if (((what_tag ==
'Output') .and. (.not. any(what_no_how == what_i)))&
556 .or. (what_tag /=
'Output'))
then
562 else if (len_trim(output_column_marker) /= 0)
then
567 if (output_interval(what_i) == 0)
then
568 call parse_variable(namespace, output_interval_tag, 50, output_interval(what_i))
570 if (how(what_i) == 0)
then
571 if (((what_tag ==
'Output') .and. (.not. any(what_no_how == what_i)))&
572 .or. (what_tag /=
'Output'))
then
591 call parse_variable(namespace, output_interval_tag, 50, output_interval(0))
595 do what_i = lbound(what, 1), ubound(what, 1)
596 if (what_tag ==
'Output')
then
597 if (what(what_i) .and. (.not. any(what_no_how == what_i)))
then
602 if (how(what_i) == 0 .and. .not.
optional_default(ignore_error, .false.))
then
603 write(
message(1),
'(a)')
'Must specify output method with variable OutputFormat.'
608 if (space%dim == 1)
then
609 if (
bitand(how(what_i), option__outputformat__axis_y) /= 0)
then
610 message(1) =
"OutputFormat = axis_y not available with Dimensions = 1."
613 if (
bitand(how(what_i), option__outputformat__plane_z) /= 0)
then
614 message(1) =
"OutputFormat = plane_z not available with Dimensions = 1."
617 if (
bitand(how(what_i), option__outputformat__xcrysden) /= 0)
then
618 message(1) =
"OutputFormat = xcrysden not available with Dimensions = 1."
623 if (space%dim <= 2)
then
624 if (
bitand(how(what_i), option__outputformat__axis_z) /= 0)
then
625 message(1) =
"OutputFormat = axis_z not available with Dimensions <= 2."
628 if (
bitand(how(what_i), option__outputformat__plane_x) /= 0)
then
629 message(1) =
"OutputFormat = plane_x not available with Dimensions <= 2."
632 if (
bitand(how(what_i), option__outputformat__plane_y) /= 0)
then
633 message(1) =
"OutputFormat = plane_y not available with Dimensions <= 2."
636 if (
bitand(how(what_i), option__outputformat__integrate_xy) /= 0)
then
637 message(1) =
"OutputFormat = integrate_xy not available with Dimensions <= 2."
640 if (
bitand(how(what_i), option__outputformat__integrate_xz) /= 0)
then
641 message(1) =
"OutputFormat = integrate_xz not available with Dimensions <= 2."
644 if (
bitand(how(what_i), option__outputformat__integrate_yz) /= 0)
then
645 message(1) =
"OutputFormat = integrate_yz not available with Dimensions <= 2."
648 if (
bitand(how(what_i), option__outputformat__dx) /= 0)
then
649 message(1) =
"OutputFormat = dx not available with Dimensions <= 2."
652 if (
bitand(how(what_i), option__outputformat__cube) /= 0)
then
653 message(1) =
"OutputFormat = cube not available with Dimensions <= 2."
658#if !defined(HAVE_NETCDF)
659 if (
bitand(how(what_i), option__outputformat__netcdf) /= 0)
then
660 message(1) =
'Octopus was compiled without NetCDF support.'
661 message(2) =
'It is not possible to write output in NetCDF format.'
665#if !defined(HAVE_ETSF_IO)
666 if (
bitand(how(what_i), option__outputformat__etsf) /= 0)
then
667 message(1) =
'Octopus was compiled without ETSF_IO support.'
668 message(2) =
'It is not possible to write output in ETSF format.'
677 if (output_interval(what_i) < 0)
then
678 message(1) =
"OutputInterval must be >= 0."
692 character(len=*),
intent(in) :: where
697 if (index(
where,
"AxisX") /= 0) how = ior(how, option__outputformat__axis_x)
698 if (index(
where,
"AxisY") /= 0) how = ior(how, option__outputformat__axis_y)
699 if (index(
where,
"AxisZ") /= 0) how = ior(how, option__outputformat__axis_z)
700 if (index(
where,
"PlaneX") /= 0) how = ior(how, option__outputformat__plane_x)
701 if (index(
where,
"PlaneY") /= 0) how = ior(how, option__outputformat__plane_y)
702 if (index(
where,
"PlaneZ") /= 0) how = ior(how, option__outputformat__plane_z)
703 if (index(
where,
"IntegrateXY") /= 0) how = ior(how, option__outputformat__integrate_xy)
704 if (index(
where,
"IntegrateXZ") /= 0) how = ior(how, option__outputformat__integrate_xz)
705 if (index(
where,
"IntegrateYZ") /= 0) how = ior(how, option__outputformat__integrate_yz)
706 if (index(
where,
"PlaneZ") /= 0) how = ior(how, option__outputformat__plane_z)
707 if (index(
where,
"DX") /= 0) how = ior(how, option__outputformat__dx)
708 if (index(
where,
"XCrySDen") /= 0) how = ior(how, option__outputformat__xcrysden)
709 if (index(
where,
"Binary") /= 0) how = ior(how, option__outputformat__binary)
710 if (index(
where,
"MeshIndex") /= 0) how = ior(how, option__outputformat__mesh_index)
711 if (index(
where,
"XYZ") /= 0) how = ior(how, option__outputformat__xyz)
712#if defined(HAVE_NETCDF)
713 if (index(
where,
"NETCDF") /= 0) how = ior(how, option__outputformat__netcdf)
715 if (index(
where,
"Cube") /= 0) how = ior(how, option__outputformat__cube)
716 if (index(
where,
"VTK") /= 0) how = ior(how, option__outputformat__vtk)
727 character(len=*),
intent(in) :: dir
728 character(len=*),
intent(in) :: fname
729 class(
space_t),
intent(in) :: space
731 real(real64),
intent(in) :: pos(:,:)
732 type(
atom_t),
intent(in) :: atoms(:)
733 class(
box_t),
intent(in) :: box
738 real(real64) :: position
744 iunit =
io_open(trim(dir)//
'/'//trim(fname)//
'.xyz', namespace, action=
'write', position=
'asis')
746 write(iunit,
'(i6)')
size(atoms)
747 write(iunit,
'(a,a,a)', advance=
'no') trim(space%short_info()),
'; ', trim(box%short_info(
unit_angstrom))
748 if (space%is_periodic())
then
749 write(iunit,
'(a,a)')
'; ', trim(latt%short_info(
unit_angstrom))
755 do iatom = 1,
size(atoms)
756 write(iunit,
'(10a)', advance=
'no') atoms(iatom)%label
758 if (idir <= space%dim)
then
759 position = pos(idir, iatom)
775 character(len=*),
intent(in) :: dir, fname
776 class(
space_t),
intent(in) :: space
778 real(real64),
intent(in) :: pos(:,:)
779 type(
atom_t),
intent(in) :: atoms(:)
780 class(
mesh_t),
intent(in) :: mesh
782 real(real64),
optional,
intent(in) :: total_forces(:,:)
785 real(real64),
allocatable :: forces(:,:)
792 iunit =
io_open(trim(dir)//
'/'//trim(fname)//
'.xsf', namespace, action=
'write', position=
'asis')
794 if (
present(total_forces))
then
795 safe_allocate(forces(1:space%dim, 1:
size(atoms)))
798 safe_deallocate_a(forces)
812 integer,
intent(in) :: iunit
813 class(
space_t),
intent(in) :: space
815 real(real64),
intent(in) :: pos(:,:)
816 type(
atom_t),
intent(in) :: atoms(:)
817 class(
mesh_t),
intent(in) :: mesh
818 real(real64),
optional,
intent(in) :: forces(:, :)
819 integer,
optional,
intent(in) :: index
821 integer :: idir, idir2, iatom, index_, natoms
822 character(len=7) :: index_str
823 real(real64) :: offset(space%dim)
824 real(real64) :: rlattice(space%dim, space%dim)
828 if (
present(index))
then
829 write(index_str,
'(a,i6)')
' ', index
832 write(index_str,
'(a)')
''
835 natoms =
size(pos, dim=2)
841 offset(1:space%dim) = latt%red_to_cart(spread(-
m_half, 1, space%dim))
843 do idir = space%periodic_dim + 1, space%dim
844 offset(idir) = -(mesh%idx%ll(idir) - 1)/2 * mesh%spacing(idir)
847 if (space%is_periodic())
then
848 if (index_ == 1)
then
849 select case (space%periodic_dim)
851 write(iunit,
'(a)')
'CRYSTAL'
853 write(iunit,
'(a)')
'SLAB'
855 write(iunit,
'(a)')
'POLYMER'
859 write(iunit,
'(a)')
'PRIMVEC'//trim(index_str)
862 rlattice = latt%rlattice
863 do idir = space%periodic_dim+1, space%dim
864 rlattice(:,idir) = rlattice(:,idir)*
m_two*mesh%box%bounding_box_l(idir)
867 do idir = 1, space%dim
871 write(iunit,
'(a)')
'PRIMCOORD'//trim(index_str)
872 write(iunit,
'(i10, a)') natoms,
' 1'
874 write(iunit,
'(a)')
'ATOMS'//trim(index_str)
879 write(iunit,
'(a10, 3f12.6)', advance=
'no') trim(atoms(iatom)%label), &
881 if (
present(forces))
then
882 write(iunit,
'(5x, 3f12.6)', advance=
'no') forces(:, iatom)
894 integer,
intent(in) :: iunit
895 class(
space_t),
intent(in) :: space
897 real(real64),
intent(in) :: pos(:,:)
898 type(
atom_t),
intent(in) :: atoms(:)
899 class(
mesh_t),
intent(in) :: mesh
900 real(real64),
intent(in) :: centers(:, :)
901 integer,
intent(in) :: supercell(:)
902 real(real64),
optional,
intent(in) :: extra_atom(:)
904 integer :: idir, idir2, iatom, index_
905 character(len=7) :: index_str
906 real(real64) :: offset(3)
907 integer :: irep, nreplica
911 write(index_str,
'(a)')
''
914 nreplica = product(supercell(1:space%dim))
919 offset(1:space%dim) = latt%red_to_cart(spread(-
m_half, 1, space%dim))
920 offset(1:3) = offset(1:3) + centers(1:3,1)
922 do idir = space%periodic_dim + 1, 3
923 offset(idir) = -(mesh%idx%ll(idir) - 1)/2 * mesh%spacing(idir)
926 if(space%is_periodic())
then
928 select case(space%periodic_dim)
930 write(iunit,
'(a)')
'CRYSTAL'
932 write(iunit,
'(a)')
'SLAB'
934 write(iunit,
'(a)')
'POLYMER'
938 write(iunit,
'(a)')
'PRIMVEC'//trim(index_str)
940 do idir = 1, space%dim
942 latt%rlattice(idir2, idir)*supercell(idir)), idir2 = 1, space%dim)
945 write(iunit,
'(a)')
'PRIMCOORD'//trim(index_str)
946 if(.not.
present(extra_atom))
then
947 write(iunit,
'(i10, a)')
size(atoms)*nreplica,
' 1'
949 write(iunit,
'(i10, a)')
size(atoms)*nreplica+1,
' 1'
952 write(iunit,
'(a)')
'ATOMS'//trim(index_str)
956 do irep = 1, nreplica
958 do iatom = 1,
size(atoms)
959 write(iunit,
'(a10, 3f12.6)', advance=
'no') trim(atoms(iatom)%label), &
961 - offset(idir)), idir = 1, space%dim)
965 write(iunit,
'(a10, 3f12.6)', advance=
'no')
'X', &
973#if defined(HAVE_NETCDF)
975 subroutine ncdf_error(func, status, filename, namespace, ierr)
976 character(len=*),
intent(in) :: func
977 integer,
intent(in) :: status
978 character(len=*),
intent(in) :: filename
980 integer,
intent(inout) :: ierr
984 if (status == nf90_noerr)
then
989 write(
message(1),
'(3a)')
"NETCDF error in function '" , trim(func) ,
"'"
990 write(
message(2),
'(3a)')
"(reading/writing ", trim(filename) ,
")"
991 write(
message(3),
'(6x,a,a)')
'Error code = ', trim(nf90_strerror(status))
1001 real(real64),
intent(in) :: in(:, :, :)
1002 real(real64),
intent(out) :: out(:, :, :)
1003 integer :: ix, iy, iz
1007 do ix = lbound(in, 1), ubound(in, 1)
1008 do iy = lbound(in, 2), ubound(in, 2)
1009 do iz = lbound(in, 3), ubound(in, 3)
1010 out(iz, iy, ix) = in(ix, iy, iz)
1020#include "io_function_inc.F90"
1023#include "complex.F90"
1024#include "io_function_inc.F90"
Fast Fourier Transform module. This module provides a single interface that works with different FFT ...
real(real64), parameter, public m_two
real(real64), parameter, public m_zero
real(real64), parameter, public m_half
This module implements the index, used for the mesh points.
subroutine zio_function_output_vector_bz(how, dir, fname, namespace, space, kpt, kpoints, ff, unit, ierr, grp, root)
subroutine, public write_canonicalized_xyz_file(dir, fname, space, latt, pos, atoms, box, namespace)
Write canonicalized xyz file with atom labels and positions in Angstroms. Includes information about ...
subroutine, public zio_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.
subroutine, public dio_function_output_global(how, dir, fname, namespace, space, mesh, ff, unit, ierr)
subroutine, public io_function_read_what_how_when(namespace, space, what, how, output_interval, what_tag_in, how_tag_in, output_interval_tag_in, ignore_error)
subroutine, public zio_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,...
subroutine, public write_xsf_geometry(iunit, space, latt, pos, atoms, mesh, forces, index)
for format specification see: http:
subroutine zio_function_output_vector(how, dir, fname, namespace, space, mesh, ff, unit, ierr, pos, atoms, grp, root)
subroutine write_xsf_geometry_supercell(iunit, space, latt, pos, atoms, mesh, centers, supercell, extra_atom)
for format specification see: http:
subroutine dio_function_output_supercell(how, dir, fname, mesh, space, latt, ff, centers, supercell, unit, ierr, namespace, pos, atoms, grp, root, is_global, extra_atom)
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.
subroutine dio_function_output_vector(how, dir, fname, namespace, space, mesh, ff, unit, ierr, pos, atoms, grp, root)
subroutine transpose3(in, out)
subroutine zio_function_output_supercell(how, dir, fname, mesh, space, latt, ff, centers, supercell, unit, ierr, namespace, pos, atoms, grp, root, is_global, extra_atom)
subroutine, public write_xsf_geometry_file(dir, fname, space, latt, pos, atoms, mesh, namespace, total_forces)
subroutine, public dout_cf_netcdf(filename, ierr, cf, cube, space, spacing, transpose, unit, namespace)
Writes a cube_function in netcdf format.
subroutine dio_function_output_global_bz(how, dir, fname, namespace, kpoints, ff, unit, ierr)
subroutine ncdf_error(func, status, filename, namespace, ierr)
integer, parameter, private zoutput_kind
subroutine, public zio_function_output_global(how, dir, fname, namespace, space, mesh, ff, unit, ierr)
subroutine, public zout_cf_netcdf(filename, ierr, cf, cube, space, spacing, transpose, unit, namespace)
Writes a cube_function in netcdf format.
subroutine dio_function_output_vector_bz(how, dir, fname, namespace, space, kpt, kpoints, ff, unit, ierr, grp, root)
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...
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,...
subroutine zio_function_output_global_bz(how, dir, fname, namespace, kpoints, ff, unit, ierr)
subroutine, public io_close(iunit, grp)
subroutine, public io_mkdir(fname, namespace, parents)
integer function, public io_open(file, namespace, action, status, form, position, die, recl, grp)
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
subroutine, public messages_variable_is_block(namespace, name)
subroutine, public messages_warning(no_lines, all_nodes, namespace)
subroutine, public messages_obsolete_variable(namespace, name, rep)
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_input_error(namespace, var, details, row, column)
logical function mpi_grp_is_root(grp)
Is the current MPI process of grpcomm, root.
type(mpi_grp_t), public mpi_world
Some general things and nomenclature:
integer function, public parse_block(namespace, name, blk, check_varinfo_)
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.
type(unit_t), public unit_angstrom
For XYZ files.
type(unit_system_t), public units_out
This module is intended to contain simple general-purpose utility functions and procedures.
class to tell whether a point is inside or outside
Describes mesh distribution to nodes.