TABLE OF CONTENTS


ABINIT/m_polynomial_term [ Functions ]

[ Top ] [ Functions ]

NAME

 m_polynomial_term

FUNCTION

 Module with the datatype polynomial terms

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

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_polynomial_term
27 
28  use defs_basis
29  use m_errors
30  use m_abicore
31 
32  implicit none
33 
34  public :: polynomial_term_init
35  public :: polynomial_term_free

m_polynomial_term/polynomial_term_free [ Functions ]

[ Top ] [ m_polynomial_term ] [ Functions ]

NAME

 polynomial_term_free

FUNCTION

 Free polynomial_term

INPUTS

 polynomial_term<type(polynomial_term)> =  datatype to free

OUTPUT

 polynomial_term<type(polynomial_term)> =  datatype to free

PARENTS

      m_effective_potential_file,m_fit_polynomial_coeff,m_polynomial_coeff
      m_polynomial_term

CHILDREN

SOURCE

317 subroutine polynomial_term_free(polynomial_term)
318 
319 
320 !This section has been created automatically by the script Abilint (TD).
321 !Do not modify the following lines by hand.
322 #undef ABI_FUNC
323 #define ABI_FUNC 'polynomial_term_free'
324 !End of the abilint section
325 
326  implicit none
327 
328 !Arguments ------------------------------------
329 !scalars
330 !arrays
331  type(polynomial_term_type), intent(inout) :: polynomial_term
332 !Local variables-------------------------------
333 !scalar
334 !arrays
335 
336 ! *************************************************************************
337 
338  polynomial_term%ndisp     = 0
339  polynomial_term%nstrain   = 0
340  polynomial_term%weight    = zero
341 
342  if(allocated(polynomial_term%atindx))then
343    polynomial_term%atindx(:,:) = 0
344    ABI_DEALLOCATE(polynomial_term%atindx)
345  end if
346 
347  if(allocated(polynomial_term%cell))then
348    polynomial_term%cell(:,:,:) = 0
349    ABI_DEALLOCATE(polynomial_term%cell)
350  end if
351 
352  if(allocated(polynomial_term%direction))then
353    polynomial_term%direction(:) = 0
354    ABI_DEALLOCATE(polynomial_term%direction)
355  end if
356 
357  if(allocated(polynomial_term%power_disp))then
358    polynomial_term%power_disp(:) = 0
359    ABI_DEALLOCATE(polynomial_term%power_disp)
360  end if
361 
362   if(allocated(polynomial_term%power_strain))then
363    polynomial_term%power_strain(:) = 0
364    ABI_DEALLOCATE(polynomial_term%power_strain)
365  end if
366 
367   if(allocated(polynomial_term%strain))then
368    polynomial_term%strain(:) = 0
369    ABI_DEALLOCATE(polynomial_term%strain)
370  end if
371 
372 
373 end subroutine polynomial_term_free

m_polynomial_term/polynomial_term_init [ Functions ]

[ Top ] [ m_polynomial_term ] [ Functions ]

NAME

 polynomial_term_init

FUNCTION

 Initialize a polynomial_term_init for given set of displacements/strain

INPUTS

 atindx(2) = Indexes of the atoms a and b in the unit cell
 cell(3,2) = Indexes of the cell of the atom a and b
 direction = direction of the perturbation => 1,2,3 for atomic displacement
 ndisp     = Number of displacement for this terms
 nstrain   = Number of strain for this terms
 power_disp   = power_disp of the displacement 2 (X_z-O_z)^2 or 1 for (X_y-O_y)^1
 power_strain = power_strain of the strain 2 (\eta)^2 or 1 for (\eta)^1
 strain(nstrain) = index of strain, 1 2 3 4 5 or 6
 weight     = Weight of the term
 check      = optional,logical => if TRUE, the term will be check, same displacement/strain
                                          are gathered in the same displacement but with an higher
                                          power_disp. For example:
                                          ((Sr_y-O1_y)^1(Sr_y-O1_y)^1 => (Sr_y-O1_y)^2)
                                 if FALSE, default, do nothing

