1(*
2 Copyright (C) 2017 M.A.L. Marques
3               2019 Susi Lehtola
4
5 This Source Code Form is subject to the terms of the Mozilla Public
6 License, v. 2.0. If a copy of the MPL was not distributed with this
7 file, You can obtain one at http://mozilla.org/MPL/2.0/.
8*)
9
10(* type: gga_exc *)
11
12$define lda_c_pw_params
13$include "lda_c_pw.mpl"
14
15(* Parameters from Table 6 *)
16gap_par0 := [
17   0.04953, 1.07924, 0.07928, (* a1, a2, a3 *)
18  -2.504e-2, 7.026e-3, -1.268e-3, 1.136e-4, -3.841e-6, (* b3, b4, b5, b6, b7 *)
19   0.031091,  (* params_a_a of lda_c_pw *)
20   0.23878    (* pre-factor of C = 0.06483*((9*Pi)/4)^(2/3) *)
21]:
22gap_par1 := [
23   0.0471985, 1.49676, 0.00179054,
24  -3.24091e-2, 9.99978e-3, -1.93483e-3, 1.79118e-4, -6.15798e-6,
25   0.015545,
26   0.064535
27]:
28
29(* Equation (20): e'(rs) *)
30gap_eps_1 := (rs, par) ->
31  par[1]*rs^(3/2)/(1 + sqrt(rs)*(par[2] + par[3]*sqrt(rs) + par[1]*rs)):
32(* Equation (21): e''(rs) *)
33gap_eps_2 := (rs, par) ->
34  add(par[i+1]*rs^i, i=3..7):
35
36(* Equation (19) *)
37gap_C := (rs, par) -> par[10]/rs^2:
38
39(* Equation (17) *)
40gap_c2 := (rs, z, par) ->
41  + (2*f_pw(rs, z)*gap_eps_1(rs, par) - gap_C(rs, par)*gap_eps_2(rs, par))
42  / (2*(gap_C(rs, par)*gap_eps_1(rs, par) - f_pw(rs, z)^2)):
43(* Equation (18) *)
44gap_c3 := (rs, z, par) ->
45  - (2*gap_eps_1(rs, par)^2 - f_pw(rs, z)*gap_eps_2(rs, par))
46  / (2*(gap_C(rs, par)*gap_eps_1(rs, par) - f_pw(rs, z)^2)):
47(* Equation (16) *)
48gap_c1 := (rs, z, par) ->
49  - gap_C(rs, par) * gap_c3(rs, z, par):
50
51(* after Equation (6): a = 30 is a parameter fixed by minimizing the
52variance of the correlation energy error for the noble gas atoms He,
53Ne, and Ar *)
54gap_par_a := 30:
55
56(* Equation (6) *)
57gap_H := (rs, t, par) ->
58  (gap_par_a  + par[9]*rs*log(rs)*t^2/beta_Hu_Langreth(rs))/(gap_par_a + t^2):
59
60gap_t := (rs, z, xt) ->
61  xt*n_total(rs)^(1/6)/(4*mphi(z)*(3/Pi)^(1/6)):
62
63(* Gap function, Equation (5) *)
64gap_G := (rs, z, xt, par) ->
65  + mphi(z)^3*beta_Hu_Langreth(rs)*gap_t(rs, z, xt)^2
66  * gap_H(rs, gap_t(rs, z, xt), par)
67  / (gap_c1(rs, z, par) - gap_c2(rs, z, par)*f_pw(rs, z)):
68
69(* Correlation energy per particle, Equation (4) *)
70gap_eps := (rs, z, xt, par) ->
71  + (f_pw(rs, z) + gap_c1(rs, z, par)*gap_G(rs, z, xt, par))
72  / (1 + gap_c2(rs, z, par)*gap_G(rs, z, xt, par) + gap_c3(rs, z, par)*gap_G(rs, z, xt, par)^2):
73
74(* Total energy, Equation (2) *)
75f_gap := (rs, z, xt) ->
76  + gap_eps(rs, 0, xt, gap_par0)
77  + f_zeta(z)*(gap_eps(rs, 1, xt, gap_par1) - gap_eps(rs, 0, xt, gap_par0)):
78
79f  := (rs, z, xt, xs0, xs1) ->
80  f_gap(rs, z, xt):
81