TABLE OF CONTENTS


ABINIT/vdw_kernelgen [ Programs ]

[ Top ] [ Programs ]

NAME

  vdw_kernelgen

FUNCTION

  Generates vdW-DF kernels from the user input.

COPYRIGHT

  Copyright (C) 2011-2018 ABINIT group (Yann Pouillon)
  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

  (main routine)

OUTPUT

  (main routine)

NOTES

  The input data must be provided in a pre-defined order and contain all
  adjustable parameters related to the generation of vdW-DF kernels.

PARENTS

CHILDREN

      abi_io_redirect,abimem_init,abinit_doctor,destroy_mpi_enreg,flush_unit
      herald,initmpi_seq,wrtout,xc_vdw_done,xc_vdw_get_params,xc_vdw_init
      xc_vdw_memcheck,xc_vdw_show,xc_vdw_write,xmpi_end,xmpi_init

SOURCE

 35 #if defined HAVE_CONFIG_H
 36 #include "config.h"
 37 #endif
 38 
 39 #include "abi_common.h"
 40 
 41 program vdw_kernelgen
 42 
 43 #if defined DEV_YP_VDWXC
 44  use defs_basis
 45  use defs_abitypes
 46  use m_build_info
 47  use m_errors
 48  use m_xc_vdw
 49  use m_mpinfo
 50  use m_xmpi
 51 #if defined HAVE_MPI2
 52  use mpi
 53 #endif
 54 
 55  use m_specialmsg,  only : specialmsg_getcount, herald
 56  use m_io_tools,    only : flush_unit
 57 
 58 !This section has been created automatically by the script Abilint (TD).
 59 !Do not modify the following lines by hand.
 60 #undef ABI_FUNC
 61 #define ABI_FUNC 'vdw_kernelgen'
 62 !End of the abilint section
 63 
 64  implicit none
 65 
 66 #if defined HAVE_MPI1
 67  include 'mpif.h'
 68 #endif
 69 
 70 !Arguments -----------------------------------
 71 
 72 !Local variables-------------------------------
 73 !no_abirules
 74 !
 75  character(len=24) :: codename
 76  character(len=500) :: message
 77  integer :: ierr
 78  type(MPI_type) :: mpi_enreg,mpi_enreg_seq
 79 
 80  type(xc_vdw_type) :: vdw_params
 81  character(len=fnlen) :: vdw_filnam
 82 
 83 #endif
 84 
 85 !******************************************************************
 86 !BEGIN EXECUTABLE SECTION
 87 
 88 #if defined DEV_YP_VDWXC
 89 
 90 !Change communicator for I/O (mandatory!)
 91  call abi_io_redirect(new_io_comm=xmpi_world)
 92 
 93 !Initialize MPI : one should write a separate routine -init_mpi_enreg-
 94 !for doing that !!
 95  call xmpi_init()
 96 
 97 !Default for sequential use
 98  call initmpi_seq(mpi_enreg)
 99 
100 !Signal MPI I/O compilation has been activated
101 #if defined HAVE_MPI_IO
102  if(xmpi_paral==0)then
103    write(message,'(3a)') &
104 &   '  In order to use MPI_IO, you must compile with the MPI flag ',ch10,&
105 &   '  Action : recompile your code with different CPP flags.'
106    MSG_ERROR(message)
107  end if
108 #endif
109 
110 !Initialize memory profiling if it is activated
111 !if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0"
112 !note that abimem.mocc files can easily be multiple GB in size so don't use this option normally
113 #ifdef HAVE_MEM_PROFILING
114  call abimem_init(0)
115 #endif
116 
117 !Other values of mpi_enreg are dataset dependent, and should NOT be initialized
118 !inside vdw_kernelgen.F90.
119 
120 !* Init fake MPI type with values for sequential case.
121  call initmpi_seq(MPI_enreg_seq)
122 
123  write(message,'(3a)') ch10,'vdW-DF functionals are not fully operational yet.',&
124 & ch10
125  MSG_ERROR(message)
126 
127 !=== Write greetings ===
128  codename='vdW_KernelGen'//repeat(' ',11)
129  call herald(codename,abinit_version,std_out)
130 !YP: calling dump_config() makes tests fail => commented
131 !call dump_config(std_out)
132 
133 !**********************************************************************
134 
135 !IMPORTANT: DO NOT TOUCH THE COMMENTS OF THE FOLLOWING BLOCK
136 
137 !Read input parameters
138 !%%% VDW-DF: BEGIN INPUT PARAMS %%%
139  read(*,*) vdw_params%functional
140  read(*,*) vdw_params%zab
141  read(*,*) vdw_params%ndpts
142  read(*,*) vdw_params%dcut
143  read(*,*) vdw_params%dratio
144  read(*,*) vdw_params%dsoft
145  read(*,*) vdw_params%phisoft
146  read(*,*) vdw_params%nqpts
147  read(*,*) vdw_params%qcut
148  read(*,*) vdw_params%qratio
149  read(*,*) vdw_params%nrpts
150  read(*,*) vdw_params%rcut
151  read(*,*) vdw_params%rsoft
152  read(*,*) vdw_params%ngpts
153  read(*,*) vdw_params%gcut
154  read(*,*) vdw_params%acutmin
155  read(*,*) vdw_params%aratio
156  read(*,*) vdw_params%damax
157  read(*,*) vdw_params%damin
158  read(*,*) vdw_params%nsmooth
159  read(*,*) vdw_params%tolerance
160  read(*,*) vdw_params%tweaks
161 !%%% VDW-DF: END INPUT PARAMS %%%
162 
163  vdw_filnam = repeat(' ',fnlen)
164  read(*,'(a)') vdw_filnam
165 
166  call xc_vdw_show(std_out,vdw_params)
167  call xc_vdw_init(vdw_params)
168  call xc_vdw_get_params(vdw_params)
169  call xc_vdw_show(std_out,vdw_params)
170  call xc_vdw_memcheck(std_out,vdw_params)
171  call xc_vdw_write(trim(vdw_filnam)//'.nc')
172  call xc_vdw_done(vdw_params)
173 
174 !**********************************************************************
175 
176  write(message,'(a,a,a)') ch10, &
177 & '+vdw_kernelgen : the run completed successfully ', &
178 & ch10
179  call wrtout(std_out,message,'COLL')
180  call flush_unit(std_out)
181 
182  call destroy_mpi_enreg(mpi_enreg)
183 
184  call abinit_doctor("__vdw_kernelgen")
185 
186  call xmpi_end()
187 #endif
188 
189  end program vdw_kernelgen