TABLE OF CONTENTS


ABINIT/m_polynomial_conf [ Functions ]

[ Top ] [ Functions ]

NAME

 m_polynomial_conf

FUNCTION

 Module for using a confinement potential
 Container type is defined, and destruction

COPYRIGHT

 Copyright (C) 2010-2018 ABINIT group (AM)
 This file is distributed under the terms of the
 GNU General Public Licence, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_polynomial_conf
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32  use m_xmpi,only   : xmpi_sum
33 
34  implicit none
35 
36  public :: polynomial_conf_init
37  public :: polynomial_conf_free
38 !AM_2017Need to debug this routine
39  public :: polynomial_conf_evaluate

m_polynomial_conf/polynomial_conf_evaluate [ Functions ]

[ Top ] [ m_polynomial_conf ] [ Functions ]

NAME

  polynomial_conf_evaluate

FUNCTION

  This fonction evaluate the energy, (soon forces and stresses) with  the confinement potential

INPUTS

   disp(3,natom_sc) = atomic displacments of a specific patern wrt to reference structure
   disp_ref(natom_uc) = Cutoff array for the atomic displacement
   factor_disp = Factor to appy to the polynomial term of the confinement (displacement)
   factor_strain = Factor to appy to the polynomial term of the confinement (strain)
   strain(6) =   strain of a specific structure wrt to reference
   strain_ref(6) = Cutoff array for the strain
   power_disp = Power of the polynome related to the displacement
   power_strain = Power of the polynome related to the strain
   natom_sc = number of atoms in the supercell
   natom_uc = number of atoms in the unit cell
   ncell   = total number of cell to treat by this cpu
   cells(ncell) = index of  the cells into the supercell (1,2,3,4,5)
   index_cells(3,ncell) = indexes  of the cells into supercell (-1 -1 -1 ,...,1 1 1)
   comm=MPI communicator

OUTPUT

   energy = contribution to the ifc to the energy
   fcart(3,natom) = contribution to the ifc to the forces
   strten(6) = contribution to the stress tensor

PARENTS

      m_effective_potential

CHILDREN

      xmpi_sum

SOURCE

261 subroutine polynomial_conf_evaluate(disp,disp_ref,energy,factor_disp,factor_strain,fcart,&
262 &                                   strain,strain_ref,strten,power_disp,power_strain,cells,&
263 &                                   natom_sc,natom_uc,ncell,index_cells,comm)
264 
265 
266 !This section has been created automatically by the script Abilint (TD).
267 !Do not modify the following lines by hand.
268 #undef ABI_FUNC
269 #define ABI_FUNC 'polynomial_conf_evaluate'
270 !End of the abilint section
271 
272  implicit none
273 
274 !Arguments -------------------------------
275 ! scalars
276   real(dp),intent(out) :: energy
277   integer,intent(in) :: natom_uc,natom_sc,ncell
278   integer,intent(in) :: power_disp,power_strain
279   integer,intent(in) :: comm
280   real(dp),intent(in) :: factor_disp,factor_strain
281 ! array
282   integer,intent(in) ::  cells(ncell),index_cells(ncell,3)
283   real(dp),intent(in) :: disp(3,natom_sc),strain(6)
284   real(dp),intent(out) :: fcart(3,natom_sc),strten(6)
285   real(dp),intent(in) :: disp_ref(natom_uc),strain_ref(6)
286 !Local variables-------------------------------
287 ! scalar
288   integer :: ia,icell,ierr,ii,kk
289   integer :: mu
290   real(dp):: diff,diff_tmp
291 ! array
292 
293 ! *************************************************************************
294 
295 ! Initialisation of variables
296   energy   = zero
297   fcart(:,:) = zero
298   strten(:) = zero
299 
300   write(std_out,*) factor_strain,index_cells,power_strain,strain,strain_ref
301   do icell = 1,ncell
302     ii = (cells(icell)-1)*natom_uc
303     do ia = 1, natom_uc
304       kk = ii + ia
305       diff_tmp = zero
306       do mu=1,3
307         diff_tmp = diff_tmp + disp(mu,kk)**2
308       end do
309       diff_tmp = diff_tmp**0.5
310 !     Compute diff between ref and curent displacement
311       diff = diff_tmp - disp_ref(ia)
312 !     Accumule energy
313       energy =  energy + (sign(half, diff)+half)*(factor_disp*((diff_tmp/disp_ref(ia))**power_disp))
314     end do
315   end do
316 
317 ! MPI_SUM
318   call xmpi_sum(energy, comm, ierr)
319 
320 end subroutine polynomial_conf_evaluate

m_polynomial_conf/polynomial_conf_free [ Functions ]

[ Top ] [ m_polynomial_conf ] [ Functions ]

NAME

 polynomial_conf_free

FUNCTION

 Free polynomial_conf

INPUTS

 polynomial_conf <type(polynomial_conf)> = polynomial_conf datatype to be free

OUTPUT

 polynomial_conf <type(polynomial_conf)> = polynomial_conf datatype to be free

PARENTS

      m_effective_potential,m_polynomial_conf

CHILDREN

      xmpi_sum

SOURCE

