gcc/testsuite/gfortran.dg/use_only_3.inc

   1     MODULE kinds
   2       INTEGER, PARAMETER :: DP = selected_real_kind(14,200)
   3       PRIVATE
   4       PUBLIC :: DP
   5     END MODULE kinds
   6
   7 MODULE constants
   8   USE kinds, ONLY : DP
   9   IMPLICIT NONE
  10   SAVE
  11   REAL(DP), PARAMETER :: pi = 3.14159265358979323846_DP
  12   REAL(DP), PARAMETER :: tpi= 2.0_DP * pi
  13   REAL(DP), PARAMETER :: fpi= 4.0_DP * pi
  14   REAL(DP), PARAMETER :: sqrtpi = 1.77245385090551602729_DP
  15   REAL(DP), PARAMETER :: sqrtpm1= 1.0_DP / sqrtpi
  16   REAL(DP), PARAMETER :: sqrt2  = 1.41421356237309504880_DP
  17   REAL(DP), PARAMETER :: H_PLANCK_SI      = 6.6260693D-34    ! J s
  18   REAL(DP), PARAMETER :: K_BOLTZMANN_SI   = 1.3806505D-23    ! J K^-1
  19   REAL(DP), PARAMETER :: ELECTRON_SI      = 1.60217653D-19   ! C
  20   REAL(DP), PARAMETER :: ELECTRONVOLT_SI  = 1.60217653D-19   ! J
  21   REAL(DP), PARAMETER :: ELECTRONMASS_SI  = 9.1093826D-31    ! Kg
  22   REAL(DP), PARAMETER :: HARTREE_SI       = 4.35974417D-18   ! J
  23   REAL(DP), PARAMETER :: RYDBERG_SI       = HARTREE_SI/2.0_DP! J
  24   REAL(DP), PARAMETER :: BOHR_RADIUS_SI   = 0.5291772108D-10 ! m
  25   REAL(DP), PARAMETER :: AMU_SI           = 1.66053886D-27   ! Kg
  26   REAL(DP), PARAMETER :: K_BOLTZMANN_AU   = K_BOLTZMANN_SI / HARTREE_SI
  27   REAL(DP), PARAMETER :: K_BOLTZMANN_RY   = K_BOLTZMANN_SI / RYDBERG_SI
  28   REAL(DP), PARAMETER :: AUTOEV           = HARTREE_SI / ELECTRONVOLT_SI
  29   REAL(DP), PARAMETER :: RYTOEV           = AUTOEV / 2.0_DP
  30   REAL(DP), PARAMETER :: AMU_AU           = AMU_SI / ELECTRONMASS_SI
  31   REAL(DP), PARAMETER :: AMU_RY           = AMU_AU / 2.0_DP
  32   REAL(DP), PARAMETER :: AU_SEC           = H_PLANCK_SI/tpi/HARTREE_SI
  33   REAL(DP), PARAMETER :: AU_PS            = AU_SEC * 1.0D+12
  34   REAL(DP), PARAMETER :: AU_GPA           = HARTREE_SI / BOHR_RADIUS_SI ** 3 &
  35                                             / 1.0D+9
  36   REAL(DP), PARAMETER :: RY_KBAR          = 10.0_dp * AU_GPA / 2.0_dp
  37   !
  38   REAL(DP), PARAMETER :: DEBYE_SI         = 3.3356409519 * 1.0D-30 ! C*m
  39   REAL(DP), PARAMETER :: AU_DEBYE         = ELECTRON_SI * BOHR_RADIUS_SI / &
  40                                             DEBYE_SI
  41   REAL(DP), PARAMETER :: eV_to_kelvin = ELECTRONVOLT_SI / K_BOLTZMANN_SI
  42   REAL(DP), PARAMETER :: ry_to_kelvin = RYDBERG_SI / K_BOLTZMANN_SI
  43   REAL(DP), PARAMETER :: eps4  = 1.0D-4
  44   REAL(DP), PARAMETER :: eps6  = 1.0D-6
  45   REAL(DP), PARAMETER :: eps8  = 1.0D-8
  46   REAL(DP), PARAMETER :: eps14 = 1.0D-14
  47   REAL(DP), PARAMETER :: eps16 = 1.0D-16
  48   REAL(DP), PARAMETER :: eps32 = 1.0D-32
  49   REAL(DP), PARAMETER :: gsmall = 1.0d-12
  50   REAL(DP), PARAMETER :: e2 = 2.D0      ! the square of the electron charge
  51   REAL(DP), PARAMETER :: degspin = 2.D0 ! the number of spins per level
  52   REAL(DP), PARAMETER :: amconv = AMU_RY
  53   REAL(DP), PARAMETER :: uakbar = RY_KBAR
  54   REAL(DP), PARAMETER :: bohr_radius_cm = bohr_radius_si * 100.0
  55   REAL(DP), PARAMETER :: BOHR_RADIUS_ANGS = bohr_radius_cm * 1.0D8
  56   REAL(DP), PARAMETER :: ANGSTROM_AU = 1.0/BOHR_RADIUS_ANGS
  57   REAL(DP), PARAMETER :: DIP_DEBYE = AU_DEBYE
  58   REAL(DP), PARAMETER :: AU_TERAHERTZ  = AU_PS
  59   REAL(DP), PARAMETER :: AU_TO_OHMCMM1 = 46000.0D0 ! (ohm cm)^-1
  60   !
  61
  62 END MODULE constants
  63
  64 !
  65 ! Copyright (C) 2001-2005 Quantum-ESPRESSO group
  66 ! This file is distributed under the terms of the
  67 ! GNU General Public License. See the file `License'
  68 ! in the root directory of the present distribution,
  69 ! or http://www.gnu.org/copyleft/gpl.txt .
  70 !
  71 !
  72 !---------------------------------------------------------------------------
  73 MODULE parameters
  74   !---------------------------------------------------------------------------
  75   !
  76   IMPLICIT NONE
  77   SAVE
  78   !
  79   INTEGER, PARAMETER :: &
  80        ntypx  = 10,     &! max number of different types of atom
  81        npsx   = ntypx,  &! max number of different PPs (obsolete)
  82        npk    = 40000,  &! max number of k-points
  83        lmaxx  = 3,      &! max non local angular momentum (l=0 to lmaxx)
  84        nchix  = 6,      &! max number of atomic wavefunctions per atom
  85        ndmx   = 2000     ! max number of points in the atomic radial mesh
  86   !
  87   INTEGER, PARAMETER :: &
  88        nbrx = 14,          &! max number of beta functions
  89        lqmax= 2*lmaxx+1,   &! max number of angular momenta of Q
  90        nqfx = 8             ! max number of coefficients in Q smoothing
  91   !
  92   INTEGER, PARAMETER :: nacx    = 10         ! max number of averaged
  93                                              ! quantities saved to the restart
  94   INTEGER, PARAMETER :: nsx     = ntypx      ! max number of species
  95   INTEGER, PARAMETER :: natx    = 5000       ! max number of atoms
  96   INTEGER, PARAMETER :: npkx    = npk        ! max number of K points
  97   INTEGER, PARAMETER :: ncnsx   = 101        ! max number of constraints
  98   INTEGER, PARAMETER :: nspinx  = 2          ! max number of spinors
  99   !
 100   INTEGER, PARAMETER :: nhclm   = 4  ! max number NH chain length, nhclm can be
 101                                      ! easily increased since the restart file
 102                                      ! should be able to handle it, perhaps
 103                                      ! better to align nhclm by 4
 104   !
 105   INTEGER, PARAMETER :: max_nconstr = 100
 106   !
 107   INTEGER, PARAMETER  ::  maxcpu = 2**17  ! Maximum number of CPU
 108   INTEGER, PARAMETER  ::  maxgrp = 128    ! Maximum number of task-groups
 109   !
 110 END MODULE parameters
 111
 112 MODULE control_flags
 113   USE kinds
 114   USE parameters
 115   IMPLICIT NONE
 116   SAVE
 117   TYPE convergence_criteria
 118      !
 119      LOGICAL  :: active
 120      INTEGER  :: nstep
 121      REAL(DP) :: ekin
 122      REAL(DP) :: derho
 123      REAL(DP) :: force
 124      !
 125   END TYPE convergence_criteria
 126   !
 127   TYPE ionic_conjugate_gradient
 128      !
 129      LOGICAL  :: active
 130      INTEGER  :: nstepix
 131      INTEGER  :: nstepex
 132      REAL(DP) :: ionthr
 133      REAL(DP) :: elethr
 134      !
 135   END TYPE ionic_conjugate_gradient
 136   !
 137   CHARACTER(LEN=4) :: program_name = ' '  !  used to control execution flow inside module
 138   !
 139   LOGICAL :: tvlocw    = .FALSE. ! write potential to unit 46 (only cp, seldom used)
 140   LOGICAL :: trhor     = .FALSE. ! read rho from      unit 47 (only cp, seldom used)
 141   LOGICAL :: trhow     = .FALSE. ! CP code, write rho to restart dir
 142   !
 143   LOGICAL :: tsde          = .FALSE. ! electronic steepest descent
 144   LOGICAL :: tzeroe        = .FALSE. ! set to zero the electronic velocities
 145   LOGICAL :: tfor          = .FALSE. ! move the ions ( calculate forces )
 146   LOGICAL :: tsdp          = .FALSE. ! ionic steepest descent
 147   LOGICAL :: tzerop        = .FALSE. ! set to zero the ionic velocities
 148   LOGICAL :: tprnfor       = .FALSE. ! print forces to standard output
 149   LOGICAL :: taurdr        = .FALSE. ! read ionic position from standard input
 150   LOGICAL :: tv0rd         = .FALSE. ! read ionic velocities from standard input
 151   LOGICAL :: tpre          = .FALSE. ! calculate stress, and (in fpmd) variable cell dynamic
 152   LOGICAL :: thdyn         = .FALSE. ! variable-cell dynamics (only cp)
 153   LOGICAL :: tsdc          = .FALSE. ! cell geometry steepest descent
 154   LOGICAL :: tzeroc        = .FALSE. ! set to zero the cell geometry velocities
 155   LOGICAL :: tstress       = .FALSE. ! print stress to standard output
 156   LOGICAL :: tortho        = .FALSE. ! use iterative orthogonalization
 157   LOGICAL :: tconjgrad     = .FALSE. ! use conjugate gradient electronic minimization
 158   LOGICAL :: timing        = .FALSE. ! print out timing information
 159   LOGICAL :: memchk        = .FALSE. ! check for memory leakage
 160   LOGICAL :: tprnsfac      = .FALSE. ! print out structure factor
 161   LOGICAL :: toptical      = .FALSE. ! print out optical properties
 162   LOGICAL :: tcarpar       = .FALSE. ! tcarpar is set TRUE for a "pure" Car Parrinello simulation
 163   LOGICAL :: tdamp         = .FALSE. ! Use damped dinamics for electrons
 164   LOGICAL :: tdampions     = .FALSE. ! Use damped dinamics for electrons
 165   LOGICAL :: tatomicwfc    = .FALSE. ! Use atomic wavefunctions as starting guess for ch. density
 166   LOGICAL :: tscreen       = .FALSE. ! Use screened coulomb potentials for cluster calculations
 167   LOGICAL :: twfcollect    = .FALSE. ! Collect wave function in the restart file at the end of run.
 168   LOGICAL :: tuspp         = .FALSE. ! Ultra-soft pseudopotential are being used
 169   INTEGER :: printwfc      = -1      ! Print wave functions, temporarely used only by ensemble-dft
 170   LOGICAL :: force_pairing = .FALSE. ! ...   Force pairing
 171   LOGICAL :: tchi2         = .FALSE. ! Compute Chi^2
 172   !
 173   TYPE (convergence_criteria) :: tconvthrs
 174                               !  thresholds used to check GS convergence
 175   !
 176   ! ... Ionic vs Electronic step frequency
 177   ! ... When "ion_nstep > 1" and "electron_dynamics = 'md' | 'sd' ", ions are
 178   ! ... propagated every "ion_nstep" electronic step only if the electronic
 179   ! ... "ekin" is lower than "ekin_conv_thr"
 180   !
 181   LOGICAL :: tionstep = .FALSE.
 182   INTEGER :: nstepe   = 1
 183                             !  parameters to control how many electronic steps
 184                             !  between ions move
 185
 186   LOGICAL :: tsteepdesc = .FALSE.
 187                             !  parameters for electronic steepest desceent
 188
 189   TYPE (ionic_conjugate_gradient) :: tconjgrad_ion
 190                             !  conjugate gradient for ionic minimization
 191
 192   INTEGER :: nbeg   = 0 ! internal code for initialization ( -1, 0, 1, 2, .. )
 193   INTEGER :: ndw    = 0 !
 194   INTEGER :: ndr    = 0 !
 195   INTEGER :: nomore = 0 !
 196   INTEGER :: iprint = 0 ! print output every iprint step
 197   INTEGER :: isave  = 0 ! write restart to ndr unit every isave step
 198   INTEGER :: nv0rd  = 0 !
 199   INTEGER :: iprsta = 0 ! output verbosity (increasing from 0 to infinity)
 200   !
 201   ! ... .TRUE. if only gamma point is used
 202   !
 203   LOGICAL :: gamma_only = .TRUE.
 204   !
 205   LOGICAL :: tnewnfi = .FALSE.
 206   INTEGER :: newnfi  = 0
 207   !
 208   ! This variable is used whenever a timestep change is requested
 209   !
 210   REAL(DP) :: dt_old = -1.0D0
 211   !
 212   ! ... Wave function randomization
 213   !
 214   LOGICAL  :: trane = .FALSE.
 215   REAL(DP) :: ampre = 0.D0
 216   !
 217   ! ... Ionic position randomization
 218   !
 219   LOGICAL  :: tranp(nsx) = .FALSE.
 220   REAL(DP) :: amprp(nsx) = 0.D0
 221   !
 222   ! ... Read the cell from standard input
 223   !
 224   LOGICAL :: tbeg = .FALSE.
 225   !
 226   ! ... This flags control the calculation of the Dipole Moments
 227   !
 228   LOGICAL :: tdipole = .FALSE.
 229   !
 230   ! ... Flags that controls DIIS electronic minimization
 231   !
 232   LOGICAL :: t_diis        = .FALSE.
 233   LOGICAL :: t_diis_simple = .FALSE.
 234   LOGICAL :: t_diis_rot    = .FALSE.
 235   !
 236   ! ... Flag controlling the Nose thermostat for electrons
 237   !
 238   LOGICAL :: tnosee = .FALSE.
 239   !
 240   ! ... Flag controlling the Nose thermostat for the cell
 241   !
 242   LOGICAL :: tnoseh = .FALSE.
 243   !
 244   ! ... Flag controlling the Nose thermostat for ions
 245   !
 246   LOGICAL  :: tnosep = .FALSE.
 247   LOGICAL  :: tcap   = .FALSE.
 248   LOGICAL  :: tcp    = .FALSE.
 249   REAL(DP) :: tolp   = 0.D0   !  tolerance for temperature variation
 250   !
 251   REAL(DP), PUBLIC :: &
 252        ekin_conv_thr = 0.D0, &!  conv. threshold for fictitious e. kinetic energy
 253        etot_conv_thr = 0.D0, &!  conv. threshold for DFT energy
 254        forc_conv_thr = 0.D0   !  conv. threshold for atomic forces
 255   INTEGER, PUBLIC :: &
 256        ekin_maxiter = 100,   &!  max number of iter. for ekin convergence
 257        etot_maxiter = 100,   &!  max number of iter. for etot convergence
 258        forc_maxiter = 100     !  max number of iter. for atomic forces conv.
 259   !
 260   ! ... Several variables controlling the run ( used mainly in PW calculations )
 261   !
 262   ! ... logical flags controlling the execution
 263   !
 264   LOGICAL, PUBLIC :: &
 265     lfixatom,           &! if .TRUE. some atom is kept fixed
 266     lscf,               &! if .TRUE. the calc. is selfconsistent
 267     lbfgs,              &! if .TRUE. the calc. is a relaxation based on new BFGS scheme
 268     lmd,                &! if .TRUE. the calc. is a dynamics
 269     lmetadyn,           &! if .TRUE. the calc. is a meta-dynamics
 270     lpath,              &! if .TRUE. the calc. is a path optimizations
 271     lneb,               &! if .TRUE. the calc. is NEB dynamics
 272     lsmd,               &! if .TRUE. the calc. is string dynamics
 273     lwf,                &! if .TRUE. the calc. is with wannier functions
 274     lphonon,            &! if .TRUE. the calc. is phonon
 275     lbands,             &! if .TRUE. the calc. is band structure
 276     lconstrain,         &! if .TRUE. the calc. is constraint
 277     ldamped,            &! if .TRUE. the calc. is a damped dynamics
 278     lrescale_t,         &! if .TRUE. the ionic temperature is rescaled
 279     langevin_rescaling, &! if .TRUE. the ionic dynamics is overdamped Langevin
 280     lcoarsegrained,     &! if .TRUE. a coarse-grained phase-space is used
 281     restart              ! if .TRUE. restart from results of a preceding run
 282   !
 283   LOGICAL, PUBLIC :: &
 284     remove_rigid_rot     ! if .TRUE. the total torque acting on the atoms is
 285                          ! removed
 286   !
 287   ! ... pw self-consistency
 288   !
 289   INTEGER, PUBLIC :: &
 290     ngm0,             &! used in mix_rho
 291     niter,            &! the maximum number of iteration
 292     nmix,             &! the number of iteration kept in the history
 293     imix               ! the type of mixing (0=plain,1=TF,2=local-TF)
 294   REAL(DP), PUBLIC  :: &
 295     mixing_beta,      &! the mixing parameter
 296     tr2                ! the convergence threshold for potential
 297   LOGICAL, PUBLIC :: &
 298     conv_elec          ! if .TRUE. electron convergence has been reached
 299   !
 300   ! ... pw diagonalization
 301   !
 302   REAL(DP), PUBLIC  :: &
 303     ethr               ! the convergence threshold for eigenvalues
 304   INTEGER, PUBLIC :: &
 305     david,            &! used on Davidson diagonalization
 306     isolve,           &! Davidson or CG or DIIS diagonalization
 307     max_cg_iter,      &! maximum number of iterations in a CG di
 308     diis_buff,        &! dimension of the buffer in diis
 309     diis_ndim          ! dimension of reduced basis in DIIS
 310   LOGICAL, PUBLIC :: &
 311     diago_full_acc     ! if true all the empty eigenvalues have the same
 312                        ! accuracy of the occupied ones
 313   !
 314   ! ... wfc and rho extrapolation
 315   !
 316   REAL(DP), PUBLIC  :: &
 317     alpha0,           &! the mixing parameters for the extrapolation
 318     beta0              ! of the starting potential
 319   INTEGER, PUBLIC :: &
 320     history,          &! number of old steps available for potential updating
 321     pot_order,        &! type of potential updating ( see update_pot )
 322     wfc_order          ! type of wavefunctions updating ( see update_pot )
 323   !
 324   ! ... ionic dynamics
 325   !
 326   INTEGER, PUBLIC :: &
 327     nstep,            &! number of ionic steps
 328     istep = 0          ! current ionic step
 329   LOGICAL, PUBLIC :: &
 330     conv_ions          ! if .TRUE. ionic convergence has been reached
 331   REAL(DP), PUBLIC  :: &
 332     upscale            ! maximum reduction of convergence threshold
 333   !
 334   ! ... system's symmetries
 335   !
 336   LOGICAL, PUBLIC :: &
 337     nosym,            &! if .TRUE. no symmetry is used
 338     noinv = .FALSE.    ! if .TRUE. eliminates inversion symmetry
 339   !
 340   ! ... phonon calculation
 341   !
 342   INTEGER, PUBLIC :: &
 343     modenum            ! for single mode phonon calculation
 344   !
 345   ! ... printout control
 346   !
 347   LOGICAL, PUBLIC :: &
 348     reduce_io          ! if .TRUE. reduce the I/O to the strict minimum
 349   INTEGER, PUBLIC :: &
 350     iverbosity         ! type of printing ( 0 few, 1 all )
 351   LOGICAL, PUBLIC :: &
 352     use_para_diago = .FALSE.  ! if .TRUE. a parallel Householder algorithm
 353   INTEGER, PUBLIC :: &
 354     para_diago_dim = 0        ! minimum matrix dimension above which a parallel
 355   INTEGER  :: ortho_max = 0    ! maximum number of iterations in routine ortho
 356   REAL(DP) :: ortho_eps = 0.D0 ! threshold for convergence in routine ortho
 357   LOGICAL, PUBLIC :: &
 358     use_task_groups = .FALSE.  ! if TRUE task groups parallelization is used
 359   INTEGER, PUBLIC :: iesr = 1
 360   LOGICAL,          PUBLIC :: tvhmean = .FALSE.
 361   REAL(DP),         PUBLIC :: vhrmin = 0.0d0
 362   REAL(DP),         PUBLIC :: vhrmax = 1.0d0
 363   CHARACTER(LEN=1), PUBLIC :: vhasse = 'Z'
 364   LOGICAL,          PUBLIC :: tprojwfc = .FALSE.
 365   CONTAINS
 366     SUBROUTINE fix_dependencies()
 367     END SUBROUTINE fix_dependencies
 368     SUBROUTINE check_flags()
 369     END SUBROUTINE check_flags
 370 END MODULE control_flags
 371
 372 !
 373 ! Copyright (C) 2002 FPMD group
 374 ! This file is distributed under the terms of the
 375 ! GNU General Public License. See the file `License'
 376 ! in the root directory of the present distribution,
 377 ! or http://www.gnu.org/copyleft/gpl.txt .
 378 !
 379
 380 !=----------------------------------------------------------------------------=!
 381    MODULE gvecw
 382 !=----------------------------------------------------------------------------=!
 383      USE kinds, ONLY: DP
 384
 385      IMPLICIT NONE
 386      SAVE
 387
 388      ! ...   G vectors less than the wave function cut-off ( ecutwfc )
 389      INTEGER :: ngw  = 0  ! local number of G vectors
 390      INTEGER :: ngwt = 0  ! in parallel execution global number of G vectors,
 391                        ! in serial execution this is equal to ngw
 392      INTEGER :: ngwl = 0  ! number of G-vector shells up to ngw
 393      INTEGER :: ngwx = 0  ! maximum local number of G vectors
 394      INTEGER :: ng0  = 0  ! first G-vector with nonzero modulus
 395                        ! needed in the parallel case (G=0 is on one node only!)
 396
 397      REAL(DP) :: ecutw = 0.0d0
 398      REAL(DP) :: gcutw = 0.0d0
 399
 400      !   values for costant cut-off computations
 401
 402      REAL(DP) :: ecfix = 0.0d0     ! value of the constant cut-off
 403      REAL(DP) :: ecutz = 0.0d0     ! height of the penalty function (above ecfix)
 404      REAL(DP) :: ecsig = 0.0d0     ! spread of the penalty function around ecfix
 405      LOGICAL   :: tecfix = .FALSE.  ! .TRUE. if constant cut-off is in use
 406
 407      ! augmented cut-off for k-point calculation
 408
 409      REAL(DP) :: ekcut = 0.0d0
 410      REAL(DP) :: gkcut = 0.0d0
 411
 412      ! array of G vectors module plus penalty function for constant cut-off
 413      ! simulation.
 414      !
 415      ! ggp = g + ( agg / tpiba**2 ) * ( 1 + erf( ( tpiba2 * g - e0gg ) / sgg ) )
 416
 417      REAL(DP), ALLOCATABLE, TARGET :: ggp(:)
 418
 419    CONTAINS
 420
 421      SUBROUTINE deallocate_gvecw
 422        IF( ALLOCATED( ggp ) ) DEALLOCATE( ggp )
 423      END SUBROUTINE deallocate_gvecw
 424
 425 !=----------------------------------------------------------------------------=!
 426    END MODULE gvecw
 427 !=----------------------------------------------------------------------------=!
 428
 429 !=----------------------------------------------------------------------------=!
 430    MODULE gvecs
 431 !=----------------------------------------------------------------------------=!
 432      USE kinds, ONLY: DP
 433
 434      IMPLICIT NONE
 435      SAVE
 436
 437      ! ...   G vectors less than the smooth grid cut-off ( ? )
 438      INTEGER :: ngs  = 0  ! local number of G vectors
 439      INTEGER :: ngst = 0  ! in parallel execution global number of G vectors,
 440                        ! in serial execution this is equal to ngw
 441      INTEGER :: ngsl = 0  ! number of G-vector shells up to ngw
 442      INTEGER :: ngsx = 0  ! maximum local number of G vectors
 443
 444      INTEGER, ALLOCATABLE :: nps(:), nms(:)
 445
 446      REAL(DP) :: ecuts = 0.0d0
 447      REAL(DP) :: gcuts = 0.0d0
 448
 449      REAL(DP) :: dual = 0.0d0
 450      LOGICAL   :: doublegrid = .FALSE.
 451
 452    CONTAINS
 453
 454      SUBROUTINE deallocate_gvecs()
 455        IF( ALLOCATED( nps ) ) DEALLOCATE( nps )
 456        IF( ALLOCATED( nms ) ) DEALLOCATE( nms )
 457      END SUBROUTINE deallocate_gvecs
 458
 459 !=----------------------------------------------------------------------------=!
 460    END MODULE gvecs
 461 !=----------------------------------------------------------------------------=!
 462
 463   MODULE electrons_base
 464       USE kinds, ONLY: DP
 465       IMPLICIT NONE
 466       SAVE
 467
 468       INTEGER :: nbnd       = 0    !  number electronic bands, each band contains
 469                                    !  two spin states
 470       INTEGER :: nbndx      = 0    !  array dimension nbndx >= nbnd
 471       INTEGER :: nspin      = 0    !  nspin = number of spins (1=no spin, 2=LSDA)
 472       INTEGER :: nel(2)     = 0    !  number of electrons (up, down)
 473       INTEGER :: nelt       = 0    !  total number of electrons ( up + down )
 474       INTEGER :: nupdwn(2)  = 0    !  number of states with spin up (1) and down (2)
 475       INTEGER :: iupdwn(2)  = 0    !  first state with spin (1) and down (2)
 476       INTEGER :: nudx       = 0    !  max (nupdw(1),nupdw(2))
 477       INTEGER :: nbsp       = 0    !  total number of electronic states
 478                                    !  (nupdwn(1)+nupdwn(2))
 479       INTEGER :: nbspx      = 0    !  array dimension nbspx >= nbsp
 480
 481       LOGICAL :: telectrons_base_initval = .FALSE.
 482       LOGICAL :: keep_occ = .FALSE.  ! if .true. when reading restart file keep
 483                                      ! the occupations calculated in initval
 484
 485       REAL(DP), ALLOCATABLE :: f(:)   ! occupation numbers ( at gamma )
 486       REAL(DP) :: qbac = 0.0d0        ! background neutralizing charge
 487       INTEGER, ALLOCATABLE :: ispin(:) ! spin of each state
 488 !
 489 !------------------------------------------------------------------------------!
 490   CONTAINS
 491 !------------------------------------------------------------------------------!
 492
 493
 494     SUBROUTINE electrons_base_initval( zv_ , na_ , nsp_ , nelec_ , nelup_ , neldw_ , nbnd_ , &
 495                nspin_ , occupations_ , f_inp, tot_charge_, multiplicity_, tot_magnetization_ )
 496       REAL(DP),         INTENT(IN) :: zv_ (:), tot_charge_
 497       REAL(DP),         INTENT(IN) :: nelec_ , nelup_ , neldw_
 498       REAL(DP),         INTENT(IN) :: f_inp(:,:)
 499       INTEGER,          INTENT(IN) :: na_ (:) , nsp_, multiplicity_, tot_magnetization_
 500       INTEGER,          INTENT(IN) :: nbnd_ , nspin_
 501       CHARACTER(LEN=*), INTENT(IN) :: occupations_
 502     END SUBROUTINE electrons_base_initval
 503
 504
 505     subroutine set_nelup_neldw ( nelec_, nelup_, neldw_, tot_magnetization_, &
 506          multiplicity_)
 507       !
 508       REAL (KIND=DP), intent(IN)    :: nelec_
 509       REAL (KIND=DP), intent(INOUT) :: nelup_, neldw_
 510       INTEGER,        intent(IN)    :: tot_magnetization_, multiplicity_
 511     end subroutine set_nelup_neldw
 512
 513 !----------------------------------------------------------------------------
 514
 515
 516     SUBROUTINE deallocate_elct()
 517       IF( ALLOCATED( f ) ) DEALLOCATE( f )
 518       IF( ALLOCATED( ispin ) ) DEALLOCATE( ispin )
 519       telectrons_base_initval = .FALSE.
 520       RETURN
 521     END SUBROUTINE deallocate_elct
 522
 523
 524 !------------------------------------------------------------------------------!
 525   END MODULE electrons_base
 526 !------------------------------------------------------------------------------!
 527
 528
 529
 530 !------------------------------------------------------------------------------!
 531   MODULE electrons_nose
 532 !------------------------------------------------------------------------------!
 533
 534       USE kinds, ONLY: DP
 535 !
 536       IMPLICIT NONE
 537       SAVE
 538
 539       REAL(DP) :: fnosee   = 0.0d0   !  frequency of the thermostat ( in THz )
 540       REAL(DP) :: qne      = 0.0d0   !  mass of teh termostat
 541       REAL(DP) :: ekincw   = 0.0d0   !  kinetic energy to be kept constant
 542
 543       REAL(DP) :: xnhe0   = 0.0d0
 544       REAL(DP) :: xnhep   = 0.0d0
 545       REAL(DP) :: xnhem   = 0.0d0
 546       REAL(DP) :: vnhe    = 0.0d0
 547   CONTAINS
 548   subroutine electrons_nose_init( ekincw_ , fnosee_ )
 549      REAL(DP), INTENT(IN) :: ekincw_, fnosee_
 550   end subroutine electrons_nose_init
 551
 552
 553   function electrons_nose_nrg( xnhe0, vnhe, qne, ekincw )
 554     real(8) :: electrons_nose_nrg
 555     real(8), intent(in) :: xnhe0, vnhe, qne, ekincw
 556     electrons_nose_nrg = 0.0
 557   end function electrons_nose_nrg
 558
 559   subroutine electrons_nose_shiftvar( xnhep, xnhe0, xnhem )
 560     implicit none
 561     real(8), intent(out) :: xnhem
 562     real(8), intent(inout) :: xnhe0
 563     real(8), intent(in) :: xnhep
 564   end subroutine electrons_nose_shiftvar
 565
 566   subroutine electrons_nosevel( vnhe, xnhe0, xnhem, delt )
 567     implicit none
 568     real(8), intent(inout) :: vnhe
 569     real(8), intent(in) :: xnhe0, xnhem, delt
 570   end subroutine electrons_nosevel
 571
 572   subroutine electrons_noseupd( xnhep, xnhe0, xnhem, delt, qne, ekinc, ekincw, vnhe )
 573     implicit none
 574     real(8), intent(out) :: xnhep, vnhe
 575     real(8), intent(in) :: xnhe0, xnhem, delt, qne, ekinc, ekincw
 576   end subroutine electrons_noseupd
 577
 578
 579   SUBROUTINE electrons_nose_info()
 580   END SUBROUTINE electrons_nose_info
 581   END MODULE electrons_nose
 582
 583 module cvan
 584   use parameters, only: nsx
 585   implicit none
 586   save
 587   integer nvb, ish(nsx)
 588   integer, allocatable:: indlm(:,:)
 589 contains
 590   subroutine allocate_cvan( nind, ns )
 591     integer, intent(in) :: nind, ns
 592   end subroutine allocate_cvan
 593
 594   subroutine deallocate_cvan( )
 595   end subroutine deallocate_cvan
 596
 597 end module cvan
 598
 599   MODULE cell_base
 600       USE kinds, ONLY : DP
 601       IMPLICIT NONE
 602       SAVE
 603         REAL(DP) :: alat = 0.0d0
 604         REAL(DP) :: celldm(6) = (/ 0.0d0, 0.0d0, 0.0d0, 0.0d0, 0.0d0, 0.0d0 /)
 605         REAL(DP) :: a1(3) = (/ 0.0d0, 0.0d0, 0.0d0 /)
 606         REAL(DP) :: a2(3) = (/ 0.0d0, 0.0d0, 0.0d0 /)
 607         REAL(DP) :: a3(3) = (/ 0.0d0, 0.0d0, 0.0d0 /)
 608         REAL(DP) :: b1(3) = (/ 0.0d0, 0.0d0, 0.0d0 /)
 609         REAL(DP) :: b2(3) = (/ 0.0d0, 0.0d0, 0.0d0 /)
 610         REAL(DP) :: b3(3) = (/ 0.0d0, 0.0d0, 0.0d0 /)
 611         REAL(DP) :: ainv(3,3) = 0.0d0
 612         REAl(DP) :: omega = 0.0d0  !  volume of the simulation cell
 613         REAL(DP) :: tpiba  = 0.0d0   !  = 2 PI / alat
 614         REAL(DP) :: tpiba2 = 0.0d0   !  = ( 2 PI / alat ) ** 2
 615         REAL(DP) :: at(3,3) = RESHAPE( (/ 0.0d0 /), (/ 3, 3 /), (/ 0.0d0 /) )
 616         REAL(DP) :: bg(3,3) = RESHAPE( (/ 0.0d0 /), (/ 3, 3 /), (/ 0.0d0 /) )
 617         INTEGER          :: ibrav      ! index of the bravais lattice
 618         CHARACTER(len=9) :: symm_type  ! 'cubic' or 'hexagonal' when ibrav=0
 619         REAL(DP) :: h(3,3)    = 0.0d0 ! simulation cell at time t
 620         REAL(DP) :: hold(3,3) = 0.0d0 ! simulation cell at time t-delt
 621         REAL(DP) :: hnew(3,3) = 0.0d0 ! simulation cell at time t+delt
 622         REAL(DP) :: velh(3,3) = 0.0d0 ! simulation cell velocity
 623         REAL(DP) :: deth      = 0.0d0 ! determinant of h ( cell volume )
 624         INTEGER   :: iforceh(3,3) = 1  ! if iforceh( i, j ) = 0 then h( i, j )
 625         LOGICAL   :: thdiag = .FALSE.  ! True if only cell diagonal elements
 626         REAL(DP) :: wmass = 0.0d0     ! cell fictitious mass
 627         REAL(DP) :: press = 0.0d0     ! external pressure
 628         REAL(DP) :: frich  = 0.0d0    ! firction parameter for cell damped dynamics
 629         REAL(DP) :: greash = 1.0d0    ! greas parameter for damped dynamics
 630         LOGICAL :: tcell_base_init = .FALSE.
 631   CONTAINS
 632         SUBROUTINE updatecell(box_tm1, box_t0, box_tp1)
 633           integer :: box_tm1, box_t0, box_tp1
 634         END SUBROUTINE updatecell
 635         SUBROUTINE dgcell( gcdot, box_tm1, box_t0, delt )
 636           REAL(DP), INTENT(OUT) :: GCDOT(3,3)
 637           REAL(DP), INTENT(IN) :: delt
 638           integer, intent(in) :: box_tm1, box_t0
 639         END SUBROUTINE dgcell
 640
 641         SUBROUTINE cell_init_ht( box, ht )
 642           integer :: box
 643           REAL(DP) :: ht(3,3)
 644         END SUBROUTINE cell_init_ht
 645
 646         SUBROUTINE cell_init_a( box, a1, a2, a3 )
 647           integer :: box
 648           REAL(DP) :: a1(3), a2(3), a3(3)
 649         END SUBROUTINE cell_init_a
 650
 651         SUBROUTINE r_to_s1 (r,s,box)
 652           REAL(DP), intent(out) ::  S(3)
 653           REAL(DP), intent(in) :: R(3)
 654           integer, intent(in) :: box
 655         END SUBROUTINE r_to_s1
 656
 657         SUBROUTINE r_to_s3 ( r, s, na, nsp, hinv )
 658           REAL(DP), intent(out) ::  S(:,:)
 659           INTEGER, intent(in) ::  na(:), nsp
 660           REAL(DP), intent(in) :: R(:,:)
 661           REAL(DP), intent(in) :: hinv(:,:)    ! hinv = TRANSPOSE( box%m1 )
 662           integer :: i, j, ia, is, isa
 663           isa = 0
 664           DO is = 1, nsp
 665             DO ia = 1, na(is)
 666               isa = isa + 1
 667               DO I=1,3
 668                 S(I,isa) = 0.D0
 669                 DO J=1,3
 670                   S(I,isa) = S(I,isa) + R(J,isa)*hinv(i,j)
 671                 END DO
 672               END DO
 673             END DO
 674           END DO
 675           RETURN
 676         END SUBROUTINE r_to_s3
 677
 678 !------------------------------------------------------------------------------!
 679
 680         SUBROUTINE r_to_s1b ( r, s, hinv )
 681           REAL(DP), intent(out) ::  S(:)
 682           REAL(DP), intent(in) :: R(:)
 683           REAL(DP), intent(in) :: hinv(:,:)    ! hinv = TRANSPOSE( box%m1 )
 684           integer :: i, j
 685           DO I=1,3
 686             S(I) = 0.D0
 687             DO J=1,3
 688               S(I) = S(I) + R(J)*hinv(i,j)
 689             END DO
 690           END DO
 691           RETURN
 692         END SUBROUTINE r_to_s1b
 693
 694
 695         SUBROUTINE s_to_r1 (S,R,box)
 696           REAL(DP), intent(in) ::  S(3)
 697           REAL(DP), intent(out) :: R(3)
 698           integer, intent(in) :: box
 699         END SUBROUTINE s_to_r1
 700
 701         SUBROUTINE s_to_r1b (S,R,h)
 702           REAL(DP), intent(in) ::  S(3)
 703           REAL(DP), intent(out) :: R(3)
 704           REAL(DP), intent(in) :: h(:,:)    ! h = TRANSPOSE( box%a )
 705         END SUBROUTINE s_to_r1b
 706
 707         SUBROUTINE s_to_r3 ( S, R, na, nsp, h )
 708           REAL(DP), intent(in) ::  S(:,:)
 709           INTEGER, intent(in) ::  na(:), nsp
 710           REAL(DP), intent(out) :: R(:,:)
 711           REAL(DP), intent(in) :: h(:,:)    ! h = TRANSPOSE( box%a )
 712         END SUBROUTINE s_to_r3
 713
 714       SUBROUTINE gethinv(box)
 715         IMPLICIT NONE
 716         integer, INTENT (INOUT) :: box
 717       END SUBROUTINE gethinv
 718
 719
 720       FUNCTION get_volume( hmat )
 721          IMPLICIT NONE
 722          REAL(DP) :: get_volume
 723          REAL(DP) :: hmat( 3, 3 )
 724           get_volume = 4.4
 725       END FUNCTION get_volume
 726
 727       FUNCTION pbc(rin,box,nl) RESULT (rout)
 728         IMPLICIT NONE
 729         integer :: box
 730         REAL (DP) :: rin(3)
 731         REAL (DP) :: rout(3), s(3)
 732         INTEGER, OPTIONAL :: nl(3)
 733         rout = 4.4
 734       END FUNCTION pbc
 735
 736           SUBROUTINE get_cell_param(box,cell,ang)
 737           IMPLICIT NONE
 738           integer, INTENT(in) :: box
 739           REAL(DP), INTENT(out), DIMENSION(3) :: cell
 740           REAL(DP), INTENT(out), DIMENSION(3), OPTIONAL :: ang
 741           END SUBROUTINE get_cell_param
 742
 743       SUBROUTINE pbcs_components(x1, y1, z1, x2, y2, z2, m)
 744         USE kinds
 745         INTEGER, INTENT(IN)  :: M
 746         REAL(DP),  INTENT(IN)  :: X1,Y1,Z1
 747         REAL(DP),  INTENT(OUT) :: X2,Y2,Z2
 748         REAL(DP) MIC
 749       END SUBROUTINE pbcs_components
 750
 751       SUBROUTINE pbcs_vectors(v, w, m)
 752         USE kinds
 753         INTEGER, INTENT(IN)  :: m
 754         REAL(DP),  INTENT(IN)  :: v(3)
 755         REAL(DP),  INTENT(OUT) :: w(3)
 756         REAL(DP) :: MIC
 757       END SUBROUTINE pbcs_vectors
 758
 759   SUBROUTINE cell_base_init( ibrav_ , celldm_ , trd_ht, cell_symmetry, rd_ht, cell_units, &
 760                a_ , b_ , c_ , cosab, cosac, cosbc, wc_ , total_ions_mass , press_ ,  &
 761                frich_ , greash_ , cell_dofree )
 762
 763     IMPLICIT NONE
 764     INTEGER, INTENT(IN) :: ibrav_
 765     REAL(DP), INTENT(IN) :: celldm_ (6)
 766     LOGICAL, INTENT(IN) :: trd_ht
 767     CHARACTER(LEN=*), INTENT(IN) :: cell_symmetry
 768     REAL(DP), INTENT(IN) :: rd_ht (3,3)
 769     CHARACTER(LEN=*), INTENT(IN) :: cell_units
 770     REAL(DP), INTENT(IN) :: a_ , b_ , c_ , cosab, cosac, cosbc
 771     CHARACTER(LEN=*), INTENT(IN) :: cell_dofree
 772     REAL(DP),  INTENT(IN) :: wc_ , frich_ , greash_ , total_ions_mass
 773     REAL(DP),  INTENT(IN) :: press_  ! external pressure from imput ( GPa )
 774   END SUBROUTINE cell_base_init
 775
 776
 777   SUBROUTINE cell_base_reinit( ht )
 778     REAL(DP), INTENT(IN) :: ht (3,3)
 779   END SUBROUTINE cell_base_reinit
 780
 781   SUBROUTINE cell_steepest( hnew, h, delt, iforceh, fcell )
 782     REAL(DP), INTENT(OUT) :: hnew(3,3)
 783     REAL(DP), INTENT(IN) :: h(3,3), fcell(3,3)
 784     INTEGER,      INTENT(IN) :: iforceh(3,3)
 785     REAL(DP), INTENT(IN) :: delt
 786   END SUBROUTINE cell_steepest
 787
 788   SUBROUTINE cell_verlet( hnew, h, hold, delt, iforceh, fcell, frich, tnoseh, hnos )
 789     REAL(DP), INTENT(OUT) :: hnew(3,3)
 790     REAL(DP), INTENT(IN) :: h(3,3), hold(3,3), hnos(3,3), fcell(3,3)
 791     INTEGER,      INTENT(IN) :: iforceh(3,3)
 792     REAL(DP), INTENT(IN) :: frich, delt
 793     LOGICAL,      INTENT(IN) :: tnoseh
 794   END SUBROUTINE cell_verlet
 795
 796   subroutine cell_hmove( h, hold, delt, iforceh, fcell )
 797     REAL(DP), intent(out) :: h(3,3)
 798     REAL(DP), intent(in) :: hold(3,3), fcell(3,3)
 799     REAL(DP), intent(in) :: delt
 800     integer, intent(in) :: iforceh(3,3)
 801   end subroutine cell_hmove
 802
 803   subroutine cell_force( fcell, ainv, stress, omega, press, wmass )
 804     REAL(DP), intent(out) :: fcell(3,3)
 805     REAL(DP), intent(in) :: stress(3,3), ainv(3,3)
 806     REAL(DP), intent(in) :: omega, press, wmass
 807   end subroutine cell_force
 808
 809   subroutine cell_move( hnew, h, hold, delt, iforceh, fcell, frich, tnoseh, vnhh, velh, tsdc )
 810     REAL(DP), intent(out) :: hnew(3,3)
 811     REAL(DP), intent(in) :: h(3,3), hold(3,3), fcell(3,3)
 812     REAL(DP), intent(in) :: vnhh(3,3), velh(3,3)
 813     integer,      intent(in) :: iforceh(3,3)
 814     REAL(DP), intent(in) :: frich, delt
 815     logical,      intent(in) :: tnoseh, tsdc
 816   end subroutine cell_move
 817
 818   subroutine cell_gamma( hgamma, ainv, h, velh )
 819     REAL(DP) :: hgamma(3,3)
 820     REAL(DP), intent(in) :: ainv(3,3), h(3,3), velh(3,3)
 821   end subroutine cell_gamma
 822
 823   subroutine cell_kinene( ekinh, temphh, velh )
 824     REAL(DP), intent(out) :: ekinh, temphh(3,3)
 825     REAL(DP), intent(in)  :: velh(3,3)
 826   end subroutine cell_kinene
 827
 828   function cell_alat( )
 829     real(DP) :: cell_alat
 830     cell_alat = 4.4
 831   end function cell_alat
 832    END MODULE cell_base
 833
 834
 835   MODULE ions_base
 836       USE kinds,      ONLY : DP
 837       USE parameters, ONLY : ntypx
 838       IMPLICIT NONE
 839       SAVE
 840       INTEGER :: nsp     = 0
 841       INTEGER :: na(5) = 0
 842       INTEGER :: nax     = 0
 843       INTEGER :: nat     = 0
 844       REAL(DP) :: zv(5)    = 0.0d0
 845       REAL(DP) :: pmass(5) = 0.0d0
 846       REAL(DP) :: amass(5) = 0.0d0
 847       REAL(DP) :: rcmax(5) = 0.0d0
 848       INTEGER,  ALLOCATABLE :: ityp(:)
 849       REAL(DP), ALLOCATABLE :: tau(:,:)     !  initial positions read from stdin (in bohr)
 850       REAL(DP), ALLOCATABLE :: vel(:,:)     !  initial velocities read from stdin (in bohr)
 851       REAL(DP), ALLOCATABLE :: tau_srt(:,:) !  tau sorted by specie in bohr
 852       REAL(DP), ALLOCATABLE :: vel_srt(:,:) !  vel sorted by specie in bohr
 853       INTEGER,  ALLOCATABLE :: ind_srt(:)   !  index of tau sorted by specie
 854       INTEGER,  ALLOCATABLE :: ind_bck(:)   !  reverse of ind_srt
 855       CHARACTER(LEN=3)      :: atm( 5 )
 856       CHARACTER(LEN=80)     :: tau_units
 857
 858
 859       INTEGER, ALLOCATABLE :: if_pos(:,:)  ! if if_pos( x, i ) = 0 then  x coordinate of
 860                                            ! the i-th atom will be kept fixed
 861       INTEGER, ALLOCATABLE :: iforce(:,:)  ! if_pos sorted by specie
 862       INTEGER :: fixatom   = -1            ! to be removed
 863       INTEGER :: ndofp     = -1            ! ionic degree of freedom
 864       INTEGER :: ndfrz     = 0             ! frozen degrees of freedom
 865
 866       REAL(DP) :: fricp   ! friction parameter for damped dynamics
 867       REAL(DP) :: greasp  ! friction parameter for damped dynamics
 868       REAL(DP), ALLOCATABLE :: taui(:,:)
 869       REAL(DP) :: cdmi(3), cdm(3)
 870       REAL(DP) :: cdms(3)
 871       LOGICAL :: tions_base_init = .FALSE.
 872   CONTAINS
 873     SUBROUTINE packtau( taup, tau, na, nsp )
 874       REAL(DP), INTENT(OUT) :: taup( :, : )
 875       REAL(DP), INTENT(IN) :: tau( :, :, : )
 876       INTEGER, INTENT(IN) :: na( : ), nsp
 877     END SUBROUTINE packtau
 878
 879     SUBROUTINE unpacktau( tau, taup, na, nsp )
 880       REAL(DP), INTENT(IN) :: taup( :, : )
 881       REAL(DP), INTENT(OUT) :: tau( :, :, : )
 882       INTEGER, INTENT(IN) :: na( : ), nsp
 883     END SUBROUTINE unpacktau
 884
 885     SUBROUTINE sort_tau( tausrt, isrt, tau, isp, nat, nsp )
 886       REAL(DP), INTENT(OUT) :: tausrt( :, : )
 887       INTEGER, INTENT(OUT) :: isrt( : )
 888       REAL(DP), INTENT(IN) :: tau( :, : )
 889       INTEGER, INTENT(IN) :: nat, nsp, isp( : )
 890       INTEGER :: ina( nsp ), na( nsp )
 891     END SUBROUTINE sort_tau
 892
 893     SUBROUTINE unsort_tau( tau, tausrt, isrt, nat )
 894       REAL(DP), INTENT(IN) :: tausrt( :, : )
 895       INTEGER, INTENT(IN) :: isrt( : )
 896       REAL(DP), INTENT(OUT) :: tau( :, : )
 897       INTEGER, INTENT(IN) :: nat
 898     END SUBROUTINE unsort_tau
 899
 900     SUBROUTINE ions_base_init( nsp_, nat_, na_, ityp_, tau_, vel_, amass_, &
 901                                atm_, if_pos_, tau_units_, alat_, a1_, a2_, &
 902                                a3_, rcmax_ )
 903       INTEGER,          INTENT(IN) :: nsp_, nat_, na_(:), ityp_(:)
 904       REAL(DP),         INTENT(IN) :: tau_(:,:)
 905       REAL(DP),         INTENT(IN) :: vel_(:,:)
 906       REAL(DP),         INTENT(IN) :: amass_(:)
 907       CHARACTER(LEN=*), INTENT(IN) :: atm_(:)
 908       CHARACTER(LEN=*), INTENT(IN) :: tau_units_
 909       INTEGER,          INTENT(IN) :: if_pos_(:,:)
 910       REAL(DP),         INTENT(IN) :: alat_, a1_(3), a2_(3), a3_(3)
 911       REAL(DP),         INTENT(IN) :: rcmax_(:)
 912     END SUBROUTINE ions_base_init
 913
 914     SUBROUTINE deallocate_ions_base()
 915     END SUBROUTINE deallocate_ions_base
 916
 917     SUBROUTINE ions_vel3( vel, taup, taum, na, nsp, dt )
 918       REAL(DP) :: vel(:,:), taup(:,:), taum(:,:)
 919       INTEGER :: na(:), nsp
 920       REAL(DP) :: dt
 921     END SUBROUTINE ions_vel3
 922
 923     SUBROUTINE ions_vel2( vel, taup, taum, nat, dt )
 924       REAL(DP) :: vel(:,:), taup(:,:), taum(:,:)
 925       INTEGER :: nat
 926       REAL(DP) :: dt
 927     END SUBROUTINE ions_vel2
 928
 929     SUBROUTINE cofmass1( tau, pmass, na, nsp, cdm )
 930       REAL(DP), INTENT(IN) :: tau(:,:,:), pmass(:)
 931       REAL(DP), INTENT(OUT) :: cdm(3)
 932       INTEGER, INTENT(IN) :: na(:), nsp
 933     END SUBROUTINE cofmass1
 934
 935     SUBROUTINE cofmass2( tau, pmass, na, nsp, cdm )
 936       REAL(DP), INTENT(IN) :: tau(:,:), pmass(:)
 937       REAL(DP), INTENT(OUT) :: cdm(3)
 938       INTEGER, INTENT(IN) :: na(:), nsp
 939     END SUBROUTINE cofmass2
 940
 941       SUBROUTINE randpos(tau, na, nsp, tranp, amprp, hinv, ifor )
 942          REAL(DP) :: hinv(3,3)
 943          REAL(DP) :: tau(:,:)
 944          INTEGER, INTENT(IN) :: ifor(:,:), na(:), nsp
 945          LOGICAL, INTENT(IN) :: tranp(:)
 946          REAL(DP), INTENT(IN) :: amprp(:)
 947          REAL(DP) :: oldp(3), rand_disp(3), rdisp(3)
 948
 949        END SUBROUTINE randpos
 950
 951   SUBROUTINE ions_kinene( ekinp, vels, na, nsp, h, pmass )
 952     REAL(DP), intent(out) :: ekinp     !  ionic kinetic energy
 953     REAL(DP), intent(in) :: vels(:,:)  !  scaled ionic velocities
 954     REAL(DP), intent(in) :: pmass(:)   !  ionic masses
 955     REAL(DP), intent(in) :: h(:,:)     !  simulation cell
 956     integer, intent(in) :: na(:), nsp
 957     integer :: i, j, is, ia, ii, isa
 958   END SUBROUTINE ions_kinene
 959
 960   subroutine ions_temp( tempp, temps, ekinpr, vels, na, nsp, h, pmass, ndega, nhpdim, atm2nhp, ekin2nhp )
 961     REAL(DP), intent(out) :: ekinpr, tempp
 962     REAL(DP), intent(out) :: temps(:)
 963     REAL(DP), intent(out) :: ekin2nhp(:)
 964     REAL(DP), intent(in)  :: vels(:,:)
 965     REAL(DP), intent(in)  :: pmass(:)
 966     REAL(DP), intent(in)  :: h(:,:)
 967     integer,        intent(in)  :: na(:), nsp, ndega, nhpdim, atm2nhp(:)
 968   end subroutine ions_temp
 969
 970   subroutine ions_thermal_stress( stress, pmass, omega, h, vels, nsp, na )
 971     REAL(DP), intent(inout) :: stress(3,3)
 972     REAL(DP), intent(in)  :: pmass(:), omega, h(3,3), vels(:,:)
 973     integer, intent(in) :: nsp, na(:)
 974     integer :: i, j, is, ia, isa
 975   end subroutine ions_thermal_stress
 976
 977   subroutine ions_vrescal( tcap, tempw, tempp, taup, tau0, taum, na, nsp, fion, iforce, &
 978                            pmass, delt )
 979     logical, intent(in) :: tcap
 980     REAL(DP), intent(inout) :: taup(:,:)
 981     REAL(DP), intent(in) :: tau0(:,:), taum(:,:), fion(:,:)
 982     REAL(DP), intent(in) :: delt, pmass(:), tempw, tempp
 983     integer, intent(in) :: na(:), nsp
 984     integer, intent(in) :: iforce(:,:)
 985   end subroutine ions_vrescal
 986   subroutine ions_shiftvar( varp, var0, varm )
 987     REAL(DP), intent(in) :: varp
 988     REAL(DP), intent(out) :: varm, var0
 989   end subroutine ions_shiftvar
 990    SUBROUTINE cdm_displacement( dis, tau )
 991       REAL(DP) :: dis
 992       REAL(DP) :: tau
 993    END SUBROUTINE cdm_displacement
 994    SUBROUTINE ions_displacement( dis, tau )
 995       REAL (DP), INTENT(OUT) :: dis
 996       REAL (DP), INTENT(IN)  :: tau
 997    END SUBROUTINE ions_displacement
 998   END MODULE ions_base