OUTPUT

 polynomial_term<type(polynomial_term)> = polynomial_term datatype is now initialized

PARENTS

      m_effective_potential_file,m_fit_polynomial_coeff,m_polynomial_coeff

CHILDREN

SOURCE

128 subroutine polynomial_term_init(atindx,cell,direction,ndisp,nstrain,polynomial_term,power_disp,&
129 &                               power_strain,strain,weight,check)
130 
131 
132 !This section has been created automatically by the script Abilint (TD).
133 !Do not modify the following lines by hand.
134 #undef ABI_FUNC
135 #define ABI_FUNC 'polynomial_term_init'
136 !End of the abilint section
137 
138  implicit none
139 
140 !Arguments ------------------------------------
141 !scalars
142  integer, intent(in) :: ndisp,nstrain
143  real(dp),intent(in) :: weight
144  logical,optional,intent(in)  :: check
145 !arrays
146  integer, intent(in) :: atindx(2,ndisp)
147  integer, intent(in) :: cell(3,2,ndisp)
148  integer, intent(in) :: direction(ndisp),power_disp(ndisp)
149  integer, intent(in) :: strain(nstrain),power_strain(nstrain)
150  type(polynomial_term_type), intent(out) :: polynomial_term
151 !Local variables-------------------------------
152 !scalar
153  integer :: idisp1,idisp2,ndisp_tmp,nstrain_tmp
154  logical :: check_in
155 !arrays
156  integer :: power_disp_tmp(ndisp),power_strain_tmp(nstrain)
157  character(500) :: msg
158 
159 ! *************************************************************************
160 
161 !Do some checks
162  if (size(atindx,2) /= ndisp) then
163    write(msg,'(a)')' atindx and ndisp have not the same size'
164    MSG_ERROR(msg)
165  end if
166 
167  if (size(cell,3) /= ndisp) then
168    write(msg,'(a)')' cell and ndisp have not the same size'
169    MSG_ERROR(msg)
170  end if
171 
172  if (size(direction) /= ndisp) then
173    write(msg,'(a)')' direction and ndisp have not the same size'
174    MSG_ERROR(msg)
175  end if
176 
177  if (size(power_disp) /= ndisp) then
178    write(msg,'(a)')' power_disp and ndisp have not the same size'
179    MSG_ERROR(msg)
180  end if
181 
182  if (size(power_strain) /= nstrain) then
183    write(msg,'(a)')' power_strain and nstrain have not the same size'
184    MSG_ERROR(msg)
185  end if
186  
187  if (size(strain) /= nstrain) then
188    write(msg,'(a)')' strain and nstrain have not the same size'
189    MSG_ERROR(msg)
190  end if
191 
192 !First free datatype before init
193  call polynomial_term_free(polynomial_term)
194  check_in = .false.
195 
196 !Copy the power array before check
197  power_disp_tmp(:) = power_disp(:)
198  power_strain_tmp(:) = power_strain(:)
199  
200  if(present(check)) check_in = check
201 
202  if(check_in)then
203 !Check if displacement are identical, in this case
204 !increase the power_disp 
205    do idisp1=1,ndisp
206      do idisp2=idisp1,ndisp
207        if (idisp1/=idisp2.and.&
208 &        atindx(1,idisp1)   == atindx(1,idisp2).and.&
209 &        atindx(2,idisp1)   == atindx(2,idisp2).and.&
210 &        direction(idisp1)  == direction(idisp2).and.&
211 &         all(cell(:,1,idisp1)==cell(:,1,idisp2)).and.&
212 &         all(cell(:,2,idisp1)==cell(:,2,idisp2)).and.&
213 &        power_disp_tmp(idisp2) > 0 )then
214          power_disp_tmp(idisp1) = power_disp_tmp(idisp1) + 1
215          power_disp_tmp(idisp2) = 0
216        end if
217      end do
218    end do
219 
220 ! Count the number of power_disp avec the previous check
221 ! or just remove the power_disp equal to zero
222    ndisp_tmp = 0
223    do idisp1=1,ndisp
224      if(power_disp_tmp(idisp1) > zero) then
225        ndisp_tmp = ndisp_tmp + 1
226      end if
227    end do
228 
229 !Check if strain are identical, in this case
230 !increase the power_strain
231    do idisp1=1,nstrain
232      do idisp2=idisp1,nstrain
233        if (idisp1/=idisp2.and.&
234 &        strain(idisp1)  == strain(idisp2).and.&
235 &        power_strain_tmp(idisp2) > 0 )then
236          power_strain_tmp(idisp1) = power_strain_tmp(idisp1) + 1
237          power_strain_tmp(idisp2) = 0
238        end if
239      end do
240    end do
241 
242 ! Count the number of power_strain avec the previous check
243 ! or just remove the power_strain equal to zero
244    nstrain_tmp = 0
245    do idisp1=1,nstrain
246      if(power_strain_tmp(idisp1) > zero) then
247        nstrain_tmp = nstrain_tmp + 1
248      end if
249    end do
250 
251  else
252    ndisp_tmp   = ndisp
253    nstrain_tmp = nstrain
254  end if!end check
255 
256 !init the values
257  polynomial_term%ndisp    = ndisp_tmp
258  polynomial_term%nstrain  = nstrain_tmp
259  polynomial_term%weight   = weight
260 
261  ABI_ALLOCATE(polynomial_term%atindx,(2,polynomial_term%ndisp))
262  ABI_ALLOCATE(polynomial_term%direction,(polynomial_term%ndisp))
263  ABI_ALLOCATE(polynomial_term%cell,(3,2,polynomial_term%ndisp))
264  ABI_ALLOCATE(polynomial_term%power_disp,(polynomial_term%ndisp))
265  ABI_ALLOCATE(polynomial_term%power_strain,(polynomial_term%nstrain))
266  ABI_ALLOCATE(polynomial_term%strain,(polynomial_term%nstrain))
267 
268 !Transfert displacement 
269  idisp2 = 0
270  do idisp1=1,ndisp
271    if(power_disp_tmp(idisp1) > zero)then
272      idisp2 =  idisp2 + 1
273      polynomial_term%direction(idisp2)  =  direction(idisp1)
274      polynomial_term%power_disp(idisp2) =  power_disp_tmp(idisp1)
275      polynomial_term%atindx(:,idisp2)   =  atindx(:,idisp1) 
276      polynomial_term%cell(:,:,idisp2)   =  cell(:,:,idisp1)
277      polynomial_term%power_disp(idisp2) =  power_disp_tmp(idisp1)
278    end if
279  end do
280 
281 !Transfert strain 
282  idisp2 = 0
283  do idisp1=1,nstrain
284    if(power_strain_tmp(idisp1) > zero)then
285      idisp2 =  idisp2 + 1
286      polynomial_term%power_strain(idisp2) = power_strain_tmp(idisp1)
287      polynomial_term%strain(idisp2) = strain(idisp1)
288    end if   
289  end do
290 
291 end subroutine polynomial_term_init

