TABLE OF CONTENTS


etsf_io_file_check_wavefunctions_data

[ Top ] [ etsf_io_file ] [ Methods ]

NAME

etsf_io_file_check_wavefunctions_data

FUNCTION

This is a high level routine to inquire a file about a specifications. It returns .true. in lstat if the file is a valid 'wavefunctions_data' file. It tests the existence of variables and their definition (type, shape. and dimension names).

INPUTS

OUTPUT

SOURCE

subroutine etsf_io_file_check_wavefunctions_data(ncid, lstat, error_data)
  integer, intent(intent)                  :: ncid
  logical, intent(out)                 :: lstat
  type(etsf_io_low_error), intent(out) :: error_data

  character(len = *), parameter        :: me = "etsf_io_file_check_wavefunctions_data"
  type(etsf_io_low_var_infos)          :: var_infos
  logical                              :: valid
  character(len = etsf_charlen)        :: string_value
  type(etsf_dims)                      :: dims
  type(etsf_split)                     :: split

  ! Read the dimensions
  call etsf_io_dims_get(ncid, dims, lstat, error_data)
  if (.not. lstat) then
     call etsf_io_low_error_update(error_data, me)
     return
  end if

  ! Allocate the split and read it (this will verify variable exist.
  call etsf_io_split_allocate(split, dims)
  call etsf_io_split_get(ncid, split, lstat, error_data)
  if (.not. lstat) then
     call etsf_io_low_error_update(error_data, me)
     return
  end if

  ! Variable primitive_vectors
  write(var_infos%name, "(A)") "primitive_vectors"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  write(var_infos%ncdimnames(2), "(A)") "number_of_vectors"
  write(var_infos%ncdimnames(1), "(A)") "number_of_cartesian_directions"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable reduced_symmetry_matrices
  write(var_infos%name, "(A)") "reduced_symmetry_matrices"
  var_infos%nctype  = etsf_io_low_integer
  var_infos%ncshape = 3
  allocate(var_infos%ncdimnames(3))
  write(var_infos%ncdimnames(3), "(A)") "number_of_symmetry_operations"
  write(var_infos%ncdimnames(2), "(A)") "number_of_reduced_dimensions"
  write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable reduced_symmetry_translations
  write(var_infos%name, "(A)") "reduced_symmetry_translations"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  write(var_infos%ncdimnames(2), "(A)") "number_of_symmetry_operations"
  write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable reduced_coordinates_of_kpoints
  write(var_infos%name, "(A)") "reduced_coordinates_of_kpoints"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  if (associated(split%my_kpoints)) then
    write(var_infos%ncdimnames(2), "(A)") "my_number_of_kpoints"
  else
    write(var_infos%ncdimnames(2), "(A)") "number_of_kpoints"
  end if
  write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable kpoint_weights
  write(var_infos%name, "(A)") "kpoint_weights"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 1
  allocate(var_infos%ncdimnames(1))
  if (associated(split%my_kpoints)) then
    write(var_infos%ncdimnames(1), "(A)") "my_number_of_kpoints"
  else
    write(var_infos%ncdimnames(1), "(A)") "number_of_kpoints"
  end if
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable number_of_states
  write(var_infos%name, "(A)") "number_of_states"
  var_infos%nctype  = etsf_io_low_integer
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  if (associated(split%my_spins)) then
    write(var_infos%ncdimnames(2), "(A)") "my_number_of_spins"
  else
    write(var_infos%ncdimnames(2), "(A)") "number_of_spins"
  end if
  if (associated(split%my_kpoints)) then
    write(var_infos%ncdimnames(1), "(A)") "my_number_of_kpoints"
  else
    write(var_infos%ncdimnames(1), "(A)") "number_of_kpoints"
  end if
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable eigenvalues
  write(var_infos%name, "(A)") "eigenvalues"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 3
  allocate(var_infos%ncdimnames(3))
  if (associated(split%my_spins)) then
    write(var_infos%ncdimnames(3), "(A)") "my_number_of_spins"
  else
    write(var_infos%ncdimnames(3), "(A)") "number_of_spins"
  end if
  if (associated(split%my_kpoints)) then
    write(var_infos%ncdimnames(2), "(A)") "my_number_of_kpoints"
  else
    write(var_infos%ncdimnames(2), "(A)") "number_of_kpoints"
  end if
  if (associated(split%my_states)) then
    write(var_infos%ncdimnames(1), "(A)") "my_max_number_of_states"
  else
    write(var_infos%ncdimnames(1), "(A)") "max_number_of_states"
  end if
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable occupations
  write(var_infos%name, "(A)") "occupations"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 3
  allocate(var_infos%ncdimnames(3))
  if (associated(split%my_spins)) then
    write(var_infos%ncdimnames(3), "(A)") "my_number_of_spins"
  else
    write(var_infos%ncdimnames(3), "(A)") "number_of_spins"
  end if
  if (associated(split%my_kpoints)) then
    write(var_infos%ncdimnames(2), "(A)") "my_number_of_kpoints"
  else
    write(var_infos%ncdimnames(2), "(A)") "number_of_kpoints"
  end if
  if (associated(split%my_states)) then
    write(var_infos%ncdimnames(1), "(A)") "my_max_number_of_states"
  else
    write(var_infos%ncdimnames(1), "(A)") "max_number_of_states"
  end if
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable basis_set
  write(var_infos%name, "(A)") "basis_set"
  var_infos%nctype  = etsf_io_low_character
  var_infos%ncshape = 1
  allocate(var_infos%ncdimnames(1))
  write(var_infos%ncdimnames(1), "(A)") "character_string_length"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Check these variables depends on the value of another.
  ! Read the condition value.
  call etsf_io_low_read_var(ncid, "basis_set", string_value, etsf_charlen, &
                          & lstat, error_data = error_data)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  call strip(string_value)
  if (trim(string_value) == "real_space") then
    ! Variable real_space_wavefunctions
    write(var_infos%name, "(A)") "real_space_wavefunctions"
    var_infos%nctype  = etsf_io_low_double
    var_infos%ncshape = 8
    allocate(var_infos%ncdimnames(8))
    if (associated(split%my_spins)) then
      write(var_infos%ncdimnames(8), "(A)") "my_number_of_spins"
    else
      write(var_infos%ncdimnames(8), "(A)") "number_of_spins"
    end if
    if (associated(split%my_kpoints)) then
      write(var_infos%ncdimnames(7), "(A)") "my_number_of_kpoints"
    else
      write(var_infos%ncdimnames(7), "(A)") "number_of_kpoints"
    end if
    if (associated(split%my_states)) then
      write(var_infos%ncdimnames(6), "(A)") "my_max_number_of_states"
    else
      write(var_infos%ncdimnames(6), "(A)") "max_number_of_states"
    end if
    write(var_infos%ncdimnames(5), "(A)") "number_of_spinor_components"
    if (associated(split%my_grid_points_vector3)) then
      write(var_infos%ncdimnames(4), "(A)") "my_number_of_grid_points_vect3"
    else
      write(var_infos%ncdimnames(4), "(A)") "number_of_grid_points_vector3"
    end if
    if (associated(split%my_grid_points_vector2)) then
      write(var_infos%ncdimnames(3), "(A)") "my_number_of_grid_points_vect2"
    else
      write(var_infos%ncdimnames(3), "(A)") "number_of_grid_points_vector2"
    end if
    if (associated(split%my_grid_points_vector1)) then
      write(var_infos%ncdimnames(2), "(A)") "my_number_of_grid_points_vect1"
    else
      write(var_infos%ncdimnames(2), "(A)") "number_of_grid_points_vector1"
    end if
    write(var_infos%ncdimnames(1), "(A)") "real_or_complex_wavefunctions"
    call test_var(ncid, var_infos, lstat, error_data)
    deallocate(var_infos%ncdimnames)
    if (.not. lstat) then
      call etsf_io_split_free(split)
      call etsf_io_low_error_update(error_data, me)
      return
    end if
    
  else if (trim(string_value) == "daubechies_wavelets") then
    ! Variable coordinates_of_basis_grid_points
    write(var_infos%name, "(A)") "coordinates_of_basis_grid_points"
    var_infos%nctype  = etsf_io_low_integer
    var_infos%ncshape = 3
    allocate(var_infos%ncdimnames(3))
    write(var_infos%ncdimnames(3), "(A)") "number_of_localization_regions"
    write(var_infos%ncdimnames(2), "(A)") "max_number_of_basis_grid_points"
    write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
    call test_var(ncid, var_infos, lstat, error_data)
    deallocate(var_infos%ncdimnames)
    if (.not. lstat) then
      call etsf_io_split_free(split)
      call etsf_io_low_error_update(error_data, me)
      return
    end if
    
    ! Variable number_of_coefficients_per_grid_point
    write(var_infos%name, "(A)") "number_of_coefficients_per_grid_point"
    var_infos%nctype  = etsf_io_low_integer
    var_infos%ncshape = 2
    allocate(var_infos%ncdimnames(2))
    write(var_infos%ncdimnames(2), "(A)") "number_of_localization_regions"
    write(var_infos%ncdimnames(1), "(A)") "max_number_of_basis_grid_points"
    call test_var(ncid, var_infos, lstat, error_data)
    deallocate(var_infos%ncdimnames)
    if (.not. lstat) then
      call etsf_io_split_free(split)
      call etsf_io_low_error_update(error_data, me)
      return
    end if
    
    ! Variable coefficients_of_wavefunctions
    write(var_infos%name, "(A)") "coefficients_of_wavefunctions"
    var_infos%nctype  = etsf_io_low_double
    var_infos%ncshape = 6
    allocate(var_infos%ncdimnames(6))
    if (associated(split%my_spins)) then
      write(var_infos%ncdimnames(6), "(A)") "my_number_of_spins"
    else
      write(var_infos%ncdimnames(6), "(A)") "number_of_spins"
    end if
    if (associated(split%my_kpoints)) then
      write(var_infos%ncdimnames(5), "(A)") "my_number_of_kpoints"
    else
      write(var_infos%ncdimnames(5), "(A)") "number_of_kpoints"
    end if
    if (associated(split%my_states)) then
      write(var_infos%ncdimnames(4), "(A)") "my_max_number_of_states"
    else
      write(var_infos%ncdimnames(4), "(A)") "max_number_of_states"
    end if
    write(var_infos%ncdimnames(3), "(A)") "number_of_spinor_components"
    if (associated(split%my_coefficients)) then
      write(var_infos%ncdimnames(2), "(A)") "my_max_number_of_coefficients"
    else
      write(var_infos%ncdimnames(2), "(A)") "max_number_of_coefficients"
    end if
    write(var_infos%ncdimnames(1), "(A)") "real_or_complex_coefficients"
    call test_var(ncid, var_infos, lstat, error_data)
    deallocate(var_infos%ncdimnames)
    if (.not. lstat) then
      call etsf_io_split_free(split)
      call etsf_io_low_error_update(error_data, me)
      return
    end if
    
  else if (trim(string_value) == "plane_waves") then
    ! Variable number_of_coefficients
    write(var_infos%name, "(A)") "number_of_coefficients"
    var_infos%nctype  = etsf_io_low_integer
    var_infos%ncshape = 1
    allocate(var_infos%ncdimnames(1))
    if (associated(split%my_kpoints)) then
      write(var_infos%ncdimnames(1), "(A)") "my_number_of_kpoints"
    else
      write(var_infos%ncdimnames(1), "(A)") "number_of_kpoints"
    end if
    call test_var(ncid, var_infos, lstat, error_data)
    deallocate(var_infos%ncdimnames)
    if (.not. lstat) then
      call etsf_io_split_free(split)
      call etsf_io_low_error_update(error_data, me)
      return
    end if
    
    ! Variable reduced_coordinates_of_plane_waves
    write(var_infos%name, "(A)") "reduced_coordinates_of_plane_waves"
    var_infos%nctype  = etsf_io_low_integer
    call etsf_io_low_read_flag(ncid, valid, "reduced_coordinates_of_plane_waves", &
                             & "k_dependent", lstat, error_data = error_data)
    if (valid) then
      var_infos%ncshape = 3
      allocate(var_infos%ncdimnames(3))
      if (associated(split%my_kpoints)) then
        write(var_infos%ncdimnames(3), "(A)") "my_number_of_kpoints"
      else
        write(var_infos%ncdimnames(3), "(A)") "number_of_kpoints"
      end if
      if (associated(split%my_coefficients)) then
        write(var_infos%ncdimnames(2), "(A)") "my_max_number_of_coefficients"
      else
        write(var_infos%ncdimnames(2), "(A)") "max_number_of_coefficients"
      end if
      write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
    else
      var_infos%ncshape = 2
      allocate(var_infos%ncdimnames(2))
      if (associated(split%my_coefficients)) then
        write(var_infos%ncdimnames(2), "(A)") "my_max_number_of_coefficients"
      else
        write(var_infos%ncdimnames(2), "(A)") "max_number_of_coefficients"
      end if
      write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
    end if
    call test_var(ncid, var_infos, lstat, error_data)
    deallocate(var_infos%ncdimnames)
    if (.not. lstat) then
      call etsf_io_split_free(split)
      call etsf_io_low_error_update(error_data, me)
      return
    end if
    
    ! Variable coefficients_of_wavefunctions
    write(var_infos%name, "(A)") "coefficients_of_wavefunctions"
    var_infos%nctype  = etsf_io_low_double
    var_infos%ncshape = 6
    allocate(var_infos%ncdimnames(6))
    if (associated(split%my_spins)) then
      write(var_infos%ncdimnames(6), "(A)") "my_number_of_spins"
    else
      write(var_infos%ncdimnames(6), "(A)") "number_of_spins"
    end if
    if (associated(split%my_kpoints)) then
      write(var_infos%ncdimnames(5), "(A)") "my_number_of_kpoints"
    else
      write(var_infos%ncdimnames(5), "(A)") "number_of_kpoints"
    end if
    if (associated(split%my_states)) then
      write(var_infos%ncdimnames(4), "(A)") "my_max_number_of_states"
    else
      write(var_infos%ncdimnames(4), "(A)") "max_number_of_states"
    end if
    write(var_infos%ncdimnames(3), "(A)") "number_of_spinor_components"
    if (associated(split%my_coefficients)) then
      write(var_infos%ncdimnames(2), "(A)") "my_max_number_of_coefficients"
    else
      write(var_infos%ncdimnames(2), "(A)") "max_number_of_coefficients"
    end if
    write(var_infos%ncdimnames(1), "(A)") "real_or_complex_coefficients"
    call test_var(ncid, var_infos, lstat, error_data)
    deallocate(var_infos%ncdimnames)
    if (.not. lstat) then
      call etsf_io_split_free(split)
      call etsf_io_low_error_update(error_data, me)
      return
    end if
    
  else
    call etsf_io_split_free(split)
    call etsf_io_low_error_set(error_data, ERROR_MODE_SPEC, &
                             & ERROR_TYPE_ARG, me, &
                             & tgtname = "basis_set", &
                             & errmess = "Empty or unknown value '"//trim(string_value)//"'.")
    lstat = .false.
    return
  end if


  ! Deallocate the split data.
  call etsf_io_split_free(split)

  lstat = .true.
end subroutine etsf_io_file_check_wavefunctions_data