TABLE OF CONTENTS


ABINIT/m_wvl_denspot [ Modules ]

[ Top ] [ Modules ]

NAME

  m_wvl_denspot

FUNCTION

COPYRIGHT

  Copyright (C) 2008-2018 ABINIT group (DC)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_wvl_denspot
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32  use m_xmpi
33 
34  use m_geometry,   only : xred2xcart
35 
36  implicit none
37 
38  private

ABINIT/wvl_denspot_free [ Functions ]

[ Top ] [ Functions ]

NAME

  wvl_denspot_free

FUNCTION

INPUTS

OUTPUT

PARENTS

      gstate,wvl_wfsinp_reformat

CHILDREN

      deallocate_denspot_distribution,deallocate_rho_descriptors
      denspot_free_history,f_free_ptr

SOURCE

196 subroutine wvl_denspot_free(den)
197 
198  use defs_wvltypes
199 #if defined HAVE_BIGDFT
200  use BigDFT_API, only: deallocate_rho_descriptors, &
201       & deallocate_denspot_distribution, denspot_free_history
202  use dynamic_memory
203 #endif
204 
205 !This section has been created automatically by the script Abilint (TD).
206 !Do not modify the following lines by hand.
207 #undef ABI_FUNC
208 #define ABI_FUNC 'wvl_denspot_free'
209 !End of the abilint section
210 
211  implicit none
212 
213 !Arguments ------------------------------------
214  type(wvl_denspot_type), intent(inout) :: den
215 
216 !Local variables-------------------------------
217 
218 ! *************************************************************************
219 
220 !DEBUG
221 !write (std_out,*) ' wvl_denspot_free : enter'
222 !ENDDEBUG
223 
224 #if defined HAVE_BIGDFT
225  if(associated(den%denspot%rhov)) then
226    call f_free_ptr(den%denspot%rhov)
227  end if
228  if(associated(den%denspot%rho_psi)) then
229    call f_free_ptr(den%denspot%rho_psi)
230  end if
231  if(associated(den%denspot%rho_C)) then
232    call f_free_ptr(den%denspot%rho_C)
233  end if
234  if(associated(den%denspot%V_ext)) then
235    call f_free_ptr(den%denspot%V_ext)
236  end if
237  if(associated(den%denspot%V_XC)) then
238    call f_free_ptr(den%denspot%V_XC)
239  end if
240  if(associated(den%denspot%Vloc_KS)) then
241    call f_free_ptr(den%denspot%Vloc_KS)
242  end if
243  if(associated(den%denspot%f_XC)) then
244    call f_free_ptr(den%denspot%f_XC)
245  end if
246  if(associated(den%denspot%rho_work)) then
247    call f_free_ptr(den%denspot%rho_work)
248  end if
249  if(associated(den%denspot%pot_work)) then
250    call f_free_ptr(den%denspot%pot_work)
251  end if
252  nullify(den%denspot%rhov)
253  nullify(den%denspot%rho_psi)
254  nullify(den%denspot%rho_C)
255  nullify(den%denspot%V_ext)
256  nullify(den%denspot%V_XC)
257  nullify(den%denspot%Vloc_KS)
258  nullify(den%denspot%f_XC)
259  nullify(den%denspot%rho_work)
260  nullify(den%denspot%pot_work)
261  !
262  call deallocate_rho_descriptors(den%denspot%rhod)
263  call deallocate_denspot_distribution(den%denspot%dpbox)
264  call denspot_free_history(den%denspot)
265 #else
266  BIGDFT_NOTENABLED_ERROR()
267  if (.false.) write(std_out,*) den%symObj
268 #endif
269 
270 end subroutine wvl_denspot_free

ABINIT/wvl_denspot_set [ Functions ]

[ Top ] [ Functions ]

NAME

  wvl_denspot_set

FUNCTION

  Fill in denspot datatype with information related
  to density and potential data.

INPUTS

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      gstate,wvl_wfsinp_reformat

CHILDREN

      allocaterhopot,density_descriptors,dpbox_set
      initialize_dft_local_fields,wrtout,xred2xcart