m_polynomial_term/polynomial_term_type [ Types ]

[ Top ] [ m_polynomial_term ] [ Types ]

NAME

 polynomial_term_type

FUNCTION

 Datatype for a terms  (displacements or strain)
 related to a polynomial coefficient

SOURCE

48  type, public :: polynomial_term_type
49 
50    integer :: ndisp = 0
51 !     Number of displacement for this terms
52 !     1 for (X_y-O_y)^3, 2 for (X_y-O_y)(X_x-O_y)^2...
53 
54    integer :: nstrain = 0
55 !     Number of strain for this terms
56 
57    integer,allocatable :: atindx(:,:)
58 !     atindx(2,ndisp)
59 !     Indexes of the atoms a and b in the unit cell
60 
61    integer,allocatable :: cell(:,:,:)
62 !     cell(3,2,ndisp)
63 !     indexes of the cell of the atom a and b
64 
65    integer,allocatable :: direction(:)
66 !     direction(ndisp)
67 !     direction of the displacement
68 
69    integer,allocatable :: strain(:)
70 !     strain(nstrain)
71 !     strain
72 
73    integer,allocatable :: power_disp(:)
74 !     power_disp(ndisp)
75 !     power of the displacement 2 (X_z-O_z)^2 or 1 for (X_y-O_y)^1
76 
77    integer,allocatable :: power_strain(:)
78 !     power_strain(nstrain)
79 !     power of the strain 2 (\eta_1)^2 or 1 for (\eta_1)^1
80 
81    real(dp) :: weight = zero
82 !     weight of the term
83 
84  end type polynomial_term_type

