TABLE OF CONTENTS


ABINIT/consist [ Functions ]

[ Top ] [ Functions ]

NAME

 consist

FUNCTION

 Checking of the consistency between the values of input variables

COPYRIGHT

 Copyright (C) 2002-2017 ABINIT group (PCasek,FF,XG)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

  aim_dtset= the structured entity containing all input variables
  tstngr= information about the test on the ngrid input variable
  tstvpt= information about the test on the vpts input variable

OUTPUT

  (only checking : print error message and stop if there is a problem)

 WARNING
 This file does not follow the ABINIT coding rules (yet)

PARENTS

      adini

CHILDREN

SOURCE

 34 #if defined HAVE_CONFIG_H
 35 #include "config.h"
 36 #endif
 37 
 38 #include "abi_common.h"
 39 
 40 
 41 subroutine consist(aim_dtset,tstngr,tstvpt)
 42 
 43  use defs_basis
 44  use defs_aimprom
 45  use defs_abitypes
 46  use m_profiling_abi
 47  use m_errors
 48 
 49 !This section has been created automatically by the script Abilint (TD).
 50 !Do not modify the following lines by hand.
 51 #undef ABI_FUNC
 52 #define ABI_FUNC 'consist'
 53 !End of the abilint section
 54 
 55  implicit none
 56 
 57 !Arguments ------------------------------------
 58 !scalars
 59  integer,intent(in) :: tstngr,tstvpt
 60  type(aim_dataset_type),intent(in) :: aim_dtset
 61 
 62 !Local variables ------------------------------
 63 
 64 ! *********************************************************************
 65 
 66 !write(std_out,*) tstngr, tstvpt
 67 
 68  if (((aim_dtset%denout/=0).or.(aim_dtset%lapout/=0)).and.((tstngr < 1).or.(tstvpt < 2))) then
 69    MSG_ERROR('in input1 - I cannot do the output !')
 70  end if
 71  if ((aim_dtset%denout > 0).and.(aim_dtset%lapout>0)) then
 72    if (aim_dtset%denout/=aim_dtset%lapout) then
 73      write(std_out,*) 'ERROR in input - when both denout and lapout are positive non-zero,'
 74      write(std_out,*) 'they must be equal.'
 75      MSG_ERROR("Aborting now")
 76    end if
 77    if ((tstvpt < aim_dtset%denout+1).or.(tstngr < aim_dtset%denout)) then
 78      write(std_out,*) 'ERROR in input2 - I cannot do the output !'
 79      MSG_ERROR("Aborting now")
 80    end if
 81  elseif (aim_dtset%denout > 0) then
 82    if ((tstvpt < aim_dtset%denout+1).or.(tstngr < aim_dtset%denout)) then
 83      write(std_out,*) 'ERROR in input - I cannot do the output !'
 84      MSG_ERROR("Aborting now")
 85    end if
 86  elseif (aim_dtset%lapout > 0) then
 87    if ((tstvpt < aim_dtset%lapout+1).or.(tstngr < aim_dtset%lapout)) then
 88      write(std_out,*) 'ERROR in input - I cannot do the output !'
 89      MSG_ERROR("Aborting now")
 90    end if
 91  end if
 92 
 93  if ((aim_dtset%isurf==1).and.(aim_dtset%crit==0)) then
 94    write(std_out,*) 'ERROR in input - must have crit/=0 for isurf==1'
 95    MSG_ERROR("Aborting now")
 96  end if
 97 
 98  if (((aim_dtset%ivol/=0).or.(aim_dtset%irho/=0)).and.(aim_dtset%isurf==0)) then
 99    MSG_ERROR('in input - I cannot integrate without surface !')
100  end if
101 
102 end subroutine consist