1!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2! Copyright 2010.  Los Alamos National Security, LLC. This material was    !
3! produced under U.S. Government contract DE-AC52-06NA25396 for Los Alamos !
4! National Laboratory (LANL), which is operated by Los Alamos National     !
5! Security, LLC for the U.S. Department of Energy. The U.S. Government has !
6! rights to use, reproduce, and distribute this software.  NEITHER THE     !
7! GOVERNMENT NOR LOS ALAMOS NATIONAL SECURITY, LLC MAKES ANY WARRANTY,     !
8! EXPRESS OR IMPLIED, OR ASSUMES ANY LIABILITY FOR THE USE OF THIS         !
9! SOFTWARE.  If software is modified to produce derivative works, such     !
10! modified software should be clearly marked, so as not to confuse it      !
11! with the version available from LANL.                                    !
12!                                                                          !
13! Additionally, this program is free software; you can redistribute it     !
14! and/or modify it under the terms of the GNU General Public License as    !
15! published by the Free Software Foundation; version 2.0 of the License.   !
16! Accordingly, this program is distributed in the hope that it will be     !
17! useful, but WITHOUT ANY WARRANTY; without even the implied warranty of   !
18! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General !
19! Public License for more details.                                         !
20!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
21
22SUBROUTINE UNIVTAILCOEF(A)
23
24  USE CONSTANTS_MOD
25  USE MYPRECISION
26
27  IMPLICIT NONE
28
29  INTEGER :: I
30  REAL(LATTEPREC) :: A(14)
31  REAL(LATTEPREC) :: R1, R1SQ, RCUT, RMOD
32  REAL(LATTEPREC) :: SCL_R1, POLYNOM, DPOLY, DDPOLY
33  REAL(LATTEPREC) :: DELTA, DELTA2, DELTA3, DELTA4
34  IF (EXISTERROR) RETURN
35
36  IF ( ABS(A(1)) .LT. 0.000000000001 ) THEN
37
38     A(9) = ZERO
39     A(10) = ZERO
40     A(11) = ZERO
41     A(12) = ZERO
42     A(13) = ZERO
43     A(14) = ZERO
44
45  ELSE
46
47     R1 = A(7)
48     RCUT = A(8)
49
50     R1SQ = R1*R1
51
52     RMOD = R1 - A(6)
53
54     POLYNOM = RMOD*(A(2) + RMOD*(A(3) + RMOD*(A(4) + A(5)*RMOD)))
55
56     SCL_R1 = EXP(POLYNOM)
57
58     !     CALL UNIVSCALE(R1, A, SCL_R1)
59
60     !     SCL_R1 = SCL_R1/A(1)
61
62     DELTA = RCUT - R1
63
64     ! Now we're using a 6th order polynomial: fitted to value, first,
65     ! and second derivatives at R1 and R_cut
66
67     A(9) = SCL_R1
68
69     RMOD = R1 - A(6)
70
71     DPOLY = A(2) + TWO*A(3)*RMOD + THREE*A(4)*RMOD*RMOD + &
72          FOUR*A(5)*RMOD*RMOD*RMOD
73
74     A(10) = DPOLY*SCL_R1
75
76     DDPOLY = TWO*A(3) + SIX*A(4)*RMOD + TWELVE*A(5)*RMOD*RMOD
77
78     A(11) = (DPOLY*DPOLY + DDPOLY)*SCL_R1/TWO
79
80     DELTA2 = DELTA*DELTA
81     DELTA3 = DELTA2*DELTA
82     DELTA4 = DELTA3*DELTA
83
84     A(12) = (MINUSONE/DELTA3)*(THREE*A(11)*DELTA2 + &
85          SIX*A(10)*DELTA + TEN*A(9))
86
87     A(13) = (ONE/DELTA4)*(THREE*A(11)*DELTA2 + &
88          EIGHT*A(10)*DELTA + FIFTEEN*A(9))
89
90     A(14) = (MINUSONE/(TEN*DELTA3)) * &
91          (SIX*A(13)*DELTA2 + THREE*A(12)*DELTA + A(11))
92
93  ENDIF
94
95  RETURN
96
97END SUBROUTINE UNIVTAILCOEF
98