MPQC
2.3.1
Main Page
Related Pages
Classes
Files
File List
src
lib
chemistry
qc
dft
integrator.h
1
//
2
// integrator.h --- definition of the electron density integrator
3
//
4
// Copyright (C) 1997 Limit Point Systems, Inc.
5
//
6
// Author: Curtis Janssen <cljanss@limitpt.com>
7
// Maintainer: LPS
8
//
9
// This file is part of the SC Toolkit.
10
//
11
// The SC Toolkit is free software; you can redistribute it and/or modify
12
// it under the terms of the GNU Library General Public License as published by
13
// the Free Software Foundation; either version 2, or (at your option)
14
// any later version.
15
//
16
// The SC Toolkit is distributed in the hope that it will be useful,
17
// but WITHOUT ANY WARRANTY; without even the implied warranty of
18
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19
// GNU Library General Public License for more details.
20
//
21
// You should have received a copy of the GNU Library General Public License
22
// along with the SC Toolkit; see the file COPYING.LIB. If not, write to
23
// the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
24
//
25
// The U.S. Government is granted a limited license as per AL 91-7.
26
//
27
28
#ifndef _chemistry_qc_dft_integrator_h
29
#define _chemistry_qc_dft_integrator_h
30
31
#ifdef __GNUC__
32
#pragma interface
33
#endif
34
35
#include <util/state/state.h>
36
#include <util/group/thread.h>
37
#include <chemistry/qc/dft/functional.h>
38
#include <chemistry/qc/basis/extent.h>
39
#include <chemistry/qc/wfn/density.h>
40
41
namespace
sc {
42
44
class
DenIntegrator
:
virtual
public
SavableState
{
45
protected
:
46
Ref<Wavefunction>
wfn_;
47
//clj Ref<ShellExtent> extent_;
48
Ref<BatchElectronDensity>
den_;
49
50
Ref<ThreadGrp>
threadgrp_;
51
Ref<MessageGrp>
messagegrp_;
52
53
double
value_;
54
double
accuracy_;
55
56
double
*alpha_vmat_;
57
double
*beta_vmat_;
58
59
//clj double *alpha_dmat_;
60
//clj double *beta_dmat_;
61
//clj double *dmat_bound_;
62
63
int
spin_polarized_;
64
65
int
need_density_;
// specialization must set to 1 if it needs density_
66
double
density_;
67
int
nbasis_;
68
int
nshell_;
69
int
n_integration_center_;
70
int
natom_;
71
int
compute_potential_integrals_;
// 1 if potential integrals are needed
72
73
int
linear_scaling_;
74
int
use_dmat_bound_;
75
76
void
init_integration(
const
Ref<DenFunctional>
&func,
77
const
RefSymmSCMatrix
& densa,
78
const
RefSymmSCMatrix
& densb,
79
double
*nuclear_gradient);
80
void
done_integration();
81
82
void
init_object();
83
public
:
85
DenIntegrator
();
87
DenIntegrator
(
const
Ref<KeyVal>
&);
89
DenIntegrator
(
StateIn
&);
90
~
DenIntegrator
();
91
void
save_data_state
(
StateOut
&);
92
94
Ref<Wavefunction>
wavefunction
()
const
{
return
wfn_; }
96
double
value
()
const
{
return
value_; }
97
99
void
set_accuracy
(
double
a);
100
double
get_accuracy(
void
) {
return
accuracy_; }
103
void
set_compute_potential_integrals
(
int
);
106
const
double
*
alpha_vmat
()
const
{
return
alpha_vmat_; }
109
const
double
*
beta_vmat
()
const
{
return
beta_vmat_; }
110
113
virtual
void
init
(
const
Ref<Wavefunction>
&);
115
virtual
void
done
();
119
virtual
void
integrate
(
const
Ref<DenFunctional>
&,
120
const
RefSymmSCMatrix
& densa =0,
121
const
RefSymmSCMatrix
& densb =0,
122
double
*nuclear_grad = 0) = 0;
123
};
124
125
127
class
IntegrationWeight
:
virtual
public
SavableState
{
128
129
protected
:
130
131
Ref<Molecule>
mol_;
132
double
tol_;
133
134
void
fd_w(
int
icenter,
SCVector3
&
point
,
double
*fd_grad_w);
135
136
public
:
137
IntegrationWeight
();
138
IntegrationWeight
(
const
Ref<KeyVal>
&);
139
IntegrationWeight
(
StateIn
&);
140
~
IntegrationWeight
();
141
void
save_data_state
(
StateOut
&);
142
143
void
test(
int
icenter,
SCVector3
&point);
144
void
test();
145
147
virtual
void
init
(
const
Ref<Molecule>
&,
double
tolerance);
149
virtual
void
done
();
154
virtual
double
w
(
int
center,
SCVector3
&point,
double
*grad_w = 0) = 0;
155
};
156
157
159
class
BeckeIntegrationWeight
:
public
IntegrationWeight
{
160
161
int
n_integration_centers;
162
SCVector3
*centers;
163
double
*atomic_radius;
164
165
double
**a_mat;
166
double
**oorab;
167
168
void
compute_grad_p(
int
gc,
int
ic,
int
wc,
SCVector3
&r,
double
p,
169
SCVector3
&g);
170
void
compute_grad_nu(
int
gc,
int
bc,
SCVector3
&
point
,
SCVector3
&grad);
171
172
double
compute_t(
int
ic,
int
jc,
SCVector3
&point);
173
double
compute_p(
int
icenter,
SCVector3
&point);
174
175
public
:
176
BeckeIntegrationWeight
();
177
BeckeIntegrationWeight
(
const
Ref<KeyVal>
&);
178
BeckeIntegrationWeight
(
StateIn
&);
179
~
BeckeIntegrationWeight
();
180
void
save_data_state
(
StateOut
&);
181
182
void
init
(
const
Ref<Molecule>
&,
double
tolerance);
183
void
done
();
184
double
w
(
int
center,
SCVector3
&point,
double
*grad_w = 0);
185
};
186
188
class
RadialIntegrator
:
virtual
public
SavableState
{
189
protected
:
190
int
nr_;
191
public
:
192
RadialIntegrator
();
193
RadialIntegrator
(
const
Ref<KeyVal>
&);
194
RadialIntegrator
(
StateIn
&);
195
~
RadialIntegrator
();
196
void
save_data_state
(
StateOut
&);
197
198
virtual
int
nr()
const
= 0;
199
virtual
double
radial_value(
int
ir,
int
nr,
double
radii,
200
double
&multiplier) = 0;
201
};
202
203
205
class
AngularIntegrator
:
virtual
public
SavableState
{
206
protected
:
207
public
:
208
AngularIntegrator
();
209
AngularIntegrator
(
const
Ref<KeyVal>
&);
210
AngularIntegrator
(
StateIn
&);
211
~
AngularIntegrator
();
212
void
save_data_state
(
StateOut
&);
213
214
virtual
int
nw(
void
)
const
= 0;
215
virtual
int
num_angular_points(
double
r_value,
int
ir) = 0;
216
virtual
double
angular_point_cartesian(
int
iangular,
double
r,
217
SCVector3
&integration_point)
const
= 0;
218
};
219
220
223
class
EulerMaclaurinRadialIntegrator
:
public
RadialIntegrator
{
224
public
:
225
EulerMaclaurinRadialIntegrator
();
226
EulerMaclaurinRadialIntegrator
(
int
i);
230
EulerMaclaurinRadialIntegrator
(
const
Ref<KeyVal>
&);
231
EulerMaclaurinRadialIntegrator
(
StateIn
&);
232
~
EulerMaclaurinRadialIntegrator
();
233
void
save_data_state
(
StateOut
&);
234
235
int
nr()
const
;
236
double
radial_value(
int
ir,
int
nr,
double
radii,
double
&multiplier);
237
238
void
print
(std::ostream & =
ExEnv::out0
())
const
;
239
};
240
282
class
LebedevLaikovIntegrator
:
public
AngularIntegrator
{
283
protected
:
284
int
npoint_;
285
double
*x_, *y_, *z_, *w_;
286
287
void
init(
int
n);
288
public
:
289
LebedevLaikovIntegrator
();
295
LebedevLaikovIntegrator
(
const
Ref<KeyVal>
&);
296
LebedevLaikovIntegrator
(
StateIn
&);
297
LebedevLaikovIntegrator
(
int
);
298
~
LebedevLaikovIntegrator
();
299
void
save_data_state
(
StateOut
&);
300
301
int
nw(
void
)
const
;
302
int
num_angular_points(
double
r_value,
int
ir);
303
double
angular_point_cartesian(
int
iangular,
double
r,
304
SCVector3
&integration_point)
const
;
305
void
print
(std::ostream & =
ExEnv::out0
())
const
;
306
};
307
310
class
GaussLegendreAngularIntegrator
:
public
AngularIntegrator
{
311
protected
:
312
int
ntheta_;
313
int
nphi_;
314
int
Ktheta_;
315
int
ntheta_r_;
316
int
nphi_r_;
317
int
Ktheta_r_;
318
double
*theta_quad_weights_;
319
double
*theta_quad_points_;
320
321
int
get_ntheta(
void
)
const
;
322
void
set_ntheta(
int
i);
323
int
get_nphi(
void
)
const
;
324
void
set_nphi(
int
i);
325
int
get_Ktheta(
void
)
const
;
326
void
set_Ktheta(
int
i);
327
int
get_ntheta_r(
void
)
const
;
328
void
set_ntheta_r(
int
i);
329
int
get_nphi_r(
void
)
const
;
330
void
set_nphi_r(
int
i);
331
int
get_Ktheta_r(
void
)
const
;
332
void
set_Ktheta_r(
int
i);
333
int
nw(
void
)
const
;
334
double
sin_theta(
SCVector3
&
point
)
const
;
335
void
gauleg(
double
x1,
double
x2,
int
n);
336
public
:
337
GaussLegendreAngularIntegrator
();
344
GaussLegendreAngularIntegrator
(
const
Ref<KeyVal>
&);
345
GaussLegendreAngularIntegrator
(
StateIn
&);
346
~
GaussLegendreAngularIntegrator
();
347
void
save_data_state
(
StateOut
&);
348
349
int
num_angular_points(
double
r_value,
int
ir);
350
double
angular_point_cartesian(
int
iangular,
double
r,
351
SCVector3
&integration_point)
const
;
352
void
print
(std::ostream & =
ExEnv::out0
())
const
;
353
};
354
357
class
RadialAngularIntegrator
:
public
DenIntegrator
{
358
private
:
359
int
prune_grid_;
360
double
**Alpha_coeffs_;
361
int
gridtype_;
362
int
**nr_points_, *xcoarse_l_;
363
int
npruned_partitions_;
364
double
*grid_accuracy_;
365
int
dynamic_grids_;
366
int
natomic_rows_;
367
int
max_gridtype_;
368
protected
:
369
Ref<IntegrationWeight>
weight_;
370
Ref<RadialIntegrator>
radial_user_;
371
Ref<AngularIntegrator>
angular_user_;
372
Ref<AngularIntegrator>
***angular_grid_;
373
Ref<RadialIntegrator>
**radial_grid_;
374
public
:
375
RadialAngularIntegrator
();
416
RadialAngularIntegrator
(
const
Ref<KeyVal>
&);
417
RadialAngularIntegrator
(
StateIn
&);
418
~
RadialAngularIntegrator
();
419
void
save_data_state
(
StateOut
&);
420
421
void
integrate
(
const
Ref<DenFunctional>
&,
422
const
RefSymmSCMatrix
& densa =0,
423
const
RefSymmSCMatrix
& densb =0,
424
double
*nuclear_gradient = 0);
425
void
print
(std::ostream & =
ExEnv::out0
())
const
;
426
AngularIntegrator
*get_angular_grid(
double
radius,
double
atomic_radius,
427
int
charge,
int
deriv_order);
428
RadialIntegrator
*get_radial_grid(
int
charge,
int
deriv_order);
429
void
init_default_grids(
void
);
430
int
angular_grid_offset(
int
i);
431
void
set_grids(
void
);
432
int
get_atomic_row(
int
i);
433
void
init_parameters(
void
);
434
void
init_parameters(
const
Ref<KeyVal>
& keyval);
435
void
init_pruning_coefficients(
const
Ref<KeyVal>
& keyval);
436
void
init_pruning_coefficients(
void
);
437
void
init_alpha_coefficients(
void
);
438
int
select_dynamic_grid(
void
);
439
Ref<IntegrationWeight>
weight() {
return
weight_; }
440
};
441
442
}
443
444
#endif
445
446
// Local Variables:
447
// mode: c++
448
// c-file-style: "CLJ"
449
// End:
Generated at Sun Oct 6 2013 05:58:42 for
MPQC
2.3.1 using the documentation package
Doxygen
1.8.3.1.