TABLE OF CONTENTS


ABINIT/glue_lotf [ Modules ]

[ Top ] [ Modules ]

NAME

 glue_lotf

FUNCTION

  Contains the GLUE procedure and parameters for Lotf

COPYRIGHT

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

NOTES

  Parameters for the glue potential
  ref: Ercolessi et al, Phil Mag A 58(1), 213 (1988)

PARENTS

CHILDREN

SOURCE

25 #if defined HAVE_CONFIG_H
26 #include "config.h"
27 #endif
28 
29 #include "abi_common.h"
30 
31 module glue_lotf
32 
33  use defs_basis
34  use m_errors
35 
36  implicit none
37  private
38 
39  !--Parameters for the glue potential (Ercolessi et al, Phil Mag A 58(1), 213 (1988))
40 
41  !--Function phi_r
42  real(dp),public :: dphi,rcphi
43  real(dp) :: a0I,a1I,a2I,a3I,a4I
44  real(dp) :: a0II,a1II,a2II,a3II,a4II,a5II,a6II
45 
46  !--Function rho_r
47  real(dp),public :: rmrho
48  real(dp) :: rbrho
49  real(dp) :: b0I,b1I,b2I,b3I
50  real(dp) :: b0II,b1II,b2II,b3II,b0III,b1III,b2III,b3III
51 
52  !--Function U_n
53  real(dp) :: n0U,nsU
54  real(dp) :: c0I,c1I,c2I,c3I,c4I
55  real(dp) :: c0II,c1II,c2II,c3II,c4II,c0III,c1III,c2III,c3III
56 
57  !--Auxiliary variables to change glue parameters
58  real(dp),parameter :: ddphi = 0.028779246_dp
59  real(dp),parameter :: scalefactor_phi = one
60 
61 
62  public ::             &
63    glue_init,          &
64    glue_pair_devs,     &
65    glue_pair,          &
66    calc_coord,         &
67    calc_rhop,          &
68    calc_coord_new_d,   &
69    rho_devs,           &
70    rhop_value,         &
71    eval_U_n,           &
72    eval_Upp_n
73 
74 contains !===========================================================

eval_lotf/eval_Upp_n [ Functions ]

[ Top ] [ eval_lotf ] [ Functions ]

NAME

 eval_Upp_n

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

700  subroutine eval_Upp_n(coordatom_i,up_dum,upp_dum)
701 
702 
703 !This section has been created automatically by the script Abilint (TD).
704 !Do not modify the following lines by hand.
705 #undef ABI_FUNC
706 #define ABI_FUNC 'eval_Upp_n'
707 !End of the abilint section
708 
709   implicit none
710 
711   !Arguments ------------------------
712   real(dp),intent(in) :: coordatom_i
713   real(dp),intent(out) :: up_dum,upp_dum
714   !Local ---------------------------
715   real(dp) :: cns,cns2,cns3,cn0,cn02,cn03
716   character(len=500) :: msg
717 
718 ! *************************************************************************
719 
720   !--Evaluate U and Up for atom i
721   up_dum = zero
722   upp_dum = zero
723 
724   if (coordatom_i < zero) then
725 
726     write(msg,'(3a,i8,2a)')&
727       &    'ERROR: negative coordination!!! ',ch10,&
728       &    'coord = ', coordatom_i,ch10,&
729       &    'The coordination cannot be negative!'
730     MSG_ERROR(msg)
731 
732   elseif (coordatom_i < nsU) then
733 
734     cns = coordatom_i - nsU
735     cns2 = cns * cns
736     cns3 = cns2 * cns
737 
738     up_dum = 4.d0*c4I*cns3 + 3.d0*c3I*cns2 + 2.d0*c2I*cns + c1I
739     upp_dum = 12.d0*c4I*cns2 + 6.d0*c3I*cns + 2.d0*c2I
740 
741   elseif (coordatom_i < n0U) then
742 
743     cn0 = coordatom_i - n0U
744     cn02 = cn0 * cn0
745     cn03 = cn02 * cn0
746 
747     up_dum = 4.d0*c4II*cn03 + 3.d0*c3II*cn02 + 2.d0*c2II*cn0 + c1II
748     upp_dum = 12.d0*c4II*cn02 + 6.d0*c3II*cn0 + 2.d0*c2II
749 
750   elseif (coordatom_i >= n0U) then
751 
752     cn0 = coordatom_i - n0U
753     cn02 = cn0 * cn0
754 
755     up_dum = 3.d0*c3III*cn02 + 2.d0*c2III*cn0 + c1III
756     upp_dum = 6.d0*c3III*cn0 + 2.d0*c2III
757 
758   endif
759  end subroutine eval_Upp_n
760 
761 end module glue_lotf