SOURCE

 75 subroutine wvl_denspot_set(den,gth_params,ixc,natom,nsppol,rprimd,wvl,&
 76 &                          wvl_crmult,wvl_frmult,wvl_mpi_comm,xred)
 77 
 78  use defs_datatypes
 79  use defs_wvltypes
 80 
 81 #if defined HAVE_BIGDFT
 82  use BigDFT_API,only: initialize_DFT_local_fields,allocateRhoPot, &
 83 &                     input_variables,dpbox_set,density_descriptors
 84 #endif
 85 
 86 !This section has been created automatically by the script Abilint (TD).
 87 !Do not modify the following lines by hand.
 88 #undef ABI_FUNC
 89 #define ABI_FUNC 'wvl_denspot_set'
 90 !End of the abilint section
 91 
 92  implicit none
 93 
 94 !Arguments ------------------------------------
 95   integer,intent(in):: ixc,natom,nsppol,wvl_mpi_comm
 96   real(dp), intent(in) :: rprimd(3, 3)
 97   real(dp), intent(in) :: wvl_frmult,wvl_crmult
 98   real(dp), intent(inout)  :: xred(3,natom)
 99   type(wvl_denspot_type), intent(out) :: den
100   type(wvl_internal_type),intent(in)  :: wvl
101   type(pseudopotential_gth_type),intent(in)::gth_params
102 
103 !Local variables-------------------------------
104 #if defined HAVE_BIGDFT
105   integer :: groupsize,me,nproc
106   real(dp), allocatable :: xcart(:,:)
107   character(len=3),parameter :: rho_commun='DBL'
108   character(len=500) :: message
109   character(len=4) :: SICapproach
110   type(local_zone_descriptors) :: Lzd
111 #endif
112 
113   ! *************************************************************************
114 
115 !DEBUG
116 !write (std_out,*) ' wvl_denspot_set : enter'
117 !ENDDEBUG
118 
119 #if defined HAVE_BIGDFT
120 
121  write(message, '(a,a)' ) ch10,&
122 & ' wvl_denspot_set: Create wavelet type denspot.'
123  call wrtout(std_out,message,'COLL')
124 
125  nproc=xmpi_comm_size(wvl_mpi_comm)
126  me=xmpi_comm_rank(wvl_mpi_comm)
127  groupsize=0
128 
129 !Store xcart for each atom
130  ABI_ALLOCATE(xcart,(3, natom))
131  call xred2xcart(natom, rprimd, xcart, xred)
132 
133  call initialize_DFT_local_fields(den%denspot, ixc, nsppol)
134 
135 !number of planes for the density
136 !dpbox%nscatterarr(jproc, 1) = ngfft3_density
137 !number of planes for the potential
138 !dpbox%nscatterarr(jproc, 2) = ngfft3_potential
139 !starting offset for the potential
140 !dpbox%nscatterarr(jproc, 3) = density_start + potential_shift - 1
141 !GGA XC shift between density and potential
142 !dpbox%nscatterarr(jproc, 4) = potential_shift
143 
144  SICapproach="NONE"
145  Lzd%hgrids(1:3)=wvl%h(1:3)
146  Lzd%Glr%d%n1i=wvl%Glr%d%n1i
147  Lzd%Glr%d%n2i=wvl%Glr%d%n2i
148  Lzd%Glr%d%n3i=wvl%Glr%d%n3i
149  call dpbox_set(den%denspot%dpbox,Lzd,den%denspot%xc,me,nproc,wvl_mpi_comm,groupsize,&
150 & SICapproach,wvl%atoms%astruct%geocode,nsppol)
151 
152 !here dpbox can be put as input
153  call density_descriptors(me,nproc,den%denspot%xc,nsppol,wvl_crmult,wvl_frmult,wvl%atoms,&
154  den%denspot%dpbox,rho_commun,xcart,gth_params%radii_cf,den%denspot%rhod)
155 
156 !Note: change allocateRhoPot
157  call allocateRhoPot(wvl%Glr,nsppol,wvl%atoms,xcart,den%denspot)
158 
159 !Aditional informations.
160  den%symObj = wvl%atoms%astruct%sym%symObj
161 
162  ABI_DEALLOCATE(xcart)
163 
164 #else
165  BIGDFT_NOTENABLED_ERROR()
166  if (.false.) write(std_out,*) ixc,natom,nsppol,wvl_mpi_comm,rprimd(1,1),wvl_frmult,wvl_crmult,&
167 & xred(1,1),den%symObj,wvl%h(1),gth_params%psppar
168 #endif
169 
170 !DEBUG
171 !write (std_out,*) ' wvl_denspot_set : exit'
172 !ENDDEBUG
173 
174 end subroutine wvl_denspot_set