m_polynomial_term/terms_compare [ Functions ]

[ Top ] [ m_polynomial_term ] [ Functions ]

NAME

  equal

FUNCTION

  Compare two polynomial_term_dot

INPUTS

 t1<type(polynomial_term)> =  datatype of the first term
 t2<type(polynomial_term)> =  datatype of the second term

OUTPUT

 res = logical 

SOURCE

392 pure function terms_compare(t1,t2) result (res)
393 !Arguments ------------------------------------
394 
395 !This section has been created automatically by the script Abilint (TD).
396 !Do not modify the following lines by hand.
397 #undef ABI_FUNC
398 #define ABI_FUNC 'terms_compare'
399 !End of the abilint section
400 
401  implicit none
402 
403 !Arguments ------------------------------------
404   type(polynomial_term_type), intent(in) :: t1,t2
405   logical :: res
406 !local
407 !variable
408   integer :: ia,idisp1,idisp2,mu
409   logical :: found
410 !array
411   integer :: blkval(2,t1%ndisp+t1%nstrain)
412 ! *************************************************************************
413   res = .true.
414   blkval(:,:) = 0
415   if(t1%ndisp==t2%ndisp.and.t1%nstrain==t2%nstrain)then
416 !   Check strain    
417     blkval(:,:) = 0
418     do idisp1=1,t1%nstrain
419       if(blkval(1,t1%ndisp+idisp1)==1)cycle!already found
420       do idisp2=1,t2%nstrain
421         if(blkval(2,t1%ndisp+idisp2)==1)cycle!already found
422         found = .false.
423         if(t1%strain(idisp1) ==  t2%strain(idisp2).and.&
424 &          t1%power_strain(idisp1) == t2%power_strain(idisp2))then
425           found=.true.
426         end if
427         if(found)then
428           blkval(1,t1%ndisp+idisp1) = 1 
429           blkval(2,t1%ndisp+idisp2) = 1
430         end if
431       end do
432     end do
433     if(any(blkval(:,t1%ndisp+1:t1%ndisp+t1%nstrain) == 0))then
434       res = .false.
435       return
436     end if
437 !   Check displacement
438     do idisp1=1,t1%ndisp
439       if(blkval(1,idisp1)==1)cycle!already found
440       do idisp2=1,t2%ndisp
441         if(blkval(2,idisp2)==1)cycle!already found
442         found = .false.
443         if(t1%atindx(1,idisp1)  ==  t2%atindx(1,idisp2).and.&
444 &          t1%atindx(2,idisp1)  ==  t2%atindx(2,idisp2).and.&
445 &          t1%direction(idisp1) ==  t2%direction(idisp2).and.&
446 &          t1%power_disp(idisp1) == t2%power_disp(idisp2))then
447           found=.true.
448           do ia=1,2
449             do mu=1,3
450               if(t1%cell(mu,ia,idisp1) /= t2%cell(mu,ia,idisp2))then
451                 found = .false.
452                 cycle
453               end if
454             end do
455           end do
456           if(found)then
457             blkval(1,idisp1) = 1 
458             blkval(2,idisp2) = 1
459           end if
460         end if
461       end do
462     end do
463     if(any(blkval(:,:)==0))res = .false.
464   else
465     res = .false.
466   end if
467 
468 end function terms_compare