TABLE OF CONTENTS


ABINIT/m_wvl_projectors [ Modules ]

[ Top ] [ Modules ]

NAME

  m_wvl_projectors

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_projectors
28 
29  use defs_basis
30  use defs_datatypes
31  use defs_wvltypes
32  use m_errors
33  use m_abicore
34  use m_atomdata
35 
36  use m_geometry,     only : xred2xcart
37 
38  implicit none
39 
40  private

ABINIT/wvl_projectors_free [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_projectors_free

FUNCTION

 Freeing routine.

INPUTS

OUTPUT

SIDE EFFECTS

  proj <type(wvl_projectors_type)>=projectors informations in a wavelet basis.

PARENTS

      gstate,wvl_wfsinp_reformat

CHILDREN

      free_dft_psp_projectors

SOURCE

193 subroutine wvl_projectors_free(proj)
194 
195 #if defined HAVE_BIGDFT
196   use BigDFT_API, only : free_DFT_PSP_projectors,deallocate_gwf
197 #endif
198 
199 !This section has been created automatically by the script Abilint (TD).
200 !Do not modify the following lines by hand.
201 #undef ABI_FUNC
202 #define ABI_FUNC 'wvl_projectors_free'
203 !End of the abilint section
204 
205  implicit none
206 
207 !Arguments ------------------------------------
208 !scalars
209  type(wvl_projectors_type),intent(inout) :: proj
210 
211 !Local variables -------------------------
212 #if defined HAVE_BIGDFT
213  integer :: ii
214 #endif
215 
216   ! *********************************************************************
217 
218 #if defined HAVE_BIGDFT
219 
220  call free_DFT_PSP_projectors(proj%nlpsp)
221 
222  if (allocated(proj%G)) then
223    do ii=1,size(proj%G)
224 !    MT dec 2014: cannot call bigdft deallocation routine
225 !    because content of proj%G datastructure was created
226 !    without f_malloc (without memory profiling).
227 !    call deallocate_gwf(proj%G(ii))
228      if (associated(proj%G(ii)%ndoc)) then
229        ABI_DEALLOCATE(proj%G(ii)%ndoc)
230      end if
231      if (associated(proj%G(ii)%nam)) then
232        ABI_DEALLOCATE(proj%G(ii)%nam)
233      end if
234      if (associated(proj%G(ii)%nshell)) then
235        ABI_DEALLOCATE(proj%G(ii)%nshell)
236      end if
237      if (associated(proj%G(ii)%psiat)) then
238        ABI_DEALLOCATE(proj%G(ii)%psiat)
239      end if
240      if (associated(proj%G(ii)%xp)) then
241        ABI_DEALLOCATE(proj%G(ii)%xp)
242      end if
243    end do
244    ABI_DATATYPE_DEALLOCATE(proj%G)
245  end if
246 
247 #else
248  BIGDFT_NOTENABLED_ERROR()
249  if (.false.) write(std_out,*) proj%nlpsp
250 #endif
251 
252 end subroutine wvl_projectors_free

ABINIT/wvl_projectors_set [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_projectors_set

FUNCTION

 Allocate and compute the access keys for the projectors when the positions
 of the atoms are given. The array to store projectors
 is also allocated, use wvl_projectors_free() to free them after use.

INPUTS

  dtset <type(dataset_type)>=internal variables used by wavelets, describing
   | wvl_internal=desciption of the wavelet box.
   | natom=number of atoms.
  mpi_enreg=informations about MPI parallelization
  psps <type(pseudopotential_type)>=variables related to pseudopotentials
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  xred(3,natom)=reduced dimensionless atomic coordinates

OUTPUT

  proj <type(wvl_projector_type)>=projectors informations for wavelets.
   | keys=its access keys for compact storage.

SIDE EFFECTS

PARENTS

      gstate,wvl_wfsinp_reformat

CHILDREN

      createprojectorsarrays,wrtout,wvl_timing,xred2xcart

SOURCE

 83 subroutine wvl_projectors_set(me, natom, proj, psps, rprimd, wfs, wvl, wvl_frmult, xred)
 84 
 85 #if defined HAVE_BIGDFT
 86  use BigDFT_API, only: createProjectorsArrays, wvl_timing => timing
 87 #endif
 88 
 89 !This section has been created automatically by the script Abilint (TD).
 90 !Do not modify the following lines by hand.
 91 #undef ABI_FUNC
 92 #define ABI_FUNC 'wvl_projectors_set'
 93 !End of the abilint section
 94 
 95  implicit none
 96 
 97 !Arguments ------------------------------------
 98 !scalars
 99  integer, intent(in) :: natom, me
100  real(dp), intent(in) :: wvl_frmult
101  type(pseudopotential_type),intent(in) :: psps
102  type(wvl_projectors_type),intent(inout) :: proj
103  type(wvl_wf_type),intent(in) :: wfs
104  type(wvl_internal_type), intent(in) :: wvl
105 !arrays
106  real(dp),intent(in) :: rprimd(3,3),xred(3,natom)
107 
108 !Local variables-------------------------------
109 !scalars
110 #if defined HAVE_BIGDFT
111  integer :: idata
112  logical,parameter :: wvl_debug=.false.
113  character(len=500) :: message
114 !arrays
115  real(dp),allocatable :: xcart(:,:)
116 #endif
117 
118 ! *********************************************************************
119 
120 #if defined HAVE_BIGDFT
121 !Consistency checks, are all pseudo true GTH pseudo with geometric informations?
122  do idata = 1, psps%npsp, 1
123    if (.not. psps%gth_params%set(idata)) then
124      write(message, '(a,a,a,a,I0,a,a,a)' ) ch10,&
125 &     ' wvl_projectors_set :  consistency checks failed,', ch10, &
126 &     '  no GTH parameters found for type number ', idata, '.', ch10, &
127 &     '  Check your input pseudo files.'
128      MSG_ERROR(message)
129    end if
130    if (.not. psps%gth_params%hasGeometry(idata)) then
131      write(message, '(a,a,a,a,a,a)' ) ch10,&
132 &     ' wvl_projectors_set :  consistency checks failed,', ch10, &
133 &     '  the given GTH parameters has no geometry informations.', ch10, &
134 &     '  Upgrade your input pseudo files to GTH with geometric informatoins.'
135      MSG_ERROR(message)
136    end if
137  end do
138 
139  if (wvl_debug) then
140    call wvl_timing(me,'CrtProjectors ','ON')
141  end if
142 
143 !Store xcart for each atom
144  ABI_ALLOCATE(xcart,(3, natom))
145  call xred2xcart(natom, rprimd, xcart, xred)
146  call createProjectorsArrays(wfs%ks%Lzd%Glr,xcart,wvl%atoms,wfs%ks%orbs,&
147  psps%gth_params%radii_cf,wvl_frmult,wvl_frmult,wvl%h(1),wvl%h(2),&
148  wvl%h(3),.false.,proj%nlpsp,proj%G)
149  write(message, '(a,a,a,a,I0)' ) ch10,&
150 & ' wvl_projectors_set : allocate projectors data,', ch10, &
151 & '  size of the compressed array: ', proj%nlpsp%nprojel
152  call wrtout(std_out,message,'COLL')
153 
154 !Deallocations
155  ABI_DEALLOCATE(xcart)
156 
157  if (wvl_debug) then
158    call wvl_timing(me,'CrtProjectors ','OF')
159  end if
160 
161 #else
162  BIGDFT_NOTENABLED_ERROR()
163  if (.false.) write(std_out,*) natom,me,wvl_frmult,psps%npsp,proj%nlpsp,wfs%ks,wvl%h(1),&
164 & rprimd(1,1),xred(1,1)
165 #endif
166 
167 end subroutine wvl_projectors_set