2 ! { dg-options "-O2 -ffast-math" }
5 INTEGER, PARAMETER :: dp
=8
9 SUBROUTINE b97_lsd_eval(rho_set
,deriv_set
,grad_deriv
,b97_params
)
10 INTEGER, INTENT(in
) :: grad_deriv
11 INTEGER :: handle
, npoints
, param
, stat
13 REAL(kind
=dp
) :: epsilon_drho
, epsilon_rho
, &
15 REAL(kind
=dp
), DIMENSION(:, :, :), POINTER :: dummy
, e_0
, e_ndra
, &
16 e_ndra_ndra
, e_ndra_ndrb
, e_ndra_ra
, e_ndra_rb
, e_ndrb
, e_ndrb_ndrb
, &
17 e_ndrb_ra
, e_ndrb_rb
, e_ra
, e_ra_ra
, e_ra_rb
, e_rb
, e_rb_rb
, &
18 norm_drhoa
, norm_drhob
, rhoa
, rhob
19 IF (.NOT
. failure
) THEN
21 rhoa
=rhoa
, rhob
=rhob
, norm_drhoa
=norm_drhoa
,&
22 norm_drhob
=norm_drhob
, e_0
=e_0
, &
23 e_ra
=e_ra
, e_rb
=e_rb
, &
24 e_ndra
=e_ndra
, e_ndrb
=e_ndrb
, &
25 e_ra_ra
=e_ra_ra
, e_ra_rb
=e_ra_rb
, e_rb_rb
=e_rb_rb
,&
26 e_ra_ndra
=e_ndra_ra
, e_ra_ndrb
=e_ndrb_ra
, &
27 e_rb_ndrb
=e_ndrb_rb
, e_rb_ndra
=e_ndra_rb
,&
28 e_ndra_ndra
=e_ndra_ndra
, e_ndrb_ndrb
=e_ndrb_ndrb
,&
29 e_ndra_ndrb
=e_ndra_ndrb
,&
30 grad_deriv
=grad_deriv
, npoints
=npoints
, &
31 epsilon_rho
=epsilon_rho
,epsilon_drho
=epsilon_drho
,&
32 param
=param
,scale_c_in
=scale_c
,scale_x_in
=scale_x
)
34 END SUBROUTINE b97_lsd_eval
35 SUBROUTINE b97_lsd_calc(rhoa
, rhob
, norm_drhoa
, norm_drhob
,&
36 e_0
, e_ra
, e_rb
, e_ndra
, e_ndrb
, &
37 e_ra_ndra
,e_ra_ndrb
, e_rb_ndra
, e_rb_ndrb
,&
38 e_ndra_ndra
, e_ndrb_ndrb
, e_ndra_ndrb
, &
39 e_ra_ra
, e_ra_rb
, e_rb_rb
,&
40 grad_deriv
,npoints
,epsilon_rho
,epsilon_drho
, &
41 param
, scale_c_in
, scale_x_in
)
42 REAL(kind
=dp
), DIMENSION(*), INTENT(in
) :: rhoa
, rhob
, norm_drhoa
, &
44 REAL(kind
=dp
), DIMENSION(*), INTENT(inout
) :: e_0
, e_ra
, e_rb
, e_ndra
, &
45 e_ndrb
, e_ra_ndra
, e_ra_ndrb
, e_rb_ndra
, e_rb_ndrb
, e_ndra_ndra
, &
46 e_ndrb_ndrb
, e_ndra_ndrb
, e_ra_ra
, e_ra_rb
, e_rb_rb
47 INTEGER, INTENT(in
) :: grad_deriv
, npoints
48 REAL(kind
=dp
), INTENT(in
) :: epsilon_rho
, epsilon_drho
49 INTEGER, INTENT(in
) :: param
50 REAL(kind
=dp
), INTENT(in
) :: scale_c_in
, scale_x_in
51 REAL(kind
=dp
) :: A_1
, A_2
, A_3
, alpha_1_1
, alpha_1_2
, alpha_1_3
, alpha_c
, &
52 rs_b
, rs_brhob
, rs_brhobrhob
, rsrhoa
, rsrhoarhoa
, rsrhoarhob
, rsrhob
, &
53 t1014
, t102
, t1047
, t1049
, t105
, t106
, t107
54 rsrhoa
= -t4
* t212
* t208
/ 0.12e2_dp
55 t235
= t224
* rsrhoa
/ 0.2e1_dp
+ beta_2_1
* rsrhoa
+ &
56 0.3e1_dp
/ 0.2e1_dp
* t228
* rsrhoa
+ t50
* t48
* rsrhoa
* t232
58 e_c_u_0rhoa
= -0.2e1_dp
* t216
* rsrhoa
* t56
+ t222
* t237
59 epsilon_c_unifrhoa
= e_c_u_0rhoa
+ t285
* t110
+ t287
* t110
- &
60 t293
+ t295
* t108
+ t297
* t108
+ t301
61 e_lsda_c_abrhoa
= epsilon_c_unifrhoa
* rho
+ epsilon_c_unif
- e_lsda_c_arhoa
62 exc_rhoa
= scale_x
* (e_lsda_x_arhoa
* gx_a
+ e_lsda_x_a
* gx_arhoa
) + &
63 scale_c
* (e_lsda_c_abrhoa
* gc_ab
+ e_lsda_c_ab
* gc_abrhoa
+ &
64 e_lsda_c_arhoa
* gc_a
+ e_lsda_c_a
* gc_arhoa
)
65 e_ra(ii
)=e_ra(ii
)+exc_rhoa
66 END SUBROUTINE b97_lsd_calc