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-2022 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.

SOURCE

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