drift_bywlk_mod.f95

drift_bywlk_mod.f95
! -*- f90 -*-
!-----------------------------------------------------------------------------
!     This file is part of QmcMol A GNU QMC program for molecules.
!    Copyright (C) 2002 CNRS - U.P.M.C. Paris 6 - FRANCE
!     R. Assaraf, M. Caffarel, F. Colonna, X. Krokidis, B. Levy,
!     P. Pernod, and P. Reinhardt - qmcmol@lct.jussieu.fr
!            http://www.lct.jussieu.fr/QmcMol
! QmcMol is free software; you can redistribute it and/or modify it
! under the terms of the GNU General Public License as published by the Free
! Software Foundation; either version 2, or (at your option) any later
! version. QmcMol is distributed in the hope that it will be useful,but
! WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
! or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
! for more details. You should have received a copy of the GNU General
! Public License along with QmcMol; see the file COPYING. If not,
! write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
! Boston, MA 02111-1307, USA.
!-----------------------------------------------------------------------------

module drift_bywlk_mod

 use basic_types_tools_mod
 use psi_bywlk_mod
! use iterator_bywlk_mod

 implicit none

! Description    : Drift (e,:,w) = Grad Psi (e,:,w) / Psi (w)
!
! Creates      : drift_bywlk (elec, dim, walk)  

 real(dp), pointer :: drift_bywlk (:,:,:) => NULL ()

! Needs       : psi_bywlk (walk)
! Needs       : psi_grad_bywlk (elec, dim, walk)

 contains 

! ========================================================================
 subroutine drift_bywlk_bld
! ------------------------------------------------------------------------
!
! Name      : drift_bywlk_bld
!
! Description  : Drift (e,:,w) = Grad Psi (e,:,w) / Psi (w)
!
! Creates    : drift_bywlk (elec, dim, walk)
!
! Needs     : psi_bywlk (walk)
! Needs     : psi_grad_bywlk (elec, dim, walk)
!
! Author: P.Reinhardt
! Date : 26 Oct 2000, 03/02/2001
! Revision: 12 Feb 2002 by F. Colonna
! ------------------------------------------------------------------------
 implicit none

! i/o:

! local:
 character(len=max_string_len_routine_name),save :: l_here
 integer(i4b)         :: dimensions_nb
 integer(i4b)         :: electrons_nb 
 integer(i4b)         :: walkers_nb
 integer(i4b)         :: walk_i, dim_i, elec_i

! header:

 if (header_exe) then

  l_here = "drift_bywlk_bld"

  call object_create ('drift_bywlk')

  call object_needed (' psi_bywlk ')
  call object_needed (' psi_grad_bywlk ')

  return

 end if

! begin:
 call enter (l_here)

! requirements:

 include 'chksiz_psi_bywlk.inc'
 include 'chksiz_psi_grad_bywlk.inc'

 electrons_nb = size (psi_grad_bywlk, 1)
 dimensions_nb = size (psi_grad_bywlk, 2)
 walkers_nb  = size (psi_grad_bywlk, 3)

! end requirements

 if(debug) then
  write (logf,*)trim(l_here),'-d: dimensions_nb   = ',dimensions_nb
  write (logf,*)trim(l_here),'-d: electrons_nb   = ',electrons_nb
  write (logf,*)trim(l_here),'-d: walkers_nb    = ',walkers_nb
 end if

! allocations:

 if(.not.associated(drift_bywlk)) then
  if(trace_depend) &
   write (logf,*)trim(l_here),': drift_bywlk not yet associated'
   call alloc ('drift_bywlk', &
            drift_bywlk, &
            electrons_nb, &
            dimensions_nb, &
            walkers_nb)
 else 
  if(trace_depend) &
   write (logf,*)trim(l_here),': drift_bywlk left in memory'
 end if

! end allocations

 drift_bywlk = zero

 do walk_i = 1, walkers_nb
  do dim_i = 1, dimensions_nb
   do elec_i = 1, electrons_nb

    drift_bywlk (elec_i, dim_i, walk_i)    &
    = psi_grad_bywlk (elec_i, dim_i, walk_i) &
    / psi_bywlk (walk_i)

   end do ! elec_i 
  end do ! dim_i 
 end do ! walk_i 

 include 'debug_drift_bywlk.inc'
!  include 'ensnb_drift_bywlk.inc'

 call exit (l_here)

 end subroutine drift_bywlk_bld 

end module drift_bywlk_mod

back to top