190 subroutine polynomial_conf_free(polynomial_conf)
191 
192 
193 !This section has been created automatically by the script Abilint (TD).
194 !Do not modify the following lines by hand.
195 #undef ABI_FUNC
196 #define ABI_FUNC 'polynomial_conf_free'
197 !End of the abilint section
198 
199  implicit none
200 
201 !Arguments ------------------------------------
202 !scalars
203 !arrays
204  type(polynomial_conf_type), intent(inout) :: polynomial_conf
205 !Local variables-------------------------------
206 !scalar
207 !arrays
208 
209 ! *************************************************************************
210 
211  if(allocated(polynomial_conf%cutoff_disp))then
212    ABI_DEALLOCATE(polynomial_conf%cutoff_disp)
213  end if
214 
215  polynomial_conf%power_disp    = 0
216  polynomial_conf%power_strain  = 0
217  polynomial_conf%factor_disp   = zero
218  polynomial_conf%factor_strain = zero
219  polynomial_conf%cutoff_strain = zero
220  polynomial_conf%need_confinement = .FALSE.
221 
222 end subroutine polynomial_conf_free

m_polynomial_conf/polynomial_conf_init [ Functions ]

[ Top ] [ m_polynomial_conf ] [ Functions ]

NAME

 polynomial_conf_init

FUNCTION

 Initialize polynomial_conf_init

INPUTS

 cutoff_disp(6) = Cutoff array for the strain
 cutoff_strain(ndisp) = Cutoff array for the atomic displacement
 factor_disp = Factor to appy to the polynomial term of the confinement (displacement)
 factor_strain = Factor to appy to the polynomial term of the confinement (strain)
 ndisp = Number of displacement (atoms) for the cut off
 power_disp = Power of the polynome related to the displacement
 power_strain = Power of the polynome related to the strain
 need_confinement = optional,Logical related to the necessity of the confinement

OUTPUT

 polynomial_conf <type(polynomial_conf)> = datatype with the information for the confinement
                                           polynomial

PARENTS

      m_effective_potential

CHILDREN

      xmpi_sum

SOURCE

115 subroutine polynomial_conf_init(cutoff_disp,cutoff_strain,factor_disp,factor_strain,ndisp,&
116 &                               polynomial_conf,power_disp,power_strain,need_confinement)
117 
118 
119 !This section has been created automatically by the script Abilint (TD).
120 !Do not modify the following lines by hand.
121 #undef ABI_FUNC
122 #define ABI_FUNC 'polynomial_conf_init'
123 !End of the abilint section
124 
125  implicit none
126 
127 !Arguments ------------------------------------
128 !scalars
129  integer, intent(in) :: ndisp,power_disp,power_strain
130  real(dp),intent(in) :: factor_disp,factor_strain
131  logical,optional,intent(in)  :: need_confinement
132 !arrays
133  real(dp),intent(in) :: cutoff_disp(ndisp),cutoff_strain(6)
134  type(polynomial_conf_type),intent(inout) :: polynomial_conf
135 !Local variables-------------------------------
136 !scalar
137 !arrays
138  character(len=500) :: msg
139 
140 ! *************************************************************************
141 
142 !Checks
143  if (ndisp <= 0) then
144    write(msg,'(a,a)')' ndisp can not be inferior or equal to zero'
145    MSG_ERROR(msg)
146  end if
147 
148 !First free the type
149  call  polynomial_conf_free(polynomial_conf)
150 
151  polynomial_conf%power_disp    = power_disp
152  polynomial_conf%power_strain  = power_strain
153  polynomial_conf%factor_disp   = factor_disp
154  polynomial_conf%factor_strain = factor_strain
155  polynomial_conf%need_confinement = .FALSE.
156 
157  polynomial_conf%ndisp   = ndisp
158  ABI_ALLOCATE(polynomial_conf%cutoff_disp,(polynomial_conf%ndisp))
159  polynomial_conf%cutoff_disp(:) = cutoff_disp(:)
160 
161  polynomial_conf%cutoff_strain = cutoff_strain(:)
162  if (present(need_confinement)) polynomial_conf%need_confinement = need_confinement
163 
164 end subroutine polynomial_conf_init

m_polynomial_conf/polynomial_conf_type [ Types ]

[ Top ] [ m_polynomial_conf ] [ Types ]

NAME

 polynomial_conf_type

FUNCTION

 datatype for specific confinement potential

SOURCE

51  type, public :: polynomial_conf_type
52 
53    integer :: ndisp = 0
54 !  Number of displacement (atoms) for the cut off
55 
56    integer :: power_disp = 0
57 !  Power of the polynome related to the displacement
58 
59    integer :: power_strain = 0
60 !  Power of the polynome related to the strain
61 
62    real(dp):: factor_disp = 0
63 !  Factor to appy to the polynomial term of the confinement (displacement)
64 
65    real(dp):: factor_strain = 0
66 !  Factor to appy to the polynomial term of the confinement (strain)
67 
68    real(dp):: cutoff_strain(6)
69 !  Cutoff array for the strain
70 
71    real(dp),allocatable :: cutoff_disp(:)
72 !  Cutoff array for the atomic displacement
73 
74    logical :: need_confinement =.FALSE.
75 !  Logical related to the necessity of the confinement
76 
77  end type polynomial_conf_type