13 double *
const opmata;
14 double *
const opmatb;
17 double *ka,
double *opa,
double *kb,
double *opb) :
18 gmata(ga), gmatb(gb), kmata(ka), kmatb(kb),
19 pmata(pa), pmatb(pb), opmata(opa), opmatb(opb) {}
22 void set_bound(
double,
double) {}
24 inline void cont1(
int ij,
int kl,
double val) {
25 gmata[ij] += val*pmata[kl];
26 gmata[kl] += val*pmata[ij];
28 gmatb[ij] += val*pmatb[kl];
29 gmatb[kl] += val*pmatb[ij];
32 inline void cont2(
int ij,
int kl,
double val) {
34 gmata[ij] -= val*pmata[kl];
35 gmata[kl] -= val*pmata[ij];
37 gmatb[ij] -= val*pmatb[kl];
38 gmatb[kl] -= val*pmatb[ij];
40 kmata[ij] += val*opmata[kl];
41 kmata[kl] += val*opmata[ij];
43 kmatb[ij] += val*opmatb[kl];
44 kmatb[kl] += val*opmatb[ij];
47 inline void cont3(
int ij,
int kl,
double val) {
49 gmata[ij] -= val*pmata[kl];
50 gmata[kl] -= val*pmata[ij];
52 gmatb[ij] -= val*pmatb[kl];
53 gmatb[kl] -= val*pmatb[ij];
55 kmata[ij] += val*opmata[kl];
56 kmata[kl] += val*opmata[ij];
58 kmatb[ij] += val*opmatb[kl];
59 kmatb[kl] += val*opmatb[ij];
62 inline void cont4(
int ij,
int kl,
double val) {
63 gmata[ij] += 0.75*val*pmata[kl];
64 gmata[kl] += 0.75*val*pmata[ij];
66 gmatb[ij] += 0.75*val*pmatb[kl];
67 gmatb[kl] += 0.75*val*pmatb[ij];
69 kmata[ij] += 0.25*val*opmata[kl];
70 kmata[kl] += 0.25*val*opmata[ij];
72 kmatb[ij] += 0.25*val*opmatb[kl];
73 kmatb[kl] += 0.25*val*opmatb[ij];
76 inline void cont5(
int ij,
int kl,
double val) {
78 gmata[ij] += val*pmata[kl];
79 gmata[kl] += val*pmata[ij];
81 gmatb[ij] += val*pmatb[kl];
82 gmatb[kl] += val*pmatb[ij];
84 kmata[ij] += val*opmata[kl];
85 kmata[kl] += val*opmata[ij];
87 kmatb[ij] += val*opmatb[kl];
88 kmatb[kl] += val*opmatb[ij];
96 double *
const opmata;
97 double *
const opmatb;
108 pmata(pa), pmatb(pb), opmata(opa), opmatb(opb) {
115 void set_bound(
double,
double) {};
117 inline void cont1(
int ij,
int kl,
double val) {
118 eca += val*pmata[ij]*pmata[kl];
119 ecb += val*pmatb[ij]*pmatb[kl];
120 ecab += val*(pmata[ij]*pmatb[kl]-pmatb[ij]*pmata[kl]);
123 inline void cont2(
int ij,
int kl,
double val) {
124 exa -= 0.25*val*pmata[ij]*pmata[kl];
125 exb -= 0.25*val*pmatb[ij]*pmatb[kl];
126 exab -= 0.25*val*(pmata[ij]*pmatb[kl]-pmatb[ij]*pmata[kl]);
129 inline void cont3(
int ij,
int kl,
double val) {
130 exa -= 0.5*val*pmata[ij]*pmata[kl];
131 exb -= 0.5*val*pmatb[ij]*pmatb[kl];
132 exab -= 0.5*val*(pmata[ij]*pmatb[kl]-pmatb[ij]*pmata[kl]);
135 inline void cont4(
int ij,
int kl,
double val) {
136 eca += val*pmata[ij]*pmata[kl];
137 ecb += val*pmatb[ij]*pmatb[kl];
138 ecab += val*(pmata[ij]*pmatb[kl]-pmatb[ij]*pmata[kl]);
140 exa -= 0.25*val*pmata[ij]*pmata[kl];
141 exb -= 0.25*val*pmatb[ij]*pmatb[kl];
142 exab -= 0.25*val*(pmata[ij]*pmatb[kl]-pmatb[ij]*pmata[kl]);
145 inline void cont5(
int ij,
int kl,
double val) {
146 eca += val*pmata[ij]*pmata[kl];
147 ecb += val*pmatb[ij]*pmatb[kl];
148 ecab += val*(pmata[ij]*pmatb[kl]-pmatb[ij]*pmata[kl]);
150 exa -= 0.5*val*pmata[ij]*pmata[kl];
151 exb -= 0.5*val*pmatb[ij]*pmatb[kl];
152 exab -= 0.5*val*(pmata[ij]*pmatb[kl]-pmatb[ij]*pmata[kl]);
159 double *
const pmata;
160 double *
const pmatb;
167 double c1,
double c2) :
168 pmat(p), pmata(pa), pmatb(pb)
176 inline double cont1(
int ij,
int kl) {
177 return pmat[ij]*pmat[kl] +
178 c1sq*(pmata[ij]*pmat[kl] + pmat[ij]*pmata[kl]) +
179 c2sq*(pmatb[ij]*pmat[kl] + pmat[ij]*pmatb[kl]) +
180 0.5*c1sq*pmata[ij]*pmata[kl] +
181 0.5*c2sq*pmatb[ij]*pmatb[kl];
184 inline double cont2(
int ij,
int kl) {
185 return pmat[ij]*pmat[kl] +
186 c1sq*(pmata[ij]*pmat[kl] + pmat[ij]*pmata[kl]) +
187 c2sq*(pmatb[ij]*pmat[kl] + pmat[ij]*pmatb[kl]) -
188 c1c2*(pmata[ij]*pmatb[kl] + pmatb[ij]*pmata[kl]);