1 /*
2
3
4 !
5 ! Dalton, a molecular electronic structure program
6 ! Copyright (C) by the authors of Dalton.
7 !
8 ! This program is free software; you can redistribute it and/or
9 ! modify it under the terms of the GNU Lesser General Public
10 ! License version 2.1 as published by the Free Software Foundation.
11 !
12 ! This program is distributed in the hope that it will be useful,
13 ! but WITHOUT ANY WARRANTY; without even the implied warranty of
14 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 ! Lesser General Public License for more details.
16 !
17 ! If a copy of the GNU LGPL v2.1 was not distributed with this
18 ! code, you can obtain one at https://www.gnu.org/licenses/old-licenses/lgpl-2.1.en.html.
19 !
20
21 !
22
23 */
24 /* fun-example4.c:
25 implementation of Example4 functional and its derivatives
26 (c) Pawel Salek, pawsa@theochem.kth.se, aug 2001
27 NOTE:
28 this file may seem unnecessarily complex but the structure really pays off
29 when implementing multiple functionals depending on different parameters.
30 */
31 /* strictly conform to XOPEN ANSI C standard */
32 #define __USE_XOPEN
33
34 #include <math.h>
35 #include <stdio.h>
36 #include "general.h"
37
38 #define __CVERSION__
39
40 #include "functionals.h"
41
42 /* INTERFACE PART */
example4_isgga(void)43 static integer example4_isgga(void) { return 1; }
44 static integer example4_read(const char* conf_line);
45 static real example4_energy(const FunDensProp* dp);
46 static void example4_first(FunFirstFuncDrv *ds, real factor, const FunDensProp* dp);
47 static void example4_second(FunSecondFuncDrv *ds, real factor, const FunDensProp* dp);
48 static void example4_third(FunThirdFuncDrv *ds, real factor, const FunDensProp* dp);
49
50 Functional Example4Functional = {
51 "Example4", /* name */
52 example4_isgga, /* gga-corrected */
53 1,
54 example4_read,
55 NULL,
56 example4_energy,
57 example4_first,
58 example4_second,
59 example4_third
60 };
61
62 /* IMPLEMENTATION PART */
63 static integer
example4_read(const char * conf_line)64 example4_read(const char* conf_line)
65 {
66 fun_set_hf_weight(0);
67 return 1;
68 }
69
70 static const real PREF= -5e-5;
71
72 static real
example4_energy(const FunDensProp * dp)73 example4_energy(const FunDensProp* dp)
74 {
75 real grad2 = dp->grada*dp->grada+dp->gradb*dp->gradb+2.0*dp->gradab;
76 return PREF*dp->rhoa*dp->rhob*grad2;
77 }
78
79 static void
example4_first(FunFirstFuncDrv * ds,real factor,const FunDensProp * dp)80 example4_first(FunFirstFuncDrv *ds, real factor, const FunDensProp* dp)
81 {
82 real grada2 = dp->grada*dp->grada;
83 real gradb2 = dp->gradb*dp->gradb;
84 real grad2 = grada2+gradb2+2.0*dp->gradab;
85 real rhoab = dp->rhoa*dp->rhob;
86 ds->df1000 += PREF*factor*dp->rhob*grad2;
87 ds->df0100 += PREF*factor*dp->rhoa*grad2;
88 ds->df0010 += 2.0*PREF*factor*rhoab*dp->grada;
89 ds->df0001 += 2.0*PREF*factor*rhoab*dp->gradb;
90 ds->df00001 += 2.0*PREF*factor*rhoab;
91 }
92 static void
example4_second(FunSecondFuncDrv * ds,real factor,const FunDensProp * dp)93 example4_second(FunSecondFuncDrv *ds, real factor, const FunDensProp* dp)
94 {
95 real grada2 = dp->grada*dp->grada;
96 real gradb2 = dp->gradb*dp->gradb;
97 real grad2 = grada2+gradb2+2.0*dp->gradab;
98 real rhoab = dp->rhoa*dp->rhob;
99 /* first derivatives */
100 ds->df1000 += PREF*factor*dp->rhob*grad2;
101 ds->df0100 += PREF*factor*dp->rhoa*grad2;
102 ds->df0010 += 2.0*PREF*factor*rhoab*dp->grada;
103 ds->df0001 += 2.0*PREF*factor*rhoab*dp->gradb;
104 ds->df00001 += 2.0*PREF*factor*rhoab;
105 /* second derivatives */
106 ds->df0020 += 2.0*PREF*factor*rhoab;
107 ds->df0002 += 2.0*PREF*factor*rhoab;
108 /* mixed derivatives */
109 ds->df1100 += PREF*factor*grad2;
110 ds->df1010 += 2.0*PREF*factor*dp->rhob*dp->grada;
111 ds->df1001 += 2.0*PREF*factor*dp->rhob*dp->gradb;
112 ds->df0101 += 2.0*PREF*factor*dp->rhoa*dp->gradb;
113 ds->df0110 += 2.0*PREF*factor*dp->rhoa*dp->grada;
114 ds->df10001 += 2.0*PREF*factor*dp->rhob;
115 ds->df01001 += 2.0*PREF*factor*dp->rhoa;
116 }
117
118 /* example4_third:
119 Example4 functional derivatives.
120 */
121 static void
example4_third(FunThirdFuncDrv * ds,real factor,const FunDensProp * dp)122 example4_third(FunThirdFuncDrv *ds, real factor, const FunDensProp* dp)
123 {
124 real r2 = dp->rhoa*dp->rhoa;
125 real denom = 0.3+r2;
126 real d2 = denom*denom;
127
128 ds->df1000 += PREF*(0.3-r2)/d2*dp->grada*dp->grada*factor;
129 ds->df0010 += PREF*dp->rhoa/denom*2*dp->grada*factor;
130 ds->df2000 += PREF*2*(r2*dp->rhoa-0.9*dp->rhoa)/(denom*d2)*dp->grada*dp->grada*factor;
131 ds->df1010 += PREF*(0.3-r2)/d2*2*dp->grada*factor;
132 ds->df0020 += PREF*dp->rhoa/denom*2*factor;
133
134 ds->df3000 += PREF*(-6*r2*r2+10.8*r2-0.54)/(d2*d2)*dp->grada*dp->grada*factor;
135 ds->df2010 += PREF*2*(r2*dp->rhoa-0.9*dp->rhoa)/(denom*d2)*2*dp->grada*factor;
136 ds->df1020 += PREF*(0.3-r2)/d2*2*factor;
137 }
138