glue_lotf/calc_coord [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 calc_coord

FUNCTION

INPUTS

PARENTS

      m_eval_lotf

CHILDREN

SOURCE

350  subroutine calc_coord(r_au,coordatom_dum)
351 
352 
353 !This section has been created automatically by the script Abilint (TD).
354 !Do not modify the following lines by hand.
355 #undef ABI_FUNC
356 #define ABI_FUNC 'calc_coord'
357 !End of the abilint section
358 
359   implicit none
360 
361   !Arguments ------------------------
362   real(dp) :: r_au,coordatom_dum
363   !Local ---------------------------
364   real(dp) :: rst,r,rho_neigh
365   ! *********************************************************************
366 !    --Adds the density of a single given neighbour
367 
368      rst = sqrt(r_au)
369      rho_neigh = zero
370 
371      if (rst < dphi) then
372        r = rst - dphi
373        rho_neigh = b3I*r*r*r + b2I*r*r + b1I*r + b0I
374      elseif (rst < rbrho) then
375        r = rst - dphi
376        rho_neigh = b3II*r*r*r + b2II*r*r + b1II*r + b0II
377      elseif (rst < rmrho) then
378        r = rst - rmrho
379        rho_neigh = b3III*r*r*r + b2III*r*r + b1III*r + b0III
380      end if
381 
382      coordatom_dum = coordatom_dum + rho_neigh
383 
384    end subroutine calc_coord

glue_lotf/calc_coord_new_d [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 calc_coord_new_d

FUNCTION

INPUTS

PARENTS

      m_eval_lotf

CHILDREN

SOURCE

456  subroutine calc_coord_new_d(r_au,alpha_d,coordatom_dum)
457 
458 
459 !This section has been created automatically by the script Abilint (TD).
460 !Do not modify the following lines by hand.
461 #undef ABI_FUNC
462 #define ABI_FUNC 'calc_coord_new_d'
463 !End of the abilint section
464 
465   implicit none
466 
467   !Arguments ------------------------
468   real(dp),intent(in) :: r_au,alpha_d
469   real(dp),intent(out) ::coordatom_dum
470   !Local ---------------------------
471   real(dp) :: rst,r,rho_neigh
472   ! *********************************************************************
473 
474 !    --Adds the density of a single given neighbour
475      rst = sqrt(r_au)
476 
477      rho_neigh = zero
478 
479      if (rst < alpha_d) then
480        r = rst - alpha_d
481        rho_neigh = b3I*r*r*r + b2I*r*r + b1I*r + b0I
482      elseif (rst < rbrho) then
483        r = rst - alpha_d
484        rho_neigh = b3II*r*r*r + b2II*r*r + b1II*r + b0II
485 
486 !      --Restriction on rmrho to keep rho and rhop continuous
487      elseif (rst < (rmrho+alpha_d-dphi)) then
488        r = rst - (rmrho+alpha_d-dphi)
489        rho_neigh = b3III*r*r*r + b2III*r*r + b1III*r + b0III
490      end if
491 
492      coordatom_dum = coordatom_dum + rho_neigh
493    end subroutine calc_coord_new_d

glue_lotf/calc_rhop [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 calc_rhop

FUNCTION

INPUTS

PARENTS

      m_eval_lotf

CHILDREN

SOURCE

404      subroutine calc_rhop(r_st,rhop_dum)
405 
406 
407 !This section has been created automatically by the script Abilint (TD).
408 !Do not modify the following lines by hand.
409 #undef ABI_FUNC
410 #define ABI_FUNC 'calc_rhop'
411 !End of the abilint section
412 
413      implicit none
414 
415 !    Arguments ------------------------
416      real(dp),intent(in) :: r_st
417      real(dp),intent(out) :: rhop_dum
418 !    Local ---------------------------
419      real(dp) :: r,r2
420 !    *********************************************************************
421 
422      rhop_dum = zero
423 
424      if (r_st < dphi) then
425        r = r_st - dphi
426        r2 = r*r
427        rhop_dum = (3.d0*b3I*r2 + 2.d0*b2I*r + b1I)
428      elseif (r_st < rbrho) then
429        r = r_st - dphi
430        r2 = r*r
431        rhop_dum = (3.d0*b3II*r2 + 2.d0*b2II*r + b1II)
432      elseif (r_st < rmrho) then
433        r = r_st - rmrho
434        r2 = r*r
435        rhop_dum = (3.d0*b3III*r2 + 2.d0*b2III*r + b1III)
436      end if
437    end subroutine calc_rhop

glue_lotf/eval_U_n [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 eval_U_n

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

622  subroutine eval_U_n(coordatom_i,epot_dum2,up_dum)
623 
624 
625 !This section has been created automatically by the script Abilint (TD).
626 !Do not modify the following lines by hand.
627 #undef ABI_FUNC
628 #define ABI_FUNC 'eval_U_n'
629 !End of the abilint section
630 
631   implicit none
632 
633   !Arguments ------------------------
634   real(dp),intent(in) :: coordatom_i
635   real(dp),intent(out) :: epot_dum2,up_dum
636   !Local ---------------------------
637   real(dp) :: cns,cns2,cns3,cns4,cn0,cn02,cn03,cn04
638   character(len=500) :: msg
639 
640   !--Evaluate U and Up for atom i
641   epot_dum2 = zero
642   up_dum = zero
643 
644   if (coordatom_i < zero) then
645 
646     write(msg,'(3a,i8,2a)')&
647       &    'ERROR: negative coordination!!! ',ch10,&
648       &    'coord = ', coordatom_i,ch10,&
649       &    'The coordination cannot be negative!'
650     MSG_ERROR(msg)
651 
652   elseif (coordatom_i < nsU) then
653 
654     cns = coordatom_i - nsU
655     cns2 = cns * cns
656     cns3 = cns2 * cns
657     cns4 = cns3 * cns
658 
659     epot_dum2 = c4I*cns4 + c3I*cns3 + c2I*cns2 + c1I*cns + c0I
660     up_dum = 4.d0*c4I*cns3 + 3.d0*c3I*cns2 + 2.d0*c2I*cns + c1I
661 
662   elseif (coordatom_i < n0U) then
663 
664     cn0 = coordatom_i - n0U
665     cn02 = cn0 * cn0
666     cn03 = cn02 * cn0
667     cn04 = cn03 * cn0
668 
669     epot_dum2 = c4II*cn04 + c3II*cn03 + c2II*cn02 + c1II*cn0 + c0II
670     up_dum = 4.d0*c4II*cn03 + 3.d0*c3II*cn02 + 2.d0*c2II*cn0 + c1II
671 
672   elseif (coordatom_i >= n0U) then
673 
674     cn0 = coordatom_i - n0U
675     cn02 = cn0 * cn0
676     cn03 = cn02 * cn0
677 
678     epot_dum2 = c3III*cn03 + c2III*cn02 + c1III*cn0 + c0III
679     up_dum = 3.d0*c3III*cn02 + 2.d0*c2III*cn0 + c1III
680 
681   endif
682  end subroutine eval_U_n

glue_lotf/glue_init [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 glue_init

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

 91  subroutine glue_init()
 92 
 93   !--Function phi_r
 94 ! *************************************************************************
 95 
 96 !This section has been created automatically by the script Abilint (TD).
 97 !Do not modify the following lines by hand.
 98 #undef ABI_FUNC
 99 #define ABI_FUNC 'glue_init'
100 !End of the abilint section
101 
102   dphi =  (0.2878207442141723d+01 + ddphi) / 0.529177d0
103   rcphi = (0.3700000000000000d+01 + ddphi) / 0.529177d0
104 
105   a0I = -0.8000000000000000d-01 / (13.6057d0 * 2.d0) *&
106     scalefactor_phi
107   a1I =  0.0000000000000000d+00 / (13.6057d0 * 2.d0) * (0.529177d0) *&
108     scalefactor_phi
109   a2I =  0.7619231375231362d+00 / (13.6057d0 * 2.d0) * (0.529177d0)**2 *&
110     scalefactor_phi
111   a3I = -0.8333333333333333d+00 / (13.6057d0 * 2.d0) * (0.529177d0)**3 *&
112     scalefactor_phi
113   a4I = -0.1211483464993159d+00 / (13.6057d0 * 2.d0) * (0.529177d0)**4 *&
114     scalefactor_phi
115 
116   !        a0II = -0.8000000000000000d-01 / (13.6057d0 * 2.d0)
117   a0II = a0I
118   !        a1II =  0.0000000000000000d+00 / (13.6057d0 * 2.d0) * (0.529177d0)
119   a1II = a1I
120   !        a2II =  0.7619231375231362d+00 / (13.6057d0 * 2.d0) * (0.529177d0)**2
121   a2II = a2I
122   a3II = -0.8333333333333333d+00 / (13.6057d0 * 2.d0) * (0.529177d0)**3 *&
123     scalefactor_phi
124   a4II = -0.1096009851140349d+01 / (13.6057d0 * 2.d0) * (0.529177d0)**4 *&
125     scalefactor_phi
126   a5II =  0.2158417178555998d+01 / (13.6057d0 * 2.d0) * (0.529177d0)**5 *&
127     scalefactor_phi
128   a6II = -0.9128915709636862d+00 / (13.6057d0 * 2.d0) * (0.529177d0)**6 *&
129     scalefactor_phi
130 
131   !--Function rho_r
132   rbrho = 0.3500000000000000d+01 / 0.529177d0
133   rmrho = (0.3900000000000000d+01 + ddphi) / 0.529177d0
134 
135   b0I =  0.1000000000000000d+01
136   b1I = -0.6800000000000000d+00 * (0.529177d0)
137   b2I =  0.7500000000000000d+00 * (0.529177d0)**2
138   b3I = -0.1333333333333333d+01 * (0.529177d0)**3
139 
140   !        b0II =  0.1000000000000000d+01
141   b0II = b0I
142   !        b1II = -0.6800000000000000d+00 * (0.529177d0)
143   b1II = b1I
144   !        b2II =  0.7500000000000000d+00 * (0.529177d0)**2
145   b2II = b2I
146   b3II = -0.1527241171296038d+01 * (0.529177d0)**3
147 
148   b0III =  0.0000000000000000d+00
149   b1III =  0.0000000000000000d+00 * (0.529177d0)
150   b2III =  0.5578188675490974d+01 * (0.529177d0)**2
151   b3III =  0.6132971688727435d+01 * (0.529177d0)**3
152 
153   !--Function U_n
154   n0U =  0.1200000000000000d+02
155   nsU =  0.9358157767784574d+01
156 
157   c0I = -0.2793388616771698d+01 / (13.6057d0 * 2.d0) * scalefactor_phi
158   c1I = -0.3419999999999999d+00 / (13.6057d0 * 2.d0) * scalefactor_phi
159   c2I =  0.3902327808424106d-01 / (13.6057d0 * 2.d0) * scalefactor_phi
160   c3I =  0.7558829951858879d-02 / (13.6057d0 * 2.d0) * scalefactor_phi
161   c4I =  0.3090472511796849d-03 / (13.6057d0 * 2.d0) * scalefactor_phi
162 
163   c0II = -0.3300000000000000d+01 / (13.6057d0 * 2.d0) * scalefactor_phi
164   c1II =  0.0000000000000000d+00 / (13.6057d0 * 2.d0) * scalefactor_phi
165   c2II =  0.8618226772941980d-01 / (13.6057d0 * 2.d0) * scalefactor_phi
166   c3II =  0.4341701445034724d-02 / (13.6057d0 * 2.d0) * scalefactor_phi
167   c4II = -0.3044398779375916d-03 / (13.6057d0 * 2.d0) * scalefactor_phi
168 
169   !        c0III = -0.3300000000000000d+01 / (13.6057d0 * 2.d0)
170   c0III = c0II
171   !        c1III =  0.0000000000000000d+00 / (13.6057d0 * 2.d0)
172   c1III = c1II
173   !        c2III =  0.8618226772941980d-01 / (13.6057d0 * 2.d0)
174   c2III = c2II
175   c3III =  0.4325981467602070d-02 / (13.6057d0 * 2.d0) * scalefactor_phi
176 
177  end subroutine glue_init

glue_lotf/glue_pair [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 glue_pair

FUNCTION

INPUTS

PARENTS

      m_eval_lotf

CHILDREN

SOURCE

282  subroutine glue_pair(RD,r_au,epot_2,fdum)
283 
284 
285 !This section has been created automatically by the script Abilint (TD).
286 !Do not modify the following lines by hand.
287 #undef ABI_FUNC
288 #define ABI_FUNC 'glue_pair'
289 !End of the abilint section
290 
291   implicit none
292 
293   !Arguments ------------------------
294   real(dp) ::  RD(3), fdum(3), r_au, epot_2
295   !Local ---------------------------
296   real(dp) ::  rst, r, r2, r3, r4, r5, r6
297   integer ::  k
298   ! *********************************************************************
299 !    --The calculation for a given pair
300      rst = sqrt(r_au)
301      r = rst - dphi
302      r2 = r*r
303      r3 = r2*r
304      r4 = r3*r
305 
306 !    --Energy
307      epot_2 = zero
308      if (rst < dphi) then
309        epot_2 = a4I*r4 + a3I*r3 + a2I*r2 + a1I*r + a0I
310      elseif (rst < rcphi) then
311        r5 = r4*r
312        r6 = r5*r
313        epot_2 = a6II*r6 + a5II*r5 + a4II*r4 + a3II*r3 + a2II*r2 +&
314        a1II*r + a0II
315      end if
316 
317 !    --Forces (without the minus sign, because we use rji)
318      fdum(:) = zero
319 
320      do k =1,3
321        if (rst < dphi) then
322          fdum(k) =  ( 4.d0*a4I*r3 + 3.d0*a3I*r2 + 2.d0*a2I*r + a1I ) * &
323          RD(k)/rst
324        elseif (rst < rcphi) then
325          fdum(k) =  ( 6.d0*a6II*r5 + 5.d0*a5II*r4 + 4.d0*a4II*r3 + &
326          3.d0*a3II*r2 + 2.d0*a2II*r + a1II ) * &
327          RD(k)/rst
328        end if
329      end do
330 
331    end subroutine glue_pair

glue_lotf/glue_pair_devs [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 glue_pair_devs

FUNCTION

INPUTS

PARENTS

      m_eval_lotf

CHILDREN

SOURCE

195  subroutine glue_pair_devs(alpha_dum,RD,r_au,epot_2,fdum,dfdum)
196 
197 
198 !This section has been created automatically by the script Abilint (TD).
199 !Do not modify the following lines by hand.
200 #undef ABI_FUNC
201 #define ABI_FUNC 'glue_pair_devs'
202 !End of the abilint section
203 
204   implicit none
205 
206   !Arguments ------------------------
207   real(dp),intent(in) ::  alpha_dum(3),RD(3)
208   real(dp),intent(in) ::  r_au
209   real(dp),intent(out) ::  epot_2
210   real(dp),intent(out) ::  fdum(3),dfdum(3,2)
211   !Local ---------------------------
212   real(dp) :: rst, r, r2, r3, r4, r5, r6
213   real(dp) :: drcphi, dphi2, rcphi2, scale_factor
214   ! *********************************************************************
215 
216 !    --The calculation for a given pair
217      drcphi = rcphi - dphi
218      dphi2 = alpha_dum(1)
219 
220 !    --Rewritten
221      rcphi2 = dphi2 + drcphi
222 
223      scale_factor = alpha_dum(2)
224 
225      rst = sqrt(r_au)
226 
227      r = rst - dphi2
228      r2 = r*r
229      r3 = r2*r
230      r4 = r3*r
231 
232 !    --Energy
233      epot_2 = zero
234 
235      if (rst < dphi2) then
236        epot_2 = scale_factor * (a4I*r4 + a3I*r3 + a2I*r2 + a1I*r + a0I)
237      elseif (rst < rcphi2) then
238        r5 = r4*r
239        r6 = r5*r
240        epot_2 = scale_factor * (a6II*r6 + a5II*r5 + a4II*r4 + a3II*r3 + a2II*r2 +&
241        a1II*r + a0II)
242      end if
243 
244 !    --Forces (without the minus sign, because we use rji) and derivatives w.r.t. parameter d
245      fdum(:) = zero
246      dfdum(:,:) = zero
247 
248      if (rst < dphi2) then
249        fdum(:) = scale_factor * (4.d0*a4I*r3 + 3.d0*a3I*r2 + 2.d0*a2I*r + a1I) * RD(:)/rst
250        dfdum(:,1) = -scale_factor * (12.d0*a4I*r2 + 6.d0*a3I*r + 2.d0*a2I) * RD(:)/rst
251        dfdum(:,2) = fdum(:) / scale_factor
252 
253      elseif (rst < rcphi2) then
254        fdum(:) = (scale_factor) * &
255 &       (6.d0*a6II*r5 + 5.d0*a5II*r4 + 4.d0*a4II*r3 + 3.d0*a3II*r2 + 2.d0*a2II*r + a1II) * &
256 &       RD(:)/rst
257 
258        dfdum(:,1) = -scale_factor * &
259 &       (30.d0*a6II*r4 + 20.d0*a5II*r3 + 12.d0*a4II*r2 + 6.d0*a3II*r + 2.d0*a2II) * RD(:)/rst
260 
261        dfdum(:,2) = fdum(:) / scale_factor
262 
263      end if
264    end subroutine glue_pair_devs

glue_lotf/rho_devs [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 rho_devs

FUNCTION

INPUTS

PARENTS

      m_eval_lotf

CHILDREN

SOURCE

511      subroutine rho_devs(r_au,alpha_d,rho_neigh_d,rho_neigh_p,rho_neigh_pd)
512 
513 
514 !This section has been created automatically by the script Abilint (TD).
515 !Do not modify the following lines by hand.
516 #undef ABI_FUNC
517 #define ABI_FUNC 'rho_devs'
518 !End of the abilint section
519 
520      real(dp),intent(in) :: r_au,alpha_d
521      real(dp),intent(out) :: rho_neigh_d,rho_neigh_p,rho_neigh_pd
522 !    Local ---------------------------
523      real(dp) :: rst,r
524 !    *********************************************************************
525 
526 !    --Adds the density of a single given neighbour
527      rst = sqrt(r_au)
528 
529      rho_neigh_p = zero
530      rho_neigh_d = zero
531      rho_neigh_pd = zero
532 
533      if (rst < alpha_d) then
534        r = rst - alpha_d
535        rho_neigh_d = - ( 3.d0*b3I*r*r + 2.d0*b2I*r + b1I )
536        rho_neigh_p = - rho_neigh_d
537        rho_neigh_pd = - ( 6.d0*b3I*r + 2.d0*b2I )
538 
539      elseif (rst < rbrho) then
540        r = rst - alpha_d
541        rho_neigh_d = - ( 3.d0*b3II*r*r + 2.d0*b2II*r + b1II )
542        rho_neigh_p = - rho_neigh_d
543        rho_neigh_pd = - ( 6.d0*b3II*r + 2.d0*b2II )
544 
545 !      --Restriction on rmrho to keep rho and rhop continuous
546      elseif (rst < (rmrho+alpha_d-dphi)) then
547        r = rst - (rmrho+alpha_d-dphi)
548        rho_neigh_d = zero
549        rho_neigh_p = 3.d0*b3III*r*r + 2.d0*b2III*r + b1III
550        rho_neigh_pd = zero
551 
552      end if
553    end subroutine rho_devs

glue_lotf/rhop_value [ Functions ]

[ Top ] [ glue_lotf ] [ Functions ]

NAME

 rhop_value

FUNCTION

INPUTS

PARENTS

      m_eval_lotf

CHILDREN

SOURCE

572  subroutine rhop_value(rst,alpha_d,rhop_dum)
573 
574 
575 !This section has been created automatically by the script Abilint (TD).
576 !Do not modify the following lines by hand.
577 #undef ABI_FUNC
578 #define ABI_FUNC 'rhop_value'
579 !End of the abilint section
580 
581   implicit none
582 
583   !Arguments ------------------------
584   real(dp),intent(in) :: rst,alpha_d
585   real(dp),intent(out) :: rhop_dum
586   !Local ---------------------------
587   real(dp) :: r
588   ! *********************************************************************
589 
590 !    --Adds the density of a single given neighbour
591      rhop_dum = zero
592 
593      if (rst < alpha_d) then
594        r = rst - alpha_d
595        rhop_dum = 3.d0*b3I*r*r + 2.d0*b2I*r + b1I
596      elseif (rst < rbrho) then
597        r = rst - alpha_d
598        rhop_dum = 3.d0*b3II*r*r + 2.d0*b2II*r + b1II
599 !      --Restriction on rmrho to keep rho and rhop continuous
600      elseif (rst < (rmrho+alpha_d-dphi)) then
601        r = rst - (rmrho+alpha_d-dphi)
602        rhop_dum = 3.d0*b3III*r*r + 2.d0*b2III*r + b1III
603      end if
604    end subroutine rhop_value