mirror of
https://github.com/TREX-CoE/qmckl.git
synced 2024-08-11 15:28:31 +02:00
Merge branch 'master' into qmckl_dgemm_integration
This commit is contained in:
vijay
GitHub
commit
eb27207bd9
@ -118,7 +118,7 @@ python-install: $(qmckl_h) $(lib_LTLIBRARIES) $(dist_python_DATA)
|
|||||||
cp src/qmckl.py . ; \
|
cp src/qmckl.py . ; \
|
||||||
export QMCKL_INCLUDEDIR="$(prefix)/include" ; \
|
export QMCKL_INCLUDEDIR="$(prefix)/include" ; \
|
||||||
export QMCKL_LIBDIR="$(prefix)/lib" ; \
|
export QMCKL_LIBDIR="$(prefix)/lib" ; \
|
||||||
pip install .
|
python3 -m pip install .
|
||||||
|
|
||||||
python-test: $(test_py)
|
python-test: $(test_py)
|
||||||
cd $(abs_srcdir)/python/test/ && \
|
cd $(abs_srcdir)/python/test/ && \
|
||||||
|
|||||||
21
README.md
21
README.md
@ -75,6 +75,27 @@ sudo make install
|
|||||||
sudo make installcheck
|
sudo make installcheck
|
||||||
```
|
```
|
||||||
|
|
||||||
|
## Python API
|
||||||
|
|
||||||
|
- [SWIG](https://www.swig.org) (>= 4.0) is required to build the Python API for maintainers
|
||||||
|
|
||||||
|
In order to install the `qmckl` Python package, first install the shared C library
|
||||||
|
`libqmckl` following the installation guide above and then run the following command:
|
||||||
|
```
|
||||||
|
make python-install
|
||||||
|
```
|
||||||
|
|
||||||
|
To test the installation, run
|
||||||
|
```
|
||||||
|
make python-test
|
||||||
|
```
|
||||||
|
|
||||||
|
Minimal example demonstrating the use of the `qmckl` Python API can be found in the
|
||||||
|
[test_api.py](https://github.com/TREX-CoE/qmckl/blob/master/python/test/test_api.py) file.
|
||||||
|
|
||||||
|
We highly recommend to use
|
||||||
|
[virtual environments](https://docs.python.org/3/tutorial/venv.html)
|
||||||
|
to avoid compatibility issues and to improve reproducibility.
|
||||||
|
|
||||||
## Installation procedure for Guix users
|
## Installation procedure for Guix users
|
||||||
|
|
||||||
|
|||||||
211
org/qmckl_ao.org
211
org/qmckl_ao.org
@ -5648,33 +5648,14 @@ end function qmckl_compute_ao_value_doc_f
|
|||||||
| ~coef_normalized~ | ~double[prim_num]~ | in | Value, gradients and Laplacian of the shells |
|
| ~coef_normalized~ | ~double[prim_num]~ | in | Value, gradients and Laplacian of the shells |
|
||||||
| ~ao_value~ | ~double[point_num][ao_num]~ | out | Values of the AOs |
|
| ~ao_value~ | ~double[point_num][ao_num]~ | out | Values of the AOs |
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
#+NAME:ao_value_constants
|
||||||
#ifdef HAVE_HPC
|
#+begin_src c :exports none
|
||||||
qmckl_exit_code
|
int32_t lstart[32] __attribute__((aligned(64)));
|
||||||
qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|
||||||
const int64_t ao_num,
|
|
||||||
const int64_t shell_num,
|
|
||||||
const int32_t* restrict prim_num_per_nucleus,
|
|
||||||
const int64_t point_num,
|
|
||||||
const int64_t nucl_num,
|
|
||||||
const double* restrict coord,
|
|
||||||
const double* restrict nucl_coord,
|
|
||||||
const int64_t* restrict nucleus_index,
|
|
||||||
const int64_t* restrict nucleus_shell_num,
|
|
||||||
const double* nucleus_range,
|
|
||||||
const int32_t* restrict nucleus_max_ang_mom,
|
|
||||||
const int32_t* restrict shell_ang_mom,
|
|
||||||
const double* restrict ao_factor,
|
|
||||||
const qmckl_matrix expo_per_nucleus,
|
|
||||||
const qmckl_tensor coef_per_nucleus,
|
|
||||||
double* restrict const ao_value )
|
|
||||||
{
|
|
||||||
int32_t lstart[32];
|
|
||||||
for (int32_t l=0 ; l<32 ; ++l) {
|
for (int32_t l=0 ; l<32 ; ++l) {
|
||||||
lstart[l] = l*(l+1)*(l+2)/6;
|
lstart[l] = l*(l+1)*(l+2)/6;
|
||||||
}
|
}
|
||||||
|
|
||||||
int64_t ao_index[shell_num+1];
|
int64_t ao_index[shell_num+1] __attribute__((aligned(64)));
|
||||||
int64_t size_max = 0;
|
int64_t size_max = 0;
|
||||||
int64_t prim_max = 0;
|
int64_t prim_max = 0;
|
||||||
int64_t shell_max = 0;
|
int64_t shell_max = 0;
|
||||||
@ -5700,17 +5681,42 @@ qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|||||||
/* Don't compute polynomials when the radial part is zero. */
|
/* Don't compute polynomials when the radial part is zero. */
|
||||||
double cutoff = 27.631021115928547; // -log(1.e-12)
|
double cutoff = 27.631021115928547; // -log(1.e-12)
|
||||||
|
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
||||||
|
#ifdef HAVE_HPC
|
||||||
|
qmckl_exit_code
|
||||||
|
qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
||||||
|
const int64_t ao_num,
|
||||||
|
const int64_t shell_num,
|
||||||
|
const int32_t* restrict prim_num_per_nucleus,
|
||||||
|
const int64_t point_num,
|
||||||
|
const int64_t nucl_num,
|
||||||
|
const double* restrict coord,
|
||||||
|
const double* restrict nucl_coord,
|
||||||
|
const int64_t* restrict nucleus_index,
|
||||||
|
const int64_t* restrict nucleus_shell_num,
|
||||||
|
const double* nucleus_range,
|
||||||
|
const int32_t* restrict nucleus_max_ang_mom,
|
||||||
|
const int32_t* restrict shell_ang_mom,
|
||||||
|
const double* restrict ao_factor,
|
||||||
|
const qmckl_matrix expo_per_nucleus,
|
||||||
|
const qmckl_tensor coef_per_nucleus,
|
||||||
|
double* restrict const ao_value )
|
||||||
|
{
|
||||||
|
<<ao_value_constants>>
|
||||||
|
|
||||||
#ifdef HAVE_OPENMP
|
#ifdef HAVE_OPENMP
|
||||||
#pragma omp parallel
|
#pragma omp parallel
|
||||||
#endif
|
#endif
|
||||||
{
|
{
|
||||||
qmckl_exit_code rc;
|
qmckl_exit_code rc;
|
||||||
double ar2[prim_max];
|
double ar2[prim_max] __attribute__((aligned(64)));
|
||||||
int32_t powers[3*size_max];
|
int32_t powers[3*size_max] __attribute__((aligned(64)));
|
||||||
double poly_vgl[5*size_max];
|
double poly_vgl[5*size_max] __attribute__((aligned(64)));
|
||||||
|
|
||||||
double exp_mat[prim_max];
|
double exp_mat[prim_max] __attribute__((aligned(64)));
|
||||||
double ce_mat[shell_max];
|
double ce_mat[shell_max] __attribute__((aligned(64)));
|
||||||
|
|
||||||
double coef_mat[nucl_num][shell_max][prim_max];
|
double coef_mat[nucl_num][shell_max][prim_max];
|
||||||
for (int i=0 ; i<nucl_num ; ++i) {
|
for (int i=0 ; i<nucl_num ; ++i) {
|
||||||
@ -5725,15 +5731,17 @@ qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|||||||
#pragma omp for
|
#pragma omp for
|
||||||
#endif
|
#endif
|
||||||
for (int64_t ipoint=0 ; ipoint < point_num ; ++ipoint) {
|
for (int64_t ipoint=0 ; ipoint < point_num ; ++ipoint) {
|
||||||
const double e_coord[3] = { coord[ipoint],
|
const double e_coord[3] __attribute__((aligned(64))) =
|
||||||
coord[ipoint + point_num],
|
{ coord[ipoint],
|
||||||
coord[ipoint + 2*point_num] };
|
coord[ipoint + point_num],
|
||||||
|
coord[ipoint + 2*point_num] };
|
||||||
|
|
||||||
for (int64_t inucl=0 ; inucl < nucl_num ; ++inucl) {
|
for (int64_t inucl=0 ; inucl < nucl_num ; ++inucl) {
|
||||||
const double n_coord[3] = { nucl_coord[inucl],
|
const double n_coord[3] __attribute__((aligned(64))) =
|
||||||
nucl_coord[inucl + nucl_num],
|
{ nucl_coord[inucl],
|
||||||
nucl_coord[inucl + 2*nucl_num] };
|
nucl_coord[inucl + nucl_num],
|
||||||
|
nucl_coord[inucl + 2*nucl_num] };
|
||||||
|
|
||||||
/* Test if the point is in the range of the nucleus */
|
/* Test if the point is in the range of the nucleus */
|
||||||
const double x = e_coord[0] - n_coord[0];
|
const double x = e_coord[0] - n_coord[0];
|
||||||
const double y = e_coord[1] - n_coord[1];
|
const double y = e_coord[1] - n_coord[1];
|
||||||
@ -5751,14 +5759,14 @@ qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|||||||
break;
|
break;
|
||||||
|
|
||||||
case 1:
|
case 1:
|
||||||
poly_vgl[0] = 0.;
|
poly_vgl[0] = 1.;
|
||||||
poly_vgl[1] = x;
|
poly_vgl[1] = x;
|
||||||
poly_vgl[2] = y;
|
poly_vgl[2] = y;
|
||||||
poly_vgl[3] = z;
|
poly_vgl[3] = z;
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 2:
|
case 2:
|
||||||
poly_vgl[0] = 0.;
|
poly_vgl[0] = 1.;
|
||||||
poly_vgl[1] = x;
|
poly_vgl[1] = x;
|
||||||
poly_vgl[2] = y;
|
poly_vgl[2] = y;
|
||||||
poly_vgl[3] = z;
|
poly_vgl[3] = z;
|
||||||
@ -5797,10 +5805,8 @@ qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|||||||
}
|
}
|
||||||
|
|
||||||
for (int i=0 ; i<nucleus_shell_num[inucl] ; ++i) {
|
for (int i=0 ; i<nucleus_shell_num[inucl] ; ++i) {
|
||||||
ce_mat[i] = 0.;
|
ce_mat[i] = 0.;
|
||||||
}
|
for (int k=0 ; k<nidx; ++k) {
|
||||||
for (int k=0 ; k<nidx; ++k) {
|
|
||||||
for (int i=0 ; i<nucleus_shell_num[inucl] ; ++i) {
|
|
||||||
ce_mat[i] += coef_mat[inucl][i][k] * exp_mat[k];
|
ce_mat[i] += coef_mat[inucl][i][k] * exp_mat[k];
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@ -5811,7 +5817,6 @@ qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|||||||
for (int64_t ishell = ishell_start ; ishell < ishell_end ; ++ishell) {
|
for (int64_t ishell = ishell_start ; ishell < ishell_end ; ++ishell) {
|
||||||
|
|
||||||
const double s1 = ce_mat[ishell-ishell_start];
|
const double s1 = ce_mat[ishell-ishell_start];
|
||||||
if (s1 == 0.0) continue;
|
|
||||||
|
|
||||||
const int64_t k = ao_index[ishell];
|
const int64_t k = ao_index[ishell];
|
||||||
double* restrict const ao_value_1 = ao_value + ipoint*ao_num + k;
|
double* restrict const ao_value_1 = ao_value + ipoint*ao_num + k;
|
||||||
@ -5819,6 +5824,13 @@ qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|||||||
const int32_t l = shell_ang_mom[ishell];
|
const int32_t l = shell_ang_mom[ishell];
|
||||||
const int32_t n = lstart[l+1]-lstart[l];
|
const int32_t n = lstart[l+1]-lstart[l];
|
||||||
|
|
||||||
|
if (s1 == 0.0) {
|
||||||
|
for (int64_t il=0 ; il<n ; ++il) {
|
||||||
|
ao_value_1[il] = 0.;
|
||||||
|
}
|
||||||
|
continue;
|
||||||
|
}
|
||||||
|
|
||||||
double* restrict poly_vgl_1 = NULL;
|
double* restrict poly_vgl_1 = NULL;
|
||||||
if (nidx > 0) {
|
if (nidx > 0) {
|
||||||
const double* restrict f = ao_factor + k;
|
const double* restrict f = ao_factor + k;
|
||||||
@ -5827,10 +5839,10 @@ qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|||||||
poly_vgl_1 = &(poly_vgl[idx]);
|
poly_vgl_1 = &(poly_vgl[idx]);
|
||||||
|
|
||||||
switch (n) {
|
switch (n) {
|
||||||
case(1):
|
case 1:
|
||||||
ao_value_1[0] = s1 * f[0];
|
ao_value_1[0] = s1 * f[0];
|
||||||
break;
|
break;
|
||||||
case (3):
|
case 3:
|
||||||
#ifdef HAVE_OPENMP
|
#ifdef HAVE_OPENMP
|
||||||
#pragma omp simd
|
#pragma omp simd
|
||||||
#endif
|
#endif
|
||||||
@ -5848,7 +5860,7 @@ qmckl_compute_ao_value_hpc_gaussian (const qmckl_context context,
|
|||||||
break;
|
break;
|
||||||
default:
|
default:
|
||||||
#ifdef HAVE_OPENMP
|
#ifdef HAVE_OPENMP
|
||||||
#pragma omp simd simdlen(8)
|
#pragma omp simd
|
||||||
#endif
|
#endif
|
||||||
for (int il=0 ; il<n ; ++il) {
|
for (int il=0 ; il<n ; ++il) {
|
||||||
ao_value_1[il] = poly_vgl_1[il] * s1 * f[il];
|
ao_value_1[il] = poly_vgl_1[il] * s1 * f[il];
|
||||||
@ -6470,36 +6482,7 @@ qmckl_compute_ao_vgl_hpc_gaussian (
|
|||||||
const qmckl_tensor coef_per_nucleus,
|
const qmckl_tensor coef_per_nucleus,
|
||||||
double* restrict const ao_vgl )
|
double* restrict const ao_vgl )
|
||||||
{
|
{
|
||||||
int32_t lstart[32];
|
<<ao_value_constants>>
|
||||||
for (int32_t l=0 ; l<32 ; ++l) {
|
|
||||||
lstart[l] = l*(l+1)*(l+2)/6;
|
|
||||||
}
|
|
||||||
|
|
||||||
int64_t ao_index[shell_num+1];
|
|
||||||
int64_t size_max = 0;
|
|
||||||
int64_t prim_max = 0;
|
|
||||||
int64_t shell_max = 0;
|
|
||||||
{
|
|
||||||
int64_t k=0;
|
|
||||||
for (int inucl=0 ; inucl < nucl_num ; ++inucl) {
|
|
||||||
prim_max = prim_num_per_nucleus[inucl] > prim_max ?
|
|
||||||
prim_num_per_nucleus[inucl] : prim_max;
|
|
||||||
shell_max = nucleus_shell_num[inucl] > shell_max ?
|
|
||||||
nucleus_shell_num[inucl] : shell_max;
|
|
||||||
const int64_t ishell_start = nucleus_index[inucl];
|
|
||||||
const int64_t ishell_end = nucleus_index[inucl] + nucleus_shell_num[inucl];
|
|
||||||
for (int64_t ishell = ishell_start ; ishell < ishell_end ; ++ishell) {
|
|
||||||
const int l = shell_ang_mom[ishell];
|
|
||||||
ao_index[ishell] = k;
|
|
||||||
k += lstart[l+1] - lstart[l];
|
|
||||||
size_max = size_max < lstart[l+1] ? lstart[l+1] : size_max;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
ao_index[shell_num] = ao_num+1;
|
|
||||||
}
|
|
||||||
|
|
||||||
/* Don't compute when the radial part is zero. */
|
|
||||||
double cutoff = 27.631021115928547; // -log(1.e-12)
|
|
||||||
|
|
||||||
#ifdef HAVE_OPENMP
|
#ifdef HAVE_OPENMP
|
||||||
#pragma omp parallel
|
#pragma omp parallel
|
||||||
@ -6508,18 +6491,6 @@ qmckl_compute_ao_vgl_hpc_gaussian (
|
|||||||
qmckl_exit_code rc;
|
qmckl_exit_code rc;
|
||||||
double ar2[prim_max] __attribute__((aligned(64)));
|
double ar2[prim_max] __attribute__((aligned(64)));
|
||||||
int32_t powers[3*size_max] __attribute__((aligned(64)));
|
int32_t powers[3*size_max] __attribute__((aligned(64)));
|
||||||
double poly_vgl_l1[4][4] __attribute__((aligned(64))) =
|
|
||||||
{{1.0, 0.0, 0.0, 0.0},
|
|
||||||
{0.0, 1.0, 0.0, 0.0},
|
|
||||||
{0.0, 0.0, 1.0, 0.0},
|
|
||||||
{0.0, 0.0, 0.0, 1.0}};
|
|
||||||
double poly_vgl_l2[5][10]__attribute__((aligned(64))) =
|
|
||||||
{{1., 0., 0., 0., 0., 0., 0., 0., 0., 0.},
|
|
||||||
{0., 1., 0., 0., 0., 0., 0., 0., 0., 0.},
|
|
||||||
{0., 0., 1., 0., 0., 0., 0., 0., 0., 0.},
|
|
||||||
{0., 0., 0., 1., 0., 0., 0., 0., 0., 0.},
|
|
||||||
{0., 0., 0., 0., 2., 0., 0., 2., 0., 2.}};
|
|
||||||
double poly_vgl[5][size_max] __attribute__((aligned(64)));
|
|
||||||
|
|
||||||
double exp_mat[prim_max][8] __attribute__((aligned(64))) ;
|
double exp_mat[prim_max][8] __attribute__((aligned(64))) ;
|
||||||
double ce_mat[shell_max][8] __attribute__((aligned(64))) ;
|
double ce_mat[shell_max][8] __attribute__((aligned(64))) ;
|
||||||
@ -6533,20 +6504,33 @@ qmckl_compute_ao_vgl_hpc_gaussian (
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
double poly_vgl_l1[4][4] __attribute__((aligned(64))) =
|
||||||
|
{{1.0, 0.0, 0.0, 0.0},
|
||||||
|
{0.0, 1.0, 0.0, 0.0},
|
||||||
|
{0.0, 0.0, 1.0, 0.0},
|
||||||
|
{0.0, 0.0, 0.0, 1.0}};
|
||||||
|
double poly_vgl_l2[5][10]__attribute__((aligned(64))) =
|
||||||
|
{{1., 0., 0., 0., 0., 0., 0., 0., 0., 0.},
|
||||||
|
{0., 1., 0., 0., 0., 0., 0., 0., 0., 0.},
|
||||||
|
{0., 0., 1., 0., 0., 0., 0., 0., 0., 0.},
|
||||||
|
{0., 0., 0., 1., 0., 0., 0., 0., 0., 0.},
|
||||||
|
{0., 0., 0., 0., 2., 0., 0., 2., 0., 2.}};
|
||||||
|
double poly_vgl[5][size_max] __attribute__((aligned(64)));
|
||||||
|
|
||||||
#ifdef HAVE_OPENMP
|
#ifdef HAVE_OPENMP
|
||||||
#pragma omp for
|
#pragma omp for
|
||||||
#endif
|
#endif
|
||||||
for (int64_t ipoint=0 ; ipoint < point_num ; ++ipoint) {
|
for (int64_t ipoint=0 ; ipoint < point_num ; ++ipoint) {
|
||||||
const double e_coord[3] __attribute__((aligned(64))) =
|
const double e_coord[3] __attribute__((aligned(64))) =
|
||||||
{ coord[ipoint],
|
{ coord[ipoint],
|
||||||
coord[ipoint + point_num],
|
coord[ipoint + point_num],
|
||||||
coord[ipoint + 2*point_num] };
|
coord[ipoint + 2*point_num] };
|
||||||
|
|
||||||
for (int64_t inucl=0 ; inucl < nucl_num ; ++inucl) {
|
for (int64_t inucl=0 ; inucl < nucl_num ; ++inucl) {
|
||||||
const double n_coord[3] __attribute__((aligned(64))) =
|
const double n_coord[3] __attribute__((aligned(64))) =
|
||||||
{ nucl_coord[inucl],
|
{ nucl_coord[inucl],
|
||||||
nucl_coord[inucl + nucl_num],
|
nucl_coord[inucl + nucl_num],
|
||||||
nucl_coord[inucl + 2*nucl_num] };
|
nucl_coord[inucl + 2*nucl_num] };
|
||||||
|
|
||||||
/* Test if the point is in the range of the nucleus */
|
/* Test if the point is in the range of the nucleus */
|
||||||
const double x = e_coord[0] - n_coord[0];
|
const double x = e_coord[0] - n_coord[0];
|
||||||
@ -6616,6 +6600,7 @@ qmckl_compute_ao_vgl_hpc_gaussian (
|
|||||||
for (int64_t iprim = 0 ; iprim < nidx ; ++iprim) {
|
for (int64_t iprim = 0 ; iprim < nidx ; ++iprim) {
|
||||||
exp_mat[iprim][0] = exp(-ar2[iprim]);
|
exp_mat[iprim][0] = exp(-ar2[iprim]);
|
||||||
}
|
}
|
||||||
|
|
||||||
for (int64_t iprim = 0 ; iprim < nidx ; ++iprim) {
|
for (int64_t iprim = 0 ; iprim < nidx ; ++iprim) {
|
||||||
double f = qmckl_mat(expo_per_nucleus, iprim, inucl) * exp_mat[iprim][0];
|
double f = qmckl_mat(expo_per_nucleus, iprim, inucl) * exp_mat[iprim][0];
|
||||||
f = -f-f;
|
f = -f-f;
|
||||||
@ -6628,21 +6613,6 @@ qmckl_compute_ao_vgl_hpc_gaussian (
|
|||||||
|
|
||||||
/* --- */
|
/* --- */
|
||||||
switch (8) {
|
switch (8) {
|
||||||
case(5):
|
|
||||||
|
|
||||||
for (int i=0 ; i<nucleus_shell_num[inucl] ; ++i) {
|
|
||||||
for (int j=0 ; j<5 ; ++j) {
|
|
||||||
ce_mat[i][j] = 0.;
|
|
||||||
}
|
|
||||||
for (int k=0 ; k<nidx; ++k) {
|
|
||||||
const double cm = coef_mat[inucl][i][k];
|
|
||||||
for (int j=0 ; j<5 ; ++j) {
|
|
||||||
ce_mat[i][j] += cm * exp_mat[k][j];
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
break;
|
|
||||||
|
|
||||||
case(8):
|
case(8):
|
||||||
|
|
||||||
for (int i=0 ; i<nucleus_shell_num[inucl] ; ++i) {
|
for (int i=0 ; i<nucleus_shell_num[inucl] ; ++i) {
|
||||||
@ -6664,6 +6634,21 @@ qmckl_compute_ao_vgl_hpc_gaussian (
|
|||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
|
case(5):
|
||||||
|
|
||||||
|
for (int i=0 ; i<nucleus_shell_num[inucl] ; ++i) {
|
||||||
|
for (int j=0 ; j<5 ; ++j) {
|
||||||
|
ce_mat[i][j] = 0.;
|
||||||
|
}
|
||||||
|
for (int k=0 ; k<nidx; ++k) {
|
||||||
|
const double cm = coef_mat[inucl][i][k];
|
||||||
|
for (int j=0 ; j<5 ; ++j) {
|
||||||
|
ce_mat[i][j] += cm * exp_mat[k][j];
|
||||||
|
}
|
||||||
|
}
|
||||||
|
}
|
||||||
|
break;
|
||||||
|
|
||||||
case(512):
|
case(512):
|
||||||
for(int i=0; i<nucleus_shell_num[inucl]; ++i){
|
for(int i=0; i<nucleus_shell_num[inucl]; ++i){
|
||||||
__m512d cemat_avx512;
|
__m512d cemat_avx512;
|
||||||
@ -6778,14 +6763,14 @@ qmckl_compute_ao_vgl_hpc_gaussian (
|
|||||||
poly_vgl_5 = &(poly_vgl[4][idx]);
|
poly_vgl_5 = &(poly_vgl[4][idx]);
|
||||||
}
|
}
|
||||||
switch (n) {
|
switch (n) {
|
||||||
case(1):
|
case 1:
|
||||||
ao_vgl_1[0] = s1 * f[0];
|
ao_vgl_1[0] = s1 * f[0];
|
||||||
ao_vgl_2[0] = s2 * f[0];
|
ao_vgl_2[0] = s2 * f[0];
|
||||||
ao_vgl_3[0] = s3 * f[0];
|
ao_vgl_3[0] = s3 * f[0];
|
||||||
ao_vgl_4[0] = s4 * f[0];
|
ao_vgl_4[0] = s4 * f[0];
|
||||||
ao_vgl_5[0] = s5;
|
ao_vgl_5[0] = s5;
|
||||||
break;
|
break;
|
||||||
case (3):
|
case 3:
|
||||||
#ifdef HAVE_OPENMP
|
#ifdef HAVE_OPENMP
|
||||||
#pragma omp simd
|
#pragma omp simd
|
||||||
#endif
|
#endif
|
||||||
@ -6800,7 +6785,7 @@ qmckl_compute_ao_vgl_hpc_gaussian (
|
|||||||
poly_vgl_4[il] * s4 )) * f[il];
|
poly_vgl_4[il] * s4 )) * f[il];
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
case(6):
|
case 6:
|
||||||
#ifdef HAVE_OPENMP
|
#ifdef HAVE_OPENMP
|
||||||
#pragma omp simd
|
#pragma omp simd
|
||||||
#endif
|
#endif
|
||||||
|
|||||||
@ -1133,6 +1133,225 @@ qmckl_exit_code test_qmckl_dgemm(qmckl_context context);
|
|||||||
assert(QMCKL_SUCCESS == test_qmckl_dgemm(context));
|
assert(QMCKL_SUCCESS == test_qmckl_dgemm(context));
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
|
** ~qmckl_dgemm_safe~
|
||||||
|
|
||||||
|
"Size-safe" proxy function with the same functionality as ~qmckl_dgemm~
|
||||||
|
but with 3 additional arguments. These arguments ~size_max_M~ (where ~M~ is a matix)
|
||||||
|
are required primarily for the Python API, where compatibility with
|
||||||
|
NumPy arrays implies that sizes of the input and output arrays are provided.
|
||||||
|
|
||||||
|
#+NAME: qmckl_dgemm_safe_args
|
||||||
|
| Variable | Type | In/Out | Description |
|
||||||
|
|--------------+-----------------+--------+---------------------------------------|
|
||||||
|
| ~context~ | ~qmckl_context~ | in | Global state |
|
||||||
|
| ~TransA~ | ~char~ | in | 'T' is transposed |
|
||||||
|
| ~TransB~ | ~char~ | in | 'T' is transposed |
|
||||||
|
| ~m~ | ~int64_t~ | in | Number of rows of the input matrix |
|
||||||
|
| ~n~ | ~int64_t~ | in | Number of columns of the input matrix |
|
||||||
|
| ~k~ | ~int64_t~ | in | Number of columns of the input matrix |
|
||||||
|
| ~alpha~ | ~double~ | in | \alpha |
|
||||||
|
| ~A~ | ~double[][lda]~ | in | Array containing the matrix $A$ |
|
||||||
|
| ~size_max_A~ | ~int64_t~ | in | Size of the matrix $A$ |
|
||||||
|
| ~lda~ | ~int64_t~ | in | Leading dimension of array ~A~ |
|
||||||
|
| ~B~ | ~double[][ldb]~ | in | Array containing the matrix $B$ |
|
||||||
|
| ~size_max_B~ | ~int64_t~ | in | Size of the matrix $B$ |
|
||||||
|
| ~ldb~ | ~int64_t~ | in | Leading dimension of array ~B~ |
|
||||||
|
| ~beta~ | ~double~ | in | \beta |
|
||||||
|
| ~C~ | ~double[][ldc]~ | out | Array containing the matrix $C$ |
|
||||||
|
| ~size_max_C~ | ~int64_t~ | in | Size of the matrix $C$ |
|
||||||
|
| ~ldc~ | ~int64_t~ | in | Leading dimension of array ~C~ |
|
||||||
|
|
||||||
|
Requirements:
|
||||||
|
|
||||||
|
- ~context~ is not ~QMCKL_NULL_CONTEXT~
|
||||||
|
- ~m > 0~
|
||||||
|
- ~n > 0~
|
||||||
|
- ~k > 0~
|
||||||
|
- ~lda >= m~
|
||||||
|
- ~ldb >= n~
|
||||||
|
- ~ldc >= n~
|
||||||
|
- ~A~ is allocated with at least $m \times k \times 8$ bytes
|
||||||
|
- ~B~ is allocated with at least $k \times n \times 8$ bytes
|
||||||
|
- ~C~ is allocated with at least $m \times n \times 8$ bytes
|
||||||
|
- ~size_max_A >= m * k~
|
||||||
|
- ~size_max_B >= k * n~
|
||||||
|
- ~size_max_C >= m * n~
|
||||||
|
|
||||||
|
#+CALL: generate_c_header(table=qmckl_dgemm_safe_args,rettyp="qmckl_exit_code",fname="qmckl_dgemm_safe")
|
||||||
|
|
||||||
|
#+RESULTS:
|
||||||
|
#+BEGIN_src c :tangle (eval h_func) :comments org
|
||||||
|
qmckl_exit_code qmckl_dgemm_safe (
|
||||||
|
const qmckl_context context,
|
||||||
|
const char TransA,
|
||||||
|
const char TransB,
|
||||||
|
const int64_t m,
|
||||||
|
const int64_t n,
|
||||||
|
const int64_t k,
|
||||||
|
const double alpha,
|
||||||
|
const double* A,
|
||||||
|
const int64_t size_max_A,
|
||||||
|
const int64_t lda,
|
||||||
|
const double* B,
|
||||||
|
const int64_t size_max_B,
|
||||||
|
const int64_t ldb,
|
||||||
|
const double beta,
|
||||||
|
double* const C,
|
||||||
|
const int64_t size_max_C,
|
||||||
|
const int64_t ldc );
|
||||||
|
#+END_src
|
||||||
|
|
||||||
|
#+begin_src f90 :tangle (eval f) :exports none
|
||||||
|
integer function qmckl_dgemm_safe_f(context, TransA, TransB, &
|
||||||
|
m, n, k, alpha, A, size_A, LDA, B, size_B, LDB, beta, C, size_C, LDC) &
|
||||||
|
result(info)
|
||||||
|
use qmckl
|
||||||
|
implicit none
|
||||||
|
integer(qmckl_context), intent(in) :: context
|
||||||
|
character , intent(in) :: TransA, TransB
|
||||||
|
integer*8 , intent(in) :: m, n, k
|
||||||
|
double precision , intent(in) :: alpha, beta
|
||||||
|
integer*8 , intent(in) :: lda
|
||||||
|
integer*8 , intent(in) :: size_A
|
||||||
|
double precision , intent(in) :: A(lda,*)
|
||||||
|
integer*8 , intent(in) :: ldb
|
||||||
|
integer*8 , intent(in) :: size_B
|
||||||
|
double precision , intent(in) :: B(ldb,*)
|
||||||
|
integer*8 , intent(in) :: ldc
|
||||||
|
integer*8 , intent(in) :: size_C
|
||||||
|
double precision , intent(out) :: C(ldc,*)
|
||||||
|
|
||||||
|
info = QMCKL_SUCCESS
|
||||||
|
|
||||||
|
if (context == QMCKL_NULL_CONTEXT) then
|
||||||
|
info = QMCKL_INVALID_CONTEXT
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (m <= 0_8) then
|
||||||
|
info = QMCKL_INVALID_ARG_4
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (n <= 0_8) then
|
||||||
|
info = QMCKL_INVALID_ARG_5
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (k <= 0_8) then
|
||||||
|
info = QMCKL_INVALID_ARG_6
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (LDA <= 0) then
|
||||||
|
info = QMCKL_INVALID_ARG_10
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (LDB <= 0) then
|
||||||
|
info = QMCKL_INVALID_ARG_13
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (LDC <= 0) then
|
||||||
|
info = QMCKL_INVALID_ARG_17
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (size_A <= 0) then
|
||||||
|
info = QMCKL_INVALID_ARG_9
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (size_B <= 0) then
|
||||||
|
info = QMCKL_INVALID_ARG_12
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (size_C <= 0) then
|
||||||
|
info = QMCKL_INVALID_ARG_16
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
call dgemm(transA, transB, int(m,4), int(n,4), int(k,4), &
|
||||||
|
alpha, A, int(LDA,4), B, int(LDB,4), beta, C, int(LDC,4))
|
||||||
|
|
||||||
|
end function qmckl_dgemm_safe_f
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
*** C interface :noexport:
|
||||||
|
|
||||||
|
#+CALL: generate_c_interface(table=qmckl_dgemm_safe_args,rettyp="qmckl_exit_code",fname="qmckl_dgemm_safe")
|
||||||
|
|
||||||
|
#+RESULTS:
|
||||||
|
#+begin_src f90 :tangle (eval f) :comments org :exports none
|
||||||
|
integer(c_int32_t) function qmckl_dgemm_safe &
|
||||||
|
(context, TransA, TransB, m, n, k, alpha, A, size_A, lda, B, size_B, ldb, beta, C, size_C, ldc) &
|
||||||
|
bind(C) result(info)
|
||||||
|
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
implicit none
|
||||||
|
|
||||||
|
integer (c_int64_t) , intent(in) , value :: context
|
||||||
|
character , intent(in) , value :: TransA
|
||||||
|
character , intent(in) , value :: TransB
|
||||||
|
integer (c_int64_t) , intent(in) , value :: m
|
||||||
|
integer (c_int64_t) , intent(in) , value :: n
|
||||||
|
integer (c_int64_t) , intent(in) , value :: k
|
||||||
|
real (c_double ) , intent(in) , value :: alpha
|
||||||
|
real (c_double ) , intent(in) :: A(lda,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_A
|
||||||
|
integer (c_int64_t) , intent(in) , value :: lda
|
||||||
|
real (c_double ) , intent(in) :: B(ldb,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_B
|
||||||
|
integer (c_int64_t) , intent(in) , value :: ldb
|
||||||
|
real (c_double ) , intent(in) , value :: beta
|
||||||
|
real (c_double ) , intent(out) :: C(ldc,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_C
|
||||||
|
integer (c_int64_t) , intent(in) , value :: ldc
|
||||||
|
|
||||||
|
integer(c_int32_t), external :: qmckl_dgemm_safe_f
|
||||||
|
info = qmckl_dgemm_safe_f &
|
||||||
|
(context, TransA, TransB, m, n, k, alpha, A, size_A, lda, B, size_B, ldb, beta, C, size_C, ldc)
|
||||||
|
|
||||||
|
end function qmckl_dgemm_safe
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
|
||||||
|
#+CALL: generate_f_interface(table=qmckl_dgemm_safe_args,rettyp="qmckl_exit_code",fname="qmckl_dgemm_safe")
|
||||||
|
|
||||||
|
#+RESULTS:
|
||||||
|
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
||||||
|
interface
|
||||||
|
integer(c_int32_t) function qmckl_dgemm_safe &
|
||||||
|
(context, TransA, TransB, m, n, k, alpha, A, size_A, lda, B, size_B, ldb, beta, C, size_C, ldc) &
|
||||||
|
bind(C)
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
import
|
||||||
|
implicit none
|
||||||
|
|
||||||
|
integer (c_int64_t) , intent(in) , value :: context
|
||||||
|
character , intent(in) , value :: TransA
|
||||||
|
character , intent(in) , value :: TransB
|
||||||
|
integer (c_int64_t) , intent(in) , value :: m
|
||||||
|
integer (c_int64_t) , intent(in) , value :: n
|
||||||
|
integer (c_int64_t) , intent(in) , value :: k
|
||||||
|
real (c_double ) , intent(in) , value :: alpha
|
||||||
|
real (c_double ) , intent(in) :: A(lda,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_A
|
||||||
|
integer (c_int64_t) , intent(in) , value :: lda
|
||||||
|
real (c_double ) , intent(in) :: B(ldb,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_B
|
||||||
|
integer (c_int64_t) , intent(in) , value :: ldb
|
||||||
|
real (c_double ) , intent(in) , value :: beta
|
||||||
|
real (c_double ) , intent(out) :: C(ldc,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_C
|
||||||
|
integer (c_int64_t) , intent(in) , value :: ldc
|
||||||
|
|
||||||
|
end function qmckl_dgemm_safe
|
||||||
|
end interface
|
||||||
|
#+end_src
|
||||||
|
|
||||||
** ~qmckl_matmul~
|
** ~qmckl_matmul~
|
||||||
|
|
||||||
Matrix multiplication using the =qmckl_matrix= data type:
|
Matrix multiplication using the =qmckl_matrix= data type:
|
||||||
@ -1452,8 +1671,6 @@ integer function qmckl_adjugate_f(context, na, A, LDA, B, ldb, det_l) &
|
|||||||
integer*8, intent(in) :: na
|
integer*8, intent(in) :: na
|
||||||
double precision, intent(inout) :: det_l
|
double precision, intent(inout) :: det_l
|
||||||
|
|
||||||
integer :: i,j
|
|
||||||
|
|
||||||
info = QMCKL_SUCCESS
|
info = QMCKL_SUCCESS
|
||||||
|
|
||||||
if (context == QMCKL_NULL_CONTEXT) then
|
if (context == QMCKL_NULL_CONTEXT) then
|
||||||
@ -2274,6 +2491,159 @@ qmckl_exit_code test_qmckl_adjugate(qmckl_context context);
|
|||||||
assert(QMCKL_SUCCESS == test_qmckl_adjugate(context));
|
assert(QMCKL_SUCCESS == test_qmckl_adjugate(context));
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
|
** ~qmckl_adjugate_safe~
|
||||||
|
|
||||||
|
"Size-safe" proxy function with the same functionality as ~qmckl_adjugate~
|
||||||
|
but with 2 additional arguments. These arguments ~size_max_M~ (where ~M~ is a matix)
|
||||||
|
are required primarily for the Python API, where compatibility with
|
||||||
|
NumPy arrays implies that sizes of the input and output arrays are provided.
|
||||||
|
|
||||||
|
|
||||||
|
#+NAME: qmckl_adjugate_safe_args
|
||||||
|
| Variable | Type | In/Out | Description |
|
||||||
|
|--------------+-----------------+--------+------------------------------------------------|
|
||||||
|
| ~context~ | ~qmckl_context~ | in | Global state |
|
||||||
|
| ~n~ | ~int64_t~ | in | Number of rows and columns of the input matrix |
|
||||||
|
| ~A~ | ~double[][lda]~ | in | Array containing the $n \times n$ matrix $A$ |
|
||||||
|
| ~size_max_A~ | ~int64_t~ | in | Size of the matrix $A$ |
|
||||||
|
| ~lda~ | ~int64_t~ | in | Leading dimension of array ~A~ |
|
||||||
|
| ~B~ | ~double[][ldb]~ | out | Adjugate of $A$ |
|
||||||
|
| ~size_max_B~ | ~int64_t~ | in | Size of the matrix $B$ |
|
||||||
|
| ~ldb~ | ~int64_t~ | in | Leading dimension of array ~B~ |
|
||||||
|
| ~det_l~ | ~double~ | inout | determinant of $A$ |
|
||||||
|
|
||||||
|
Requirements:
|
||||||
|
|
||||||
|
- ~context~ is not ~QMCKL_NULL_CONTEXT~
|
||||||
|
- ~n > 0~
|
||||||
|
- ~lda >= m~
|
||||||
|
- ~A~ is allocated with at least $m \times m \times 8$ bytes
|
||||||
|
- ~ldb >= m~
|
||||||
|
- ~B~ is allocated with at least $m \times m \times 8$ bytes
|
||||||
|
- ~size_max_A >= m * m~
|
||||||
|
- ~size_max_B >= m * m~
|
||||||
|
|
||||||
|
#+CALL: generate_c_header(table=qmckl_adjugate_safe_args,rettyp="qmckl_exit_code",fname="qmckl_adjugate_safe")
|
||||||
|
|
||||||
|
#+RESULTS:
|
||||||
|
#+BEGIN_src c :tangle (eval h_func) :comments org
|
||||||
|
qmckl_exit_code qmckl_adjugate_safe (
|
||||||
|
const qmckl_context context,
|
||||||
|
const int64_t n,
|
||||||
|
const double* A,
|
||||||
|
const int64_t size_max_A,
|
||||||
|
const int64_t lda,
|
||||||
|
double* const B,
|
||||||
|
const int64_t size_max_B,
|
||||||
|
const int64_t ldb,
|
||||||
|
double* det_l );
|
||||||
|
#+END_src
|
||||||
|
|
||||||
|
For small matrices (\le 5\times 5), we use specialized routines
|
||||||
|
for performance motivations. For larger sizes, we rely on the
|
||||||
|
LAPACK library.
|
||||||
|
|
||||||
|
#+begin_src f90 :tangle (eval f) :exports none
|
||||||
|
integer function qmckl_adjugate_safe_f(context, &
|
||||||
|
na, A, size_A, LDA, B, size_B, LDB, det_l) &
|
||||||
|
result(info)
|
||||||
|
use qmckl
|
||||||
|
implicit none
|
||||||
|
integer(qmckl_context) , intent(in) :: context
|
||||||
|
double precision, intent(in) :: A (LDA,*)
|
||||||
|
integer*8, intent(in) :: size_A
|
||||||
|
integer*8, intent(in) :: LDA
|
||||||
|
double precision, intent(out) :: B (LDB,*)
|
||||||
|
integer*8, intent(in) :: size_B
|
||||||
|
integer*8, intent(in) :: LDB
|
||||||
|
integer*8, intent(in) :: na
|
||||||
|
double precision, intent(inout) :: det_l
|
||||||
|
|
||||||
|
integer, external :: qmckl_adjugate_f
|
||||||
|
|
||||||
|
info = QMCKL_SUCCESS
|
||||||
|
|
||||||
|
if (size_A < na) then
|
||||||
|
info = QMCKL_INVALID_ARG_4
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (size_B <= 0_8) then
|
||||||
|
info = QMCKL_INVALID_ARG_7
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
info = qmckl_adjugate_f(context, na, A, LDA, B, LDB, det_l)
|
||||||
|
|
||||||
|
if (info == QMCKL_INVALID_ARG_4) then
|
||||||
|
info = QMCKL_INVALID_ARG_5
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
if (info == QMCKL_INVALID_ARG_6) then
|
||||||
|
info = QMCKL_INVALID_ARG_8
|
||||||
|
return
|
||||||
|
endif
|
||||||
|
|
||||||
|
end function qmckl_adjugate_safe_f
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
*** C interface
|
||||||
|
|
||||||
|
#+CALL: generate_c_interface(table=qmckl_adjugate_safe_args,rettyp="qmckl_exit_code",fname="qmckl_adjugate_safe")
|
||||||
|
|
||||||
|
#+RESULTS:
|
||||||
|
#+begin_src f90 :tangle (eval f) :comments org :exports none
|
||||||
|
integer(c_int32_t) function qmckl_adjugate_safe &
|
||||||
|
(context, n, A, size_A, lda, B, size_B, ldb, det_l) &
|
||||||
|
bind(C) result(info)
|
||||||
|
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
implicit none
|
||||||
|
|
||||||
|
integer (c_int64_t) , intent(in) , value :: context
|
||||||
|
integer (c_int64_t) , intent(in) , value :: n
|
||||||
|
real (c_double ) , intent(in) :: A(lda,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: lda
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_A
|
||||||
|
real (c_double ) , intent(out) :: B(ldb,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: ldb
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_B
|
||||||
|
real (c_double ) , intent(inout) :: det_l
|
||||||
|
|
||||||
|
integer(c_int32_t), external :: qmckl_adjugate_safe_f
|
||||||
|
info = qmckl_adjugate_safe_f &
|
||||||
|
(context, n, A, size_A, lda, B, size_B, ldb, det_l)
|
||||||
|
|
||||||
|
end function qmckl_adjugate_safe
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
#+CALL: generate_f_interface(table=qmckl_adjugate_safe_args,rettyp="qmckl_exit_code",fname="qmckl_adjugate_safe")
|
||||||
|
|
||||||
|
#+RESULTS:
|
||||||
|
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
||||||
|
interface
|
||||||
|
integer(c_int32_t) function qmckl_adjugate_safe &
|
||||||
|
(context, n, A, size_A, lda, B, size_B, ldb, det_l) &
|
||||||
|
bind(C)
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
import
|
||||||
|
implicit none
|
||||||
|
|
||||||
|
integer (c_int64_t) , intent(in) , value :: context
|
||||||
|
integer (c_int64_t) , intent(in) , value :: n
|
||||||
|
real (c_double ) , intent(in) :: A(lda,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: lda
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_A
|
||||||
|
real (c_double ) , intent(out) :: B(ldb,*)
|
||||||
|
integer (c_int64_t) , intent(in) , value :: ldb
|
||||||
|
integer (c_int64_t) , intent(in) , value :: size_B
|
||||||
|
real (c_double ) , intent(inout) :: det_l
|
||||||
|
|
||||||
|
end function qmckl_adjugate_safe
|
||||||
|
end interface
|
||||||
|
#+end_src
|
||||||
|
|
||||||
** ~qmckl_transpose~
|
** ~qmckl_transpose~
|
||||||
|
|
||||||
Transposes a matrix: $A^\dagger_{ji} = A_{ij}$.
|
Transposes a matrix: $A^\dagger_{ji} = A_{ij}$.
|
||||||
@ -2373,9 +2743,7 @@ qmckl_transpose (qmckl_context context,
|
|||||||
assert (qmckl_context_destroy(context) == QMCKL_SUCCESS);
|
assert (qmckl_context_destroy(context) == QMCKL_SUCCESS);
|
||||||
return 0;
|
return 0;
|
||||||
}
|
}
|
||||||
|
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
|
|
||||||
# -*- mode: org -*-
|
# -*- mode: org -*-
|
||||||
# vim: syntax=c
|
# vim: syntax=c
|
||||||
|
|||||||
@ -292,6 +292,9 @@ qmckl_context qmckl_context_create() {
|
|||||||
|
|
||||||
rc = qmckl_init_determinant(context);
|
rc = qmckl_init_determinant(context);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
assert (rc == QMCKL_SUCCESS);
|
||||||
|
|
||||||
|
rc = qmckl_init_jastrow(context);
|
||||||
|
assert (rc == QMCKL_SUCCESS);
|
||||||
}
|
}
|
||||||
|
|
||||||
/* Allocate qmckl_memory_struct */
|
/* Allocate qmckl_memory_struct */
|
||||||
|
|||||||
@ -1111,21 +1111,21 @@ const double* nucl_charge = chbrclf_charge;
|
|||||||
const double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
|
const double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
|
||||||
|
|
||||||
rc = qmckl_set_electron_num (context, chbrclf_elec_up_num, chbrclf_elec_dn_num);
|
rc = qmckl_set_electron_num (context, chbrclf_elec_up_num, chbrclf_elec_dn_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
|
rc = qmckl_set_electron_coord (context, 'N', chbrclf_walk_num, elec_coord, chbrclf_walk_num*chbrclf_elec_num*3);
|
||||||
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
assert(qmckl_electron_provided(context));
|
assert(qmckl_electron_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_electron_coord (context, 'N', chbrclf_walk_num, elec_coord, chbrclf_walk_num*chbrclf_elec_num*3);
|
|
||||||
assert(rc == QMCKL_SUCCESS);
|
|
||||||
|
|
||||||
rc = qmckl_set_nucleus_num (context, chbrclf_nucl_num);
|
rc = qmckl_set_nucleus_num (context, chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), chbrclf_nucl_num*3);
|
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), chbrclf_nucl_num*3);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_set_nucleus_charge(context, nucl_charge, chbrclf_nucl_num);
|
rc = qmckl_set_nucleus_charge(context, nucl_charge, chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
assert(qmckl_nucleus_provided(context));
|
assert(qmckl_nucleus_provided(context));
|
||||||
|
|
||||||
@ -1145,57 +1145,57 @@ const char typ = 'G';
|
|||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_type (context, typ);
|
rc = qmckl_set_ao_basis_type (context, typ);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_shell_num (context, chbrclf_shell_num);
|
rc = qmckl_set_ao_basis_shell_num (context, chbrclf_shell_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_prim_num (context, chbrclf_prim_num);
|
rc = qmckl_set_ao_basis_prim_num (context, chbrclf_prim_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_nucleus_index (context, nucleus_index, chbrclf_nucl_num);
|
rc = qmckl_set_ao_basis_nucleus_index (context, nucleus_index, chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_nucleus_shell_num (context, nucleus_shell_num, chbrclf_nucl_num);
|
rc = qmckl_set_ao_basis_nucleus_shell_num (context, nucleus_shell_num, chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_shell_ang_mom (context, shell_ang_mom, chbrclf_shell_num);
|
rc = qmckl_set_ao_basis_shell_ang_mom (context, shell_ang_mom, chbrclf_shell_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_shell_factor (context, shell_factor, chbrclf_shell_num);
|
rc = qmckl_set_ao_basis_shell_factor (context, shell_factor, chbrclf_shell_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_shell_prim_num (context, shell_prim_num, chbrclf_shell_num);
|
rc = qmckl_set_ao_basis_shell_prim_num (context, shell_prim_num, chbrclf_shell_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_shell_prim_index (context, shell_prim_index, chbrclf_shell_num);
|
rc = qmckl_set_ao_basis_shell_prim_index (context, shell_prim_index, chbrclf_shell_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_exponent (context, exponent, chbrclf_prim_num);
|
rc = qmckl_set_ao_basis_exponent (context, exponent, chbrclf_prim_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_coefficient (context, coefficient, chbrclf_prim_num);
|
rc = qmckl_set_ao_basis_coefficient (context, coefficient, chbrclf_prim_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_ao_basis_provided(context));
|
assert(!qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_prim_factor (context, prim_factor, chbrclf_prim_num);
|
rc = qmckl_set_ao_basis_prim_factor (context, prim_factor, chbrclf_prim_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_ao_num(context, chbrclf_ao_num);
|
rc = qmckl_set_ao_basis_ao_num(context, chbrclf_ao_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_set_ao_basis_ao_factor (context, ao_factor, chbrclf_ao_num);
|
rc = qmckl_set_ao_basis_ao_factor (context, ao_factor, chbrclf_ao_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
assert(qmckl_ao_basis_provided(context));
|
assert(qmckl_ao_basis_provided(context));
|
||||||
|
|
||||||
@ -1203,22 +1203,22 @@ assert(qmckl_ao_basis_provided(context));
|
|||||||
double ao_vgl[chbrclf_walk_num*chbrclf_elec_num][5][chbrclf_ao_num];
|
double ao_vgl[chbrclf_walk_num*chbrclf_elec_num][5][chbrclf_ao_num];
|
||||||
|
|
||||||
rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0]), (int64_t) 5*chbrclf_walk_num*chbrclf_elec_num*chbrclf_ao_num);
|
rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0]), (int64_t) 5*chbrclf_walk_num*chbrclf_elec_num*chbrclf_ao_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
/* Set up MO data */
|
/* Set up MO data */
|
||||||
rc = qmckl_set_mo_basis_mo_num(context, chbrclf_mo_num);
|
rc = qmckl_set_mo_basis_mo_num(context, chbrclf_mo_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
const double * mo_coefficient = &(chbrclf_mo_coef[0]);
|
const double * mo_coefficient = &(chbrclf_mo_coef[0]);
|
||||||
|
|
||||||
rc = qmckl_set_mo_basis_coefficient(context, mo_coefficient);
|
rc = qmckl_set_mo_basis_coefficient(context, mo_coefficient);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
assert(qmckl_mo_basis_provided(context));
|
assert(qmckl_mo_basis_provided(context));
|
||||||
|
|
||||||
double mo_vgl[chbrclf_walk_num*chbrclf_elec_num][5][chbrclf_mo_num];
|
double mo_vgl[chbrclf_walk_num*chbrclf_elec_num][5][chbrclf_mo_num];
|
||||||
rc = qmckl_get_mo_basis_mo_vgl(context, &(mo_vgl[0][0][0]), 5*chbrclf_walk_num*chbrclf_elec_num*chbrclf_mo_num);
|
rc = qmckl_get_mo_basis_mo_vgl(context, &(mo_vgl[0][0][0]), 5*chbrclf_walk_num*chbrclf_elec_num*chbrclf_mo_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
/* Set up determinant data */
|
/* Set up determinant data */
|
||||||
|
|
||||||
@ -1238,19 +1238,19 @@ for(k = 0; k < det_num_beta; ++k)
|
|||||||
mo_index_beta[k][i][j] = j + 1;
|
mo_index_beta[k][i][j] = j + 1;
|
||||||
|
|
||||||
rc = qmckl_set_determinant_type (context, typ);
|
rc = qmckl_set_determinant_type (context, typ);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_set_determinant_det_num_alpha (context, det_num_alpha);
|
rc = qmckl_set_determinant_det_num_alpha (context, det_num_alpha);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_set_determinant_det_num_beta (context, det_num_beta);
|
rc = qmckl_set_determinant_det_num_beta (context, det_num_beta);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_set_determinant_mo_index_alpha (context, &(mo_index_alpha[0][0][0]));
|
rc = qmckl_set_determinant_mo_index_alpha (context, &(mo_index_alpha[0][0][0]));
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_set_determinant_mo_index_beta (context, &(mo_index_beta[0][0][0]));
|
rc = qmckl_set_determinant_mo_index_beta (context, &(mo_index_beta[0][0][0]));
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
// Get slater-determinant
|
// Get slater-determinant
|
||||||
|
|
||||||
@ -1258,10 +1258,10 @@ double det_vgl_alpha[det_num_alpha][chbrclf_walk_num][5][chbrclf_elec_up_num][ch
|
|||||||
double det_vgl_beta[det_num_beta][chbrclf_walk_num][5][chbrclf_elec_dn_num][chbrclf_elec_dn_num];
|
double det_vgl_beta[det_num_beta][chbrclf_walk_num][5][chbrclf_elec_dn_num][chbrclf_elec_dn_num];
|
||||||
|
|
||||||
rc = qmckl_get_det_vgl_alpha(context, &(det_vgl_alpha[0][0][0][0][0]));
|
rc = qmckl_get_det_vgl_alpha(context, &(det_vgl_alpha[0][0][0][0][0]));
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_get_det_vgl_beta(context, &(det_vgl_beta[0][0][0][0][0]));
|
rc = qmckl_get_det_vgl_beta(context, &(det_vgl_beta[0][0][0][0][0]));
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
@ -2018,10 +2018,10 @@ double det_inv_matrix_alpha[det_num_alpha][chbrclf_walk_num][chbrclf_elec_up_num
|
|||||||
double det_inv_matrix_beta[det_num_beta][chbrclf_walk_num][chbrclf_elec_dn_num][chbrclf_elec_dn_num];
|
double det_inv_matrix_beta[det_num_beta][chbrclf_walk_num][chbrclf_elec_dn_num][chbrclf_elec_dn_num];
|
||||||
|
|
||||||
rc = qmckl_get_det_inv_matrix_alpha(context, &(det_inv_matrix_alpha[0][0][0][0]));
|
rc = qmckl_get_det_inv_matrix_alpha(context, &(det_inv_matrix_alpha[0][0][0][0]));
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_get_det_inv_matrix_beta(context, &(det_inv_matrix_beta[0][0][0][0]));
|
rc = qmckl_get_det_inv_matrix_beta(context, &(det_inv_matrix_beta[0][0][0][0]));
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
@ -2034,7 +2034,7 @@ assert (rc == QMCKL_SUCCESS);
|
|||||||
*** Test
|
*** Test
|
||||||
#+begin_src c :tangle (eval c_test)
|
#+begin_src c :tangle (eval c_test)
|
||||||
rc = qmckl_context_destroy(context);
|
rc = qmckl_context_destroy(context);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
return 0;
|
return 0;
|
||||||
}
|
}
|
||||||
|
|||||||
File diff suppressed because it is too large
Load Diff
@ -608,6 +608,62 @@ qmckl_last_error(qmckl_context context, char* buffer) {
|
|||||||
end interface
|
end interface
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
|
* Helper functions for debugging
|
||||||
|
|
||||||
|
The following function prints to ~stderr~ an error message is the return code is
|
||||||
|
not ~QMCKL_SUCCESS~.
|
||||||
|
|
||||||
|
# Header
|
||||||
|
#+begin_src c :comments org :tangle (eval h_func)
|
||||||
|
qmckl_exit_code
|
||||||
|
qmckl_check(qmckl_context context, qmckl_exit_code rc);
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
# Source
|
||||||
|
#+begin_src c :tangle (eval c) :exports none
|
||||||
|
#include <stdio.h>
|
||||||
|
|
||||||
|
qmckl_exit_code
|
||||||
|
qmckl_check(qmckl_context context, qmckl_exit_code rc)
|
||||||
|
{
|
||||||
|
|
||||||
|
char fname[QMCKL_MAX_FUN_LEN];
|
||||||
|
char message[QMCKL_MAX_MSG_LEN];
|
||||||
|
|
||||||
|
if (rc != QMCKL_SUCCESS) {
|
||||||
|
fprintf(stderr, "===========\nQMCKL ERROR\n%s\n", qmckl_string_of_error(rc));
|
||||||
|
qmckl_get_error(context, &rc, fname, message);
|
||||||
|
fprintf(stderr, "Function: %s\nMessage: %s\n===========\n", fname, message);
|
||||||
|
}
|
||||||
|
|
||||||
|
return rc;
|
||||||
|
}
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
It should be used as:
|
||||||
|
|
||||||
|
#+begin_src c
|
||||||
|
rc = qmckl_check(context,
|
||||||
|
qmckl_...(context, ...)
|
||||||
|
);
|
||||||
|
assert (rc == QMCKL_SUCCESS);
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
** Fortran inteface
|
||||||
|
|
||||||
|
#+begin_src f90 :tangle (eval fh_func) :exports none :noweb yes
|
||||||
|
interface
|
||||||
|
function qmckl_check (context, rc) bind(C, name='qmckl_check')
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
import
|
||||||
|
implicit none
|
||||||
|
integer(qmckl_exit_code) :: qmckl_check
|
||||||
|
integer (c_int64_t) , intent(in), value :: context
|
||||||
|
integer(qmckl_exit_code), intent(in) :: rc
|
||||||
|
end function qmckl_check
|
||||||
|
end interface
|
||||||
|
#+end_src
|
||||||
|
|
||||||
* End of files :noexport:
|
* End of files :noexport:
|
||||||
|
|
||||||
#+begin_src c :comments link :tangle (eval h_private_type)
|
#+begin_src c :comments link :tangle (eval h_private_type)
|
||||||
@ -623,12 +679,13 @@ qmckl_last_error(qmckl_context context, char* buffer) {
|
|||||||
qmckl_exit_code exit_code;
|
qmckl_exit_code exit_code;
|
||||||
exit_code = 1;
|
exit_code = 1;
|
||||||
|
|
||||||
assert (qmckl_set_error(context, exit_code, "qmckl_transpose", "Success") == QMCKL_SUCCESS);
|
assert (qmckl_set_error(context, exit_code, "my_function", "Message") == QMCKL_SUCCESS);
|
||||||
|
|
||||||
assert (qmckl_get_error(context, &exit_code, function_name, message) == QMCKL_SUCCESS);
|
assert (qmckl_get_error(context, &exit_code, function_name, message) == QMCKL_SUCCESS);
|
||||||
assert (exit_code == 1);
|
assert (exit_code == 1);
|
||||||
assert (strcmp(function_name,"qmckl_transpose") == 0);
|
|
||||||
assert (strcmp(message,"Success") == 0);
|
assert (strcmp(function_name,"my_function") == 0);
|
||||||
|
assert (strcmp(message,"Message") == 0);
|
||||||
|
|
||||||
exit_code = qmckl_context_destroy(context);
|
exit_code = qmckl_context_destroy(context);
|
||||||
assert(exit_code == QMCKL_SUCCESS);
|
assert(exit_code == QMCKL_SUCCESS);
|
||||||
|
|||||||
File diff suppressed because it is too large
Load Diff
507
org/qmckl_mo.org
507
org/qmckl_mo.org
@ -3,16 +3,16 @@
|
|||||||
#+INCLUDE: ../tools/lib.org
|
#+INCLUDE: ../tools/lib.org
|
||||||
|
|
||||||
The molecular orbitals (MOs) are defined in the basis of AOs along with a AO to MO
|
The molecular orbitals (MOs) are defined in the basis of AOs along with a AO to MO
|
||||||
coefficient matrix \[C\]. Using these coefficients (e.g. from Hartree Fock SCF method)
|
coefficient matrix $C$. Using these coefficients (e.g. from Hartree Fock method)
|
||||||
the MOs are defined as follows:
|
the MOs are defined as follows:
|
||||||
|
|
||||||
\[
|
\[
|
||||||
\phi_i(\mathbf{r}) = C_i * \chi_i (\mathbf{r})
|
\phi_i(\mathbf{r}) = C_i * \chi_i (\mathbf{r})
|
||||||
\]
|
\]
|
||||||
|
|
||||||
|
By default, all the MOs are computed. A subset of MOs can be selected
|
||||||
In this section we demonstrate how to use the QMCkl specific DGEMM
|
for computation, for example to remove computation of the useless
|
||||||
function to calculate the MOs.
|
virtual MOs present in a MO coefficient matrix.
|
||||||
|
|
||||||
|
|
||||||
* Headers :noexport:
|
* Headers :noexport:
|
||||||
@ -84,11 +84,11 @@ int main() {
|
|||||||
The following arrays are stored in the context:
|
The following arrays are stored in the context:
|
||||||
|
|
||||||
|
|
||||||
|-----------------+--------------------+----------------------------------------|
|
|-------------------+--------------------+----------------------------------------|
|
||||||
| ~mo_num~ | | Number of MOs |
|
| ~mo_num~ | | Number of MOs |
|
||||||
| ~coefficient~ | ~[mo_num][ao_num]~ | Orbital coefficients |
|
| ~coefficient~ | ~[mo_num][ao_num]~ | MO coefficients |
|
||||||
| ~coefficient_t~ | ~[ao_num][mo_num]~ | Transposed of the Orbital coefficients |
|
| ~coefficient_t~ | ~[ao_num][mo_num]~ | Transposed of the Orbital coefficients |
|
||||||
|-----------------+--------------------+----------------------------------------|
|
|-------------------+--------------------+----------------------------------------|
|
||||||
|
|
||||||
Computed data:
|
Computed data:
|
||||||
|
|
||||||
@ -103,7 +103,7 @@ int main() {
|
|||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_private_type)
|
#+begin_src c :comments org :tangle (eval h_private_type)
|
||||||
typedef struct qmckl_mo_basis_struct {
|
typedef struct qmckl_mo_basis_struct {
|
||||||
int64_t mo_num;
|
int64_t mo_num;
|
||||||
double * restrict coefficient;
|
double * restrict coefficient;
|
||||||
double * restrict coefficient_t;
|
double * restrict coefficient_t;
|
||||||
|
|
||||||
@ -144,150 +144,6 @@ qmckl_exit_code qmckl_init_mo_basis(qmckl_context context) {
|
|||||||
}
|
}
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
** Access functions
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_get_mo_basis_mo_num (const qmckl_context context,
|
|
||||||
int64_t* mo_num);
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_get_mo_basis_mo_num (const qmckl_context context,
|
|
||||||
int64_t* mo_num)
|
|
||||||
{
|
|
||||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_INVALID_CONTEXT,
|
|
||||||
"qmckl_get_mo_basis_mo_num",
|
|
||||||
NULL);
|
|
||||||
}
|
|
||||||
|
|
||||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
|
||||||
assert (ctx != NULL);
|
|
||||||
|
|
||||||
int32_t mask = 1;
|
|
||||||
|
|
||||||
if ( (ctx->mo_basis.uninitialized & mask) != 0) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_NOT_PROVIDED,
|
|
||||||
"qmckl_get_mo_basis_mo_num",
|
|
||||||
NULL);
|
|
||||||
}
|
|
||||||
|
|
||||||
assert (ctx->mo_basis.mo_num > (int64_t) 0);
|
|
||||||
,*mo_num = ctx->mo_basis.mo_num;
|
|
||||||
return QMCKL_SUCCESS;
|
|
||||||
}
|
|
||||||
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_get_mo_basis_coefficient (const qmckl_context context,
|
|
||||||
double* const coefficient,
|
|
||||||
const int64_t size_max);
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :exports none
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_get_mo_basis_coefficient (const qmckl_context context,
|
|
||||||
double* const coefficient,
|
|
||||||
const int64_t size_max)
|
|
||||||
{
|
|
||||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_INVALID_CONTEXT,
|
|
||||||
"qmckl_get_mo_basis_coefficient",
|
|
||||||
NULL);
|
|
||||||
}
|
|
||||||
|
|
||||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
|
||||||
assert (ctx != NULL);
|
|
||||||
|
|
||||||
int32_t mask = 1 << 1;
|
|
||||||
|
|
||||||
if ( (ctx->ao_basis.uninitialized & mask) != 0) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_NOT_PROVIDED,
|
|
||||||
"qmckl_get_mo_basis_coefficient",
|
|
||||||
NULL);
|
|
||||||
}
|
|
||||||
|
|
||||||
if (coefficient == NULL) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_INVALID_ARG_2,
|
|
||||||
"qmckl_get_mo_basis_coefficient",
|
|
||||||
"NULL pointer");
|
|
||||||
}
|
|
||||||
|
|
||||||
if (size_max < ctx->ao_basis.ao_num * ctx->mo_basis.mo_num) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_INVALID_ARG_3,
|
|
||||||
"qmckl_get_mo_basis_coefficient",
|
|
||||||
"Array too small. Expected mo_num * ao_num");
|
|
||||||
}
|
|
||||||
|
|
||||||
assert (ctx->mo_basis.coefficient != NULL);
|
|
||||||
memcpy(coefficient, ctx->mo_basis.coefficient,
|
|
||||||
ctx->ao_basis.ao_num * ctx->mo_basis.mo_num * sizeof(double));
|
|
||||||
|
|
||||||
return QMCKL_SUCCESS;
|
|
||||||
}
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
When all the data for the AOs have been provided, the following
|
|
||||||
function returns ~true~.
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_func)
|
|
||||||
bool qmckl_mo_basis_provided (const qmckl_context context);
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :exports none
|
|
||||||
bool qmckl_mo_basis_provided(const qmckl_context context) {
|
|
||||||
|
|
||||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
|
||||||
return false;
|
|
||||||
}
|
|
||||||
|
|
||||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
|
||||||
assert (ctx != NULL);
|
|
||||||
|
|
||||||
return ctx->mo_basis.provided;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
*** Fortran interfaces
|
|
||||||
|
|
||||||
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
|
||||||
interface
|
|
||||||
integer(c_int32_t) function qmckl_get_mo_basis_mo_num (context, &
|
|
||||||
mo_num) bind(C)
|
|
||||||
use, intrinsic :: iso_c_binding
|
|
||||||
import
|
|
||||||
implicit none
|
|
||||||
integer (c_int64_t) , intent(in) , value :: context
|
|
||||||
integer (c_int64_t) , intent(out) :: mo_num
|
|
||||||
end function qmckl_get_mo_basis_mo_num
|
|
||||||
end interface
|
|
||||||
|
|
||||||
interface
|
|
||||||
integer(c_int32_t) function qmckl_get_mo_basis_coefficient(context, &
|
|
||||||
coefficient, size_max) bind(C)
|
|
||||||
use, intrinsic :: iso_c_binding
|
|
||||||
import
|
|
||||||
implicit none
|
|
||||||
integer (c_int64_t) , intent(in) , value :: context
|
|
||||||
double precision, intent(out) :: coefficient(*)
|
|
||||||
integer (c_int64_t) , intent(in), value :: size_max
|
|
||||||
end function qmckl_get_mo_basis_coefficient
|
|
||||||
end interface
|
|
||||||
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
** Initialization functions
|
** Initialization functions
|
||||||
|
|
||||||
To set the basis set, all the following functions need to be
|
To set the basis set, all the following functions need to be
|
||||||
@ -398,6 +254,15 @@ qmckl_exit_code qmckl_finalize_mo_basis(qmckl_context context) {
|
|||||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||||
assert (ctx != NULL);
|
assert (ctx != NULL);
|
||||||
|
|
||||||
|
if (ctx->mo_basis.coefficient_t != NULL) {
|
||||||
|
qmckl_exit_code rc = qmckl_free(context, ctx->mo_basis.coefficient_t);
|
||||||
|
if (rc != QMCKL_SUCCESS) {
|
||||||
|
return qmckl_failwith( context, rc,
|
||||||
|
"qmckl_finalize_mo_basis",
|
||||||
|
NULL);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
|
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
|
||||||
mem_info.size = ctx->ao_basis.ao_num * ctx->mo_basis.mo_num * sizeof(double);
|
mem_info.size = ctx->ao_basis.ao_num * ctx->mo_basis.mo_num * sizeof(double);
|
||||||
double* new_array = (double*) qmckl_malloc(context, mem_info);
|
double* new_array = (double*) qmckl_malloc(context, mem_info);
|
||||||
@ -410,15 +275,6 @@ qmckl_exit_code qmckl_finalize_mo_basis(qmckl_context context) {
|
|||||||
|
|
||||||
assert (ctx->mo_basis.coefficient != NULL);
|
assert (ctx->mo_basis.coefficient != NULL);
|
||||||
|
|
||||||
if (ctx->mo_basis.coefficient_t != NULL) {
|
|
||||||
qmckl_exit_code rc = qmckl_free(context, ctx->mo_basis.coefficient);
|
|
||||||
if (rc != QMCKL_SUCCESS) {
|
|
||||||
return qmckl_failwith( context, rc,
|
|
||||||
"qmckl_finalize_mo_basis",
|
|
||||||
NULL);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
for (int64_t i=0 ; i<ctx->ao_basis.ao_num ; ++i) {
|
for (int64_t i=0 ; i<ctx->ao_basis.ao_num ; ++i) {
|
||||||
for (int64_t j=0 ; j<ctx->mo_basis.mo_num ; ++j) {
|
for (int64_t j=0 ; j<ctx->mo_basis.mo_num ; ++j) {
|
||||||
new_array[i*ctx->mo_basis.mo_num + j] = ctx->mo_basis.coefficient[j*ctx->ao_basis.ao_num + i];
|
new_array[i*ctx->mo_basis.mo_num + j] = ctx->mo_basis.coefficient[j*ctx->ao_basis.ao_num + i];
|
||||||
@ -426,11 +282,269 @@ qmckl_exit_code qmckl_finalize_mo_basis(qmckl_context context) {
|
|||||||
}
|
}
|
||||||
|
|
||||||
ctx->mo_basis.coefficient_t = new_array;
|
ctx->mo_basis.coefficient_t = new_array;
|
||||||
qmckl_exit_code rc = QMCKL_SUCCESS;
|
|
||||||
return rc;
|
qmckl_exit_code rc;
|
||||||
|
if (ctx->mo_basis.mo_vgl != NULL) {
|
||||||
|
rc = qmckl_free(context, ctx->mo_basis.mo_vgl);
|
||||||
|
if (rc != QMCKL_SUCCESS) return rc;
|
||||||
|
ctx->mo_basis.mo_vgl = NULL;
|
||||||
|
ctx->mo_basis.mo_vgl_date = 0;
|
||||||
|
}
|
||||||
|
|
||||||
|
if (ctx->mo_basis.mo_value != NULL) {
|
||||||
|
rc = qmckl_free(context, ctx->mo_basis.mo_value);
|
||||||
|
if (rc != QMCKL_SUCCESS) return rc;
|
||||||
|
ctx->mo_basis.mo_value = NULL;
|
||||||
|
ctx->mo_basis.mo_value_date = 0;
|
||||||
|
}
|
||||||
|
|
||||||
|
return QMCKL_SUCCESS;
|
||||||
}
|
}
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
|
** Access functions
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
||||||
|
qmckl_exit_code
|
||||||
|
qmckl_get_mo_basis_mo_num (const qmckl_context context,
|
||||||
|
int64_t* mo_num);
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
||||||
|
qmckl_exit_code
|
||||||
|
qmckl_get_mo_basis_mo_num (const qmckl_context context,
|
||||||
|
int64_t* mo_num)
|
||||||
|
{
|
||||||
|
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_INVALID_CONTEXT,
|
||||||
|
"qmckl_get_mo_basis_mo_num",
|
||||||
|
NULL);
|
||||||
|
}
|
||||||
|
|
||||||
|
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||||
|
assert (ctx != NULL);
|
||||||
|
|
||||||
|
int32_t mask = 1;
|
||||||
|
|
||||||
|
if ( (ctx->mo_basis.uninitialized & mask) != 0) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_NOT_PROVIDED,
|
||||||
|
"qmckl_get_mo_basis_mo_num",
|
||||||
|
NULL);
|
||||||
|
}
|
||||||
|
|
||||||
|
assert (ctx->mo_basis.mo_num > (int64_t) 0);
|
||||||
|
,*mo_num = ctx->mo_basis.mo_num;
|
||||||
|
return QMCKL_SUCCESS;
|
||||||
|
}
|
||||||
|
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
||||||
|
qmckl_exit_code
|
||||||
|
qmckl_get_mo_basis_coefficient (const qmckl_context context,
|
||||||
|
double* const coefficient,
|
||||||
|
const int64_t size_max);
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval c) :exports none
|
||||||
|
qmckl_exit_code
|
||||||
|
qmckl_get_mo_basis_coefficient (const qmckl_context context,
|
||||||
|
double* const coefficient,
|
||||||
|
const int64_t size_max)
|
||||||
|
{
|
||||||
|
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_INVALID_CONTEXT,
|
||||||
|
"qmckl_get_mo_basis_coefficient",
|
||||||
|
NULL);
|
||||||
|
}
|
||||||
|
|
||||||
|
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||||
|
assert (ctx != NULL);
|
||||||
|
|
||||||
|
int32_t mask = 1 << 1;
|
||||||
|
|
||||||
|
if ( (ctx->mo_basis.uninitialized & mask) != 0) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_NOT_PROVIDED,
|
||||||
|
"qmckl_get_mo_basis_coefficient",
|
||||||
|
NULL);
|
||||||
|
}
|
||||||
|
|
||||||
|
if (coefficient == NULL) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_INVALID_ARG_2,
|
||||||
|
"qmckl_get_mo_basis_coefficient",
|
||||||
|
"NULL pointer");
|
||||||
|
}
|
||||||
|
|
||||||
|
if (size_max < ctx->ao_basis.ao_num * ctx->mo_basis.mo_num) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_INVALID_ARG_3,
|
||||||
|
"qmckl_get_mo_basis_coefficient",
|
||||||
|
"Array too small: expected mo_num * ao_num.");
|
||||||
|
}
|
||||||
|
|
||||||
|
assert (ctx->mo_basis.coefficient != NULL);
|
||||||
|
memcpy(coefficient, ctx->mo_basis.coefficient,
|
||||||
|
ctx->ao_basis.ao_num * ctx->mo_basis.mo_num * sizeof(double));
|
||||||
|
|
||||||
|
return QMCKL_SUCCESS;
|
||||||
|
}
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
When all the data for the AOs have been provided, the following
|
||||||
|
function returns ~true~.
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval h_func)
|
||||||
|
bool qmckl_mo_basis_provided (const qmckl_context context);
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval c) :exports none
|
||||||
|
bool qmckl_mo_basis_provided(const qmckl_context context) {
|
||||||
|
|
||||||
|
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||||
|
return false;
|
||||||
|
}
|
||||||
|
|
||||||
|
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||||
|
assert (ctx != NULL);
|
||||||
|
|
||||||
|
return ctx->mo_basis.provided;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
*** Fortran interfaces
|
||||||
|
|
||||||
|
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
||||||
|
interface
|
||||||
|
integer(c_int32_t) function qmckl_get_mo_basis_mo_num (context, &
|
||||||
|
mo_num) bind(C)
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
import
|
||||||
|
implicit none
|
||||||
|
integer (c_int64_t) , intent(in) , value :: context
|
||||||
|
integer (c_int64_t) , intent(out) :: mo_num
|
||||||
|
end function qmckl_get_mo_basis_mo_num
|
||||||
|
end interface
|
||||||
|
|
||||||
|
interface
|
||||||
|
integer(c_int32_t) function qmckl_get_mo_basis_coefficient(context, &
|
||||||
|
coefficient, size_max) bind(C)
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
import
|
||||||
|
implicit none
|
||||||
|
integer (c_int64_t) , intent(in) , value :: context
|
||||||
|
double precision, intent(out) :: coefficient(*)
|
||||||
|
integer (c_int64_t) , intent(in), value :: size_max
|
||||||
|
end function qmckl_get_mo_basis_coefficient
|
||||||
|
end interface
|
||||||
|
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
** Update
|
||||||
|
|
||||||
|
Useless MOs can be removed, for instance virtual MOs in a single
|
||||||
|
determinant calculation.
|
||||||
|
|
||||||
|
To select a subset of MOs that will be kept, create an array of
|
||||||
|
integers of size =mo_num=. If the integer is zero, the MO is dropped,
|
||||||
|
otherwise it is kept.
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval h_func)
|
||||||
|
bool qmckl_mo_basis_select_mo (const qmckl_context context,
|
||||||
|
const int32_t* keep,
|
||||||
|
const int64_t size_max);
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
#+begin_src c :comments org :tangle (eval c) :exports none
|
||||||
|
bool qmckl_mo_basis_select_mo (const qmckl_context context,
|
||||||
|
const int32_t* keep,
|
||||||
|
const int64_t size_max)
|
||||||
|
{
|
||||||
|
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_INVALID_CONTEXT,
|
||||||
|
"qmckl_get_mo_basis_select_mo",
|
||||||
|
NULL);
|
||||||
|
}
|
||||||
|
|
||||||
|
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
||||||
|
assert (ctx != NULL);
|
||||||
|
|
||||||
|
if ( !(qmckl_mo_basis_provided(context)) ) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_NOT_PROVIDED,
|
||||||
|
"qmckl_get_mo_basis_select_mo",
|
||||||
|
NULL);
|
||||||
|
}
|
||||||
|
|
||||||
|
if (keep == NULL) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_INVALID_ARG_2,
|
||||||
|
"qmckl_get_mo_basis_select_mo",
|
||||||
|
"NULL pointer");
|
||||||
|
}
|
||||||
|
|
||||||
|
const int64_t mo_num = ctx->mo_basis.mo_num;
|
||||||
|
const int64_t ao_num = ctx->ao_basis.ao_num;
|
||||||
|
|
||||||
|
if (size_max < mo_num) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_INVALID_ARG_3,
|
||||||
|
"qmckl_get_mo_basis_select_mo",
|
||||||
|
"Array too small: expected mo_num.");
|
||||||
|
}
|
||||||
|
|
||||||
|
int64_t mo_num_new = 0;
|
||||||
|
for (int64_t i=0 ; i<mo_num ; ++i) {
|
||||||
|
if (keep[i] != 0) ++mo_num_new;
|
||||||
|
}
|
||||||
|
|
||||||
|
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
|
||||||
|
mem_info.size = ao_num * mo_num_new * sizeof(double);
|
||||||
|
double* restrict coefficient = (double*) qmckl_malloc(context, mem_info);
|
||||||
|
|
||||||
|
int64_t k = 0;
|
||||||
|
for (int64_t i=0 ; i<mo_num ; ++i) {
|
||||||
|
if (keep[i] != 0) {
|
||||||
|
memcpy(&(coefficient[k*ao_num]), &(ctx->mo_basis.coefficient[i*ao_num]), ao_num*sizeof(double));
|
||||||
|
++k;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
qmckl_exit_code rc = qmckl_free(context, ctx->mo_basis.coefficient);
|
||||||
|
if (rc != QMCKL_SUCCESS) return rc;
|
||||||
|
|
||||||
|
ctx->mo_basis.coefficient = coefficient;
|
||||||
|
ctx->mo_basis.mo_num = mo_num_new;
|
||||||
|
|
||||||
|
rc = qmckl_finalize_mo_basis(context);
|
||||||
|
return rc;
|
||||||
|
}
|
||||||
|
|
||||||
|
#+end_src
|
||||||
|
|
||||||
|
*** Fortran interface
|
||||||
|
|
||||||
|
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
||||||
|
interface
|
||||||
|
integer(c_int32_t) function qmckl_mo_basis_select_mo (context, &
|
||||||
|
keep, size_max) bind(C)
|
||||||
|
use, intrinsic :: iso_c_binding
|
||||||
|
import
|
||||||
|
implicit none
|
||||||
|
integer (c_int64_t) , intent(in), value :: context
|
||||||
|
integer (c_int32_t) , intent(in) :: keep(*)
|
||||||
|
integer (c_int64_t) , intent(in), value :: size_max
|
||||||
|
end function qmckl_mo_basis_select_mo
|
||||||
|
end interface
|
||||||
|
#+end_src
|
||||||
|
|
||||||
* Computation
|
* Computation
|
||||||
|
|
||||||
** Computation of MOs: values only
|
** Computation of MOs: values only
|
||||||
@ -628,15 +742,8 @@ qmckl_exit_code qmckl_provide_mo_basis_mo_value(qmckl_context context)
|
|||||||
|
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
||||||
rc = qmckl_provide_ao_basis_ao_value(context);
|
|
||||||
if (rc != QMCKL_SUCCESS) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_NOT_PROVIDED,
|
|
||||||
"qmckl_ao_value",
|
|
||||||
NULL);
|
|
||||||
}
|
|
||||||
|
|
||||||
if (ctx->mo_basis.mo_vgl_date == ctx->point.date) {
|
if (ctx->mo_basis.mo_vgl_date == ctx->point.date) {
|
||||||
|
|
||||||
// mo_vgl has been computed at this step: Just copy the data.
|
// mo_vgl has been computed at this step: Just copy the data.
|
||||||
@ -653,6 +760,14 @@ qmckl_exit_code qmckl_provide_mo_basis_mo_value(qmckl_context context)
|
|||||||
|
|
||||||
} else {
|
} else {
|
||||||
|
|
||||||
|
rc = qmckl_provide_ao_basis_ao_value(context);
|
||||||
|
if (rc != QMCKL_SUCCESS) {
|
||||||
|
return qmckl_failwith( context,
|
||||||
|
QMCKL_NOT_PROVIDED,
|
||||||
|
"qmckl_ao_value",
|
||||||
|
NULL);
|
||||||
|
}
|
||||||
|
|
||||||
rc = qmckl_compute_mo_basis_mo_value(context,
|
rc = qmckl_compute_mo_basis_mo_value(context,
|
||||||
ctx->ao_basis.ao_num,
|
ctx->ao_basis.ao_num,
|
||||||
ctx->mo_basis.mo_num,
|
ctx->mo_basis.mo_num,
|
||||||
@ -664,6 +779,8 @@ qmckl_exit_code qmckl_provide_mo_basis_mo_value(qmckl_context context)
|
|||||||
}
|
}
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
#+CALL: write_provider_post( group="mo_basis", data="mo_value" )
|
#+CALL: write_provider_post( group="mo_basis", data="mo_value" )
|
||||||
|
|
||||||
#+RESULTS:
|
#+RESULTS:
|
||||||
@ -727,8 +844,8 @@ integer function qmckl_compute_mo_basis_mo_value_doc_f(context, &
|
|||||||
do j=1,point_num
|
do j=1,point_num
|
||||||
mo_value(:,j) = 0.d0
|
mo_value(:,j) = 0.d0
|
||||||
do k=1,ao_num
|
do k=1,ao_num
|
||||||
if (ao_value(k,j) /= 0.d0) then
|
c1 = ao_value(k,j)
|
||||||
c1 = ao_value(k,j)
|
if (c1 /= 0.d0) then
|
||||||
do i=1,mo_num
|
do i=1,mo_num
|
||||||
mo_value(i,j) = mo_value(i,j) + coefficient_t(i,k) * c1
|
mo_value(i,j) = mo_value(i,j) + coefficient_t(i,k) * c1
|
||||||
end do
|
end do
|
||||||
@ -736,7 +853,7 @@ integer function qmckl_compute_mo_basis_mo_value_doc_f(context, &
|
|||||||
end do
|
end do
|
||||||
end do
|
end do
|
||||||
|
|
||||||
else ! dgemm
|
else ! dgemm for checking
|
||||||
|
|
||||||
LDA = size(coefficient_t,1)
|
LDA = size(coefficient_t,1)
|
||||||
LDB = size(ao_value,1)
|
LDB = size(ao_value,1)
|
||||||
@ -874,7 +991,7 @@ qmckl_compute_mo_basis_mo_value_hpc (const qmckl_context context,
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
int64_t n;
|
int64_t n=0;
|
||||||
|
|
||||||
for (n=0 ; n < nidx-4 ; n+=4) {
|
for (n=0 ; n < nidx-4 ; n+=4) {
|
||||||
const double* restrict ck1 = coefficient_t + idx[n ]*mo_num;
|
const double* restrict ck1 = coefficient_t + idx[n ]*mo_num;
|
||||||
@ -895,8 +1012,7 @@ qmckl_compute_mo_basis_mo_value_hpc (const qmckl_context context,
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
const int64_t n0 = n < 0 ? 0 : n;
|
for (int64_t m=n ; m < nidx ; m+=1) {
|
||||||
for (int64_t m=n0 ; m < nidx ; m+=1) {
|
|
||||||
const double* restrict ck = coefficient_t + idx[m]*mo_num;
|
const double* restrict ck = coefficient_t + idx[m]*mo_num;
|
||||||
const double a1 = av1[m];
|
const double a1 = av1[m];
|
||||||
|
|
||||||
@ -1355,7 +1471,7 @@ qmckl_compute_mo_basis_mo_vgl_hpc (const qmckl_context context,
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
int64_t n;
|
int64_t n=0;
|
||||||
|
|
||||||
for (n=0 ; n < nidx-4 ; n+=4) {
|
for (n=0 ; n < nidx-4 ; n+=4) {
|
||||||
const double* restrict ck1 = coefficient_t + idx[n ]*mo_num;
|
const double* restrict ck1 = coefficient_t + idx[n ]*mo_num;
|
||||||
@ -1400,8 +1516,7 @@ qmckl_compute_mo_basis_mo_vgl_hpc (const qmckl_context context,
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
const int64_t n0 = n < 0 ? 0 : n;
|
for (int64_t m=n ; m < nidx ; m+=1) {
|
||||||
for (int64_t m=n0 ; m < nidx ; m+=1) {
|
|
||||||
const double* restrict ck = coefficient_t + idx[m]*mo_num;
|
const double* restrict ck = coefficient_t + idx[m]*mo_num;
|
||||||
const double a1 = av1[m];
|
const double a1 = av1[m];
|
||||||
const double a2 = av2[m];
|
const double a2 = av2[m];
|
||||||
@ -1603,7 +1718,7 @@ rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0]),
|
|||||||
assert (rc == QMCKL_SUCCESS);
|
assert (rc == QMCKL_SUCCESS);
|
||||||
|
|
||||||
/* Set up MO data */
|
/* Set up MO data */
|
||||||
const int64_t mo_num = chbrclf_mo_num;
|
int64_t mo_num = chbrclf_mo_num;
|
||||||
rc = qmckl_set_mo_basis_mo_num(context, mo_num);
|
rc = qmckl_set_mo_basis_mo_num(context, mo_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
assert (rc == QMCKL_SUCCESS);
|
||||||
|
|
||||||
@ -1699,6 +1814,32 @@ printf(" mo_vgl mo_vgl[1][26][223] %25.15e\n", mo_vgl[2][1][3]);
|
|||||||
printf(" mo_vgl mo_vgl[0][26][224] %25.15e\n", mo_vgl[2][0][3]);
|
printf(" mo_vgl mo_vgl[0][26][224] %25.15e\n", mo_vgl[2][0][3]);
|
||||||
printf(" mo_vgl mo_vgl[1][26][224] %25.15e\n", mo_vgl[2][1][3]);
|
printf(" mo_vgl mo_vgl[1][26][224] %25.15e\n", mo_vgl[2][1][3]);
|
||||||
printf("\n");
|
printf("\n");
|
||||||
|
|
||||||
|
|
||||||
|
/* Check selection of MOs */
|
||||||
|
|
||||||
|
int32_t keep[mo_num];
|
||||||
|
for (int i=0 ; i<mo_num ; ++i) {
|
||||||
|
keep[i] = 0;
|
||||||
|
}
|
||||||
|
keep[2] = 1;
|
||||||
|
keep[5] = 1;
|
||||||
|
|
||||||
|
rc = qmckl_mo_basis_select_mo(context, &(keep[0]), mo_num);
|
||||||
|
assert (rc == QMCKL_SUCCESS);
|
||||||
|
|
||||||
|
rc = qmckl_get_mo_basis_mo_num(context, &mo_num);
|
||||||
|
printf(" mo_num: %ld\n", mo_num);
|
||||||
|
assert(mo_num == 2);
|
||||||
|
|
||||||
|
double mo_coefficient_new[mo_num][ao_num];
|
||||||
|
rc = qmckl_get_mo_basis_coefficient (context, &(mo_coefficient_new[0][0]), mo_num*ao_num);
|
||||||
|
for (int i=0 ; i<ao_num ; ++i) {
|
||||||
|
assert(mo_coefficient_new[0][i] == mo_coefficient[i + ao_num*2]);
|
||||||
|
assert(mo_coefficient_new[1][i] == mo_coefficient[i + ao_num*5]);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
}
|
}
|
||||||
|
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|||||||
@ -74,15 +74,12 @@ int main() {
|
|||||||
| ~charge~ | qmckl_vector | Nuclear charges |
|
| ~charge~ | qmckl_vector | Nuclear charges |
|
||||||
| ~coord~ | qmckl_matrix | Nuclear coordinates, in transposed format |
|
| ~coord~ | qmckl_matrix | Nuclear coordinates, in transposed format |
|
||||||
| ~coord_date~ | int64_t | Nuclear coordinates, date if modified |
|
| ~coord_date~ | int64_t | Nuclear coordinates, date if modified |
|
||||||
| ~rescale_factor_kappa~ | double | The distance scaling factor |
|
|
||||||
|
|
||||||
Computed data:
|
Computed data:
|
||||||
|
|
||||||
|-----------------------------+--------------+------------------------------------------------------------|
|
|-----------------------------+--------------+------------------------------------------------------------|
|
||||||
| ~nn_distance~ | qmckl_matrix | Nucleus-nucleus distances |
|
| ~nn_distance~ | qmckl_matrix | Nucleus-nucleus distances |
|
||||||
| ~nn_distance_date~ | int64_t | Date when Nucleus-nucleus distances were computed |
|
| ~nn_distance_date~ | int64_t | Date when Nucleus-nucleus distances were computed |
|
||||||
| ~nn_distance_rescaled~ | qmckl_matrix | Nucleus-nucleus rescaled distances |
|
|
||||||
| ~nn_distance_rescaled_date~ | int64_t | Date when Nucleus-nucleus rescaled distances were computed |
|
|
||||||
| ~repulsion~ | double | Nuclear repulsion energy |
|
| ~repulsion~ | double | Nuclear repulsion energy |
|
||||||
| ~repulsion_date~ | int64_t | Date when the nuclear repulsion energy was computed |
|
| ~repulsion_date~ | int64_t | Date when the nuclear repulsion energy was computed |
|
||||||
|
|
||||||
@ -93,14 +90,11 @@ typedef struct qmckl_nucleus_struct {
|
|||||||
int64_t num;
|
int64_t num;
|
||||||
int64_t repulsion_date;
|
int64_t repulsion_date;
|
||||||
int64_t nn_distance_date;
|
int64_t nn_distance_date;
|
||||||
int64_t nn_distance_rescaled_date;
|
|
||||||
int64_t coord_date;
|
int64_t coord_date;
|
||||||
qmckl_vector charge;
|
qmckl_vector charge;
|
||||||
qmckl_matrix coord;
|
qmckl_matrix coord;
|
||||||
qmckl_matrix nn_distance;
|
qmckl_matrix nn_distance;
|
||||||
qmckl_matrix nn_distance_rescaled;
|
|
||||||
double repulsion;
|
double repulsion;
|
||||||
double rescale_factor_kappa;
|
|
||||||
int32_t uninitialized;
|
int32_t uninitialized;
|
||||||
bool provided;
|
bool provided;
|
||||||
} qmckl_nucleus_struct;
|
} qmckl_nucleus_struct;
|
||||||
@ -131,20 +125,12 @@ qmckl_exit_code qmckl_init_nucleus(qmckl_context context) {
|
|||||||
ctx->nucleus.uninitialized = (1 << 3) - 1;
|
ctx->nucleus.uninitialized = (1 << 3) - 1;
|
||||||
|
|
||||||
/* Default values */
|
/* Default values */
|
||||||
ctx->nucleus.rescale_factor_kappa = 1.0;
|
|
||||||
|
|
||||||
return QMCKL_SUCCESS;
|
return QMCKL_SUCCESS;
|
||||||
}
|
}
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
** Access functions
|
** Access functions
|
||||||
#+NAME:post
|
|
||||||
#+begin_src c :exports none
|
|
||||||
if ( (ctx->nucleus.uninitialized & mask) != 0) {
|
|
||||||
return NULL;
|
|
||||||
}
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
||||||
qmckl_exit_code
|
qmckl_exit_code
|
||||||
@ -272,53 +258,6 @@ end interface
|
|||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
||||||
qmckl_exit_code
|
qmckl_exit_code
|
||||||
qmckl_get_nucleus_rescale_factor(const qmckl_context context,
|
|
||||||
double* const rescale_factor_kappa);
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_get_nucleus_rescale_factor (const qmckl_context context,
|
|
||||||
double* const rescale_factor_kappa)
|
|
||||||
{
|
|
||||||
|
|
||||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
|
||||||
return QMCKL_INVALID_CONTEXT;
|
|
||||||
}
|
|
||||||
|
|
||||||
if (rescale_factor_kappa == NULL) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_INVALID_ARG_2,
|
|
||||||
"qmckl_get_nucleus_rescale_factor",
|
|
||||||
"rescale_factor_kappa is a null pointer");
|
|
||||||
}
|
|
||||||
|
|
||||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
|
||||||
assert (ctx != NULL);
|
|
||||||
|
|
||||||
assert (ctx->nucleus.rescale_factor_kappa > 0.0);
|
|
||||||
|
|
||||||
(*rescale_factor_kappa) = ctx->nucleus.rescale_factor_kappa;
|
|
||||||
|
|
||||||
return QMCKL_SUCCESS;
|
|
||||||
}
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
|
||||||
interface
|
|
||||||
integer(c_int32_t) function qmckl_get_nucleus_rescale_factor(context, kappa) &
|
|
||||||
bind(C)
|
|
||||||
use, intrinsic :: iso_c_binding
|
|
||||||
import
|
|
||||||
implicit none
|
|
||||||
integer (c_int64_t) , intent(in) , value :: context
|
|
||||||
real (c_double) , intent(out) :: kappa
|
|
||||||
end function qmckl_get_nucleus_rescale_factor
|
|
||||||
end interface
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_func) :exports none
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_get_nucleus_coord(const qmckl_context context,
|
qmckl_get_nucleus_coord(const qmckl_context context,
|
||||||
const char transp,
|
const char transp,
|
||||||
double* const coord,
|
double* const coord,
|
||||||
@ -438,7 +377,7 @@ if (mask != 0 && !(ctx->nucleus.uninitialized & mask)) {
|
|||||||
}
|
}
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
#+NAME:post2
|
#+NAME:post
|
||||||
#+begin_src c :exports none
|
#+begin_src c :exports none
|
||||||
ctx->nucleus.uninitialized &= ~mask;
|
ctx->nucleus.uninitialized &= ~mask;
|
||||||
ctx->nucleus.provided = (ctx->nucleus.uninitialized == 0);
|
ctx->nucleus.provided = (ctx->nucleus.uninitialized == 0);
|
||||||
@ -475,7 +414,7 @@ qmckl_set_nucleus_num(qmckl_context context,
|
|||||||
|
|
||||||
ctx->nucleus.num = num;
|
ctx->nucleus.num = num;
|
||||||
|
|
||||||
<<post2>>
|
<<post>>
|
||||||
}
|
}
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
@ -541,7 +480,7 @@ qmckl_set_nucleus_charge(qmckl_context context,
|
|||||||
"Error in vector->double* conversion");
|
"Error in vector->double* conversion");
|
||||||
}
|
}
|
||||||
|
|
||||||
<<post2>>
|
<<post>>
|
||||||
}
|
}
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
@ -619,7 +558,7 @@ qmckl_set_nucleus_coord(qmckl_context context,
|
|||||||
}
|
}
|
||||||
if (rc != QMCKL_SUCCESS) return rc;
|
if (rc != QMCKL_SUCCESS) return rc;
|
||||||
|
|
||||||
<<post2>>
|
<<post>>
|
||||||
}
|
}
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
@ -638,55 +577,12 @@ interface
|
|||||||
end interface
|
end interface
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_func)
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_set_nucleus_rescale_factor(qmckl_context context,
|
|
||||||
const double kappa);
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
Sets the rescale parameter for the nuclear distances.
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_set_nucleus_rescale_factor(qmckl_context context,
|
|
||||||
const double rescale_factor_kappa)
|
|
||||||
{
|
|
||||||
int32_t mask = 0; // Can be updated
|
|
||||||
|
|
||||||
<<pre2>>
|
|
||||||
|
|
||||||
if (rescale_factor_kappa <= 0.0) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_INVALID_ARG_2,
|
|
||||||
"qmckl_set_nucleus_rescale_factor",
|
|
||||||
"rescale_factor_kappa cannot be <= 0.");
|
|
||||||
}
|
|
||||||
|
|
||||||
ctx->nucleus.rescale_factor_kappa = rescale_factor_kappa;
|
|
||||||
|
|
||||||
return QMCKL_SUCCESS;
|
|
||||||
}
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
|
||||||
interface
|
|
||||||
integer(c_int32_t) function qmckl_set_nucleus_rescale_factor(context, kappa) &
|
|
||||||
bind(C)
|
|
||||||
use, intrinsic :: iso_c_binding
|
|
||||||
import
|
|
||||||
implicit none
|
|
||||||
integer (c_int64_t) , intent(in) , value :: context
|
|
||||||
real (c_double) , intent(in) , value :: kappa
|
|
||||||
end function qmckl_set_nucleus_rescale_factor
|
|
||||||
end interface
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
** Test
|
** Test
|
||||||
|
|
||||||
#+begin_src c :tangle (eval c_test)
|
#+begin_src c :tangle (eval c_test)
|
||||||
const double* nucl_charge = chbrclf_charge;
|
const double* nucl_charge = chbrclf_charge;
|
||||||
const double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
|
const double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
|
||||||
const double nucl_rescale_factor_kappa = 2.0;
|
|
||||||
|
|
||||||
/* --- */
|
/* --- */
|
||||||
|
|
||||||
@ -700,38 +596,25 @@ assert(rc == QMCKL_NOT_PROVIDED);
|
|||||||
|
|
||||||
|
|
||||||
rc = qmckl_set_nucleus_num (context, chbrclf_nucl_num);
|
rc = qmckl_set_nucleus_num (context, chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(!qmckl_nucleus_provided(context));
|
assert(!qmckl_nucleus_provided(context));
|
||||||
|
|
||||||
rc = qmckl_get_nucleus_num (context, &n);
|
rc = qmckl_get_nucleus_num (context, &n);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert(n == chbrclf_nucl_num);
|
assert(n == chbrclf_nucl_num);
|
||||||
|
|
||||||
double k;
|
|
||||||
rc = qmckl_get_nucleus_rescale_factor (context, &k);
|
|
||||||
assert(rc == QMCKL_SUCCESS);
|
|
||||||
assert(k == 1.0);
|
|
||||||
|
|
||||||
|
|
||||||
rc = qmckl_set_nucleus_rescale_factor (context, nucl_rescale_factor_kappa);
|
|
||||||
assert(rc == QMCKL_SUCCESS);
|
|
||||||
|
|
||||||
rc = qmckl_get_nucleus_rescale_factor (context, &k);
|
|
||||||
assert(rc == QMCKL_SUCCESS);
|
|
||||||
assert(k == nucl_rescale_factor_kappa);
|
|
||||||
|
|
||||||
double nucl_coord2[3*chbrclf_nucl_num];
|
double nucl_coord2[3*chbrclf_nucl_num];
|
||||||
|
|
||||||
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2, 3*chbrclf_nucl_num);
|
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2, 3*chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_NOT_PROVIDED);
|
assert(rc == QMCKL_NOT_PROVIDED);
|
||||||
|
|
||||||
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), 3*chbrclf_nucl_num);
|
rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0]), 3*chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
assert(!qmckl_nucleus_provided(context));
|
assert(!qmckl_nucleus_provided(context));
|
||||||
|
|
||||||
rc = qmckl_get_nucleus_coord (context, 'N', nucl_coord2, 3*chbrclf_nucl_num);
|
rc = qmckl_get_nucleus_coord (context, 'N', nucl_coord2, 3*chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (size_t k=0 ; k<3 ; ++k) {
|
for (size_t k=0 ; k<3 ; ++k) {
|
||||||
for (int64_t i=0 ; i<chbrclf_nucl_num ; ++i) {
|
for (int64_t i=0 ; i<chbrclf_nucl_num ; ++i) {
|
||||||
assert( nucl_coord[chbrclf_nucl_num*k+i] == nucl_coord2[3*i+k] );
|
assert( nucl_coord[chbrclf_nucl_num*k+i] == nucl_coord2[3*i+k] );
|
||||||
@ -739,7 +622,7 @@ for (size_t k=0 ; k<3 ; ++k) {
|
|||||||
}
|
}
|
||||||
|
|
||||||
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2, 3*chbrclf_nucl_num);
|
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2, 3*chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int64_t i=0 ; i<3*chbrclf_nucl_num ; ++i) {
|
for (int64_t i=0 ; i<3*chbrclf_nucl_num ; ++i) {
|
||||||
assert( nucl_coord[i] == nucl_coord2[i] );
|
assert( nucl_coord[i] == nucl_coord2[i] );
|
||||||
}
|
}
|
||||||
@ -750,10 +633,10 @@ rc = qmckl_get_nucleus_charge(context, nucl_charge2, chbrclf_nucl_num);
|
|||||||
assert(rc == QMCKL_NOT_PROVIDED);
|
assert(rc == QMCKL_NOT_PROVIDED);
|
||||||
|
|
||||||
rc = qmckl_set_nucleus_charge(context, nucl_charge, chbrclf_nucl_num);
|
rc = qmckl_set_nucleus_charge(context, nucl_charge, chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
rc = qmckl_get_nucleus_charge(context, nucl_charge2, chbrclf_nucl_num);
|
rc = qmckl_get_nucleus_charge(context, nucl_charge2, chbrclf_nucl_num);
|
||||||
assert(rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int64_t i=0 ; i<chbrclf_nucl_num ; ++i) {
|
for (int64_t i=0 ; i<chbrclf_nucl_num ; ++i) {
|
||||||
assert( nucl_charge[i] == nucl_charge2[i] );
|
assert( nucl_charge[i] == nucl_charge2[i] );
|
||||||
}
|
}
|
||||||
@ -959,202 +842,6 @@ assert(fabs(distance[1]-2.070304721365169) < 1.e-12);
|
|||||||
|
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
** Nucleus-nucleus rescaled distances
|
|
||||||
|
|
||||||
*** Get
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_func) :noweb yes
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_get_nucleus_nn_distance_rescaled(qmckl_context context,
|
|
||||||
double* distance_rescaled,
|
|
||||||
const int64_t size_max);
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
|
||||||
qmckl_exit_code
|
|
||||||
qmckl_get_nucleus_nn_distance_rescaled(qmckl_context context,
|
|
||||||
double* distance_rescaled,
|
|
||||||
const int64_t size_max)
|
|
||||||
{
|
|
||||||
/* Check input parameters */
|
|
||||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
|
||||||
return (char) 0;
|
|
||||||
}
|
|
||||||
|
|
||||||
qmckl_exit_code rc = qmckl_provide_nn_distance_rescaled(context);
|
|
||||||
if (rc != QMCKL_SUCCESS) return rc;
|
|
||||||
|
|
||||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
|
||||||
assert (ctx != NULL);
|
|
||||||
|
|
||||||
const int64_t sze = ctx->nucleus.num * ctx->nucleus.num;
|
|
||||||
if (sze > size_max) {
|
|
||||||
return qmckl_failwith(context,
|
|
||||||
QMCKL_INVALID_ARG_3,
|
|
||||||
"qmckl_get_nucleus_nn_distance_rescaled",
|
|
||||||
"Array too small");
|
|
||||||
}
|
|
||||||
rc = qmckl_double_of_matrix(context,
|
|
||||||
ctx->nucleus.nn_distance_rescaled,
|
|
||||||
distance_rescaled,
|
|
||||||
size_max);
|
|
||||||
return rc;
|
|
||||||
}
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
|
|
||||||
interface
|
|
||||||
integer(c_int32_t) function qmckl_get_nucleus_nn_distance_rescaled(context, distance_rescaled, size_max) &
|
|
||||||
bind(C)
|
|
||||||
use, intrinsic :: iso_c_binding
|
|
||||||
import
|
|
||||||
implicit none
|
|
||||||
integer (c_int64_t) , intent(in) , value :: context
|
|
||||||
real (c_double ) , intent(out) :: distance_rescaled(*)
|
|
||||||
integer (c_int64_t) , intent(in) , value :: size_max
|
|
||||||
end function
|
|
||||||
end interface
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
*** Provide :noexport:
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval h_private_func) :noweb yes :exports none
|
|
||||||
qmckl_exit_code qmckl_provide_nn_distance_rescaled(qmckl_context context);
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
|
|
||||||
qmckl_exit_code qmckl_provide_nn_distance_rescaled(qmckl_context context)
|
|
||||||
{
|
|
||||||
/* Check input parameters */
|
|
||||||
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
|
|
||||||
return (char) 0;
|
|
||||||
}
|
|
||||||
|
|
||||||
qmckl_context_struct* const ctx = (qmckl_context_struct*) context;
|
|
||||||
assert (ctx != NULL);
|
|
||||||
|
|
||||||
if (!ctx->nucleus.provided) return QMCKL_NOT_PROVIDED;
|
|
||||||
|
|
||||||
/* Allocate array */
|
|
||||||
if (ctx->nucleus.nn_distance_rescaled.data == NULL) {
|
|
||||||
ctx->nucleus.nn_distance_rescaled =
|
|
||||||
qmckl_matrix_alloc(context, ctx->nucleus.num, ctx->nucleus.num);
|
|
||||||
|
|
||||||
if (ctx->nucleus.nn_distance_rescaled.data == NULL) {
|
|
||||||
return qmckl_failwith( context,
|
|
||||||
QMCKL_ALLOCATION_FAILED,
|
|
||||||
"qmckl_nn_distance_rescaled",
|
|
||||||
NULL);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
qmckl_exit_code rc =
|
|
||||||
qmckl_compute_nn_distance_rescaled(context,
|
|
||||||
ctx->nucleus.num,
|
|
||||||
ctx->nucleus.rescale_factor_kappa,
|
|
||||||
ctx->nucleus.coord.data,
|
|
||||||
ctx->nucleus.nn_distance_rescaled.data);
|
|
||||||
if (rc != QMCKL_SUCCESS) {
|
|
||||||
return rc;
|
|
||||||
}
|
|
||||||
|
|
||||||
ctx->nucleus.nn_distance_rescaled_date = ctx->date;
|
|
||||||
|
|
||||||
return QMCKL_SUCCESS;
|
|
||||||
}
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
*** Compute
|
|
||||||
|
|
||||||
#+NAME: qmckl_nn_distance_rescaled_args
|
|
||||||
| qmckl_context | context | in | Global state |
|
|
||||||
| int64_t | nucl_num | in | Number of nuclei |
|
|
||||||
| double | coord[3][nucl_num] | in | Nuclear coordinates (au) |
|
|
||||||
| double | nn_distance_rescaled[nucl_num][nucl_num] | out | Nucleus-nucleus rescaled distances (au) |
|
|
||||||
|
|
||||||
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
|
|
||||||
integer function qmckl_compute_nn_distance_rescaled_f(context, nucl_num, rescale_factor_kappa, coord, nn_distance_rescaled) &
|
|
||||||
result(info)
|
|
||||||
use qmckl
|
|
||||||
implicit none
|
|
||||||
integer(qmckl_context), intent(in) :: context
|
|
||||||
integer*8 , intent(in) :: nucl_num
|
|
||||||
double precision , intent(in) :: rescale_factor_kappa
|
|
||||||
double precision , intent(in) :: coord(nucl_num,3)
|
|
||||||
double precision , intent(out) :: nn_distance_rescaled(nucl_num,nucl_num)
|
|
||||||
|
|
||||||
integer*8 :: k
|
|
||||||
|
|
||||||
info = QMCKL_SUCCESS
|
|
||||||
|
|
||||||
if (context == QMCKL_NULL_CONTEXT) then
|
|
||||||
info = QMCKL_INVALID_CONTEXT
|
|
||||||
return
|
|
||||||
endif
|
|
||||||
|
|
||||||
if (nucl_num <= 0) then
|
|
||||||
info = QMCKL_INVALID_ARG_2
|
|
||||||
return
|
|
||||||
endif
|
|
||||||
|
|
||||||
info = qmckl_distance_rescaled(context, 'T', 'T', nucl_num, nucl_num, &
|
|
||||||
coord, nucl_num, &
|
|
||||||
coord, nucl_num, &
|
|
||||||
nn_distance_rescaled, nucl_num, rescale_factor_kappa)
|
|
||||||
|
|
||||||
end function qmckl_compute_nn_distance_rescaled_f
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
#+begin_src c :tangle (eval h_private_func) :comments org :exports none
|
|
||||||
qmckl_exit_code qmckl_compute_nn_distance_rescaled (
|
|
||||||
const qmckl_context context,
|
|
||||||
const int64_t nucl_num,
|
|
||||||
const double rescale_factor_kappa,
|
|
||||||
const double* coord,
|
|
||||||
double* const nn_distance_rescaled );
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
|
|
||||||
#+CALL: generate_c_interface(table=qmckl_nn_distance_rescaled_args,rettyp="qmckl_exit_code",fname="qmckl_compute_nn_distance")
|
|
||||||
|
|
||||||
#+RESULTS:
|
|
||||||
#+begin_src f90 :tangle (eval f) :comments org :exports none
|
|
||||||
integer(c_int32_t) function qmckl_compute_nn_distance_rescaled &
|
|
||||||
(context, nucl_num, rescale_factor_kappa, coord, nn_distance_rescaled) &
|
|
||||||
bind(C) result(info)
|
|
||||||
|
|
||||||
use, intrinsic :: iso_c_binding
|
|
||||||
implicit none
|
|
||||||
|
|
||||||
integer (c_int64_t) , intent(in) , value :: context
|
|
||||||
integer (c_int64_t) , intent(in) , value :: nucl_num
|
|
||||||
real (c_double ) , intent(in) , value :: rescale_factor_kappa
|
|
||||||
real (c_double ) , intent(in) :: coord(nucl_num,3)
|
|
||||||
real (c_double ) , intent(out) :: nn_distance_rescaled(nucl_num,nucl_num)
|
|
||||||
|
|
||||||
integer(c_int32_t), external :: qmckl_compute_nn_distance_rescaled_f
|
|
||||||
info = qmckl_compute_nn_distance_rescaled_f &
|
|
||||||
(context, nucl_num, rescale_factor_kappa, coord, nn_distance_rescaled)
|
|
||||||
|
|
||||||
end function qmckl_compute_nn_distance_rescaled
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
*** Test
|
|
||||||
|
|
||||||
#+begin_src c :tangle (eval c_test)
|
|
||||||
/* Reference input data */
|
|
||||||
/* TODO */
|
|
||||||
|
|
||||||
//assert(qmckl_nucleus_provided(context));
|
|
||||||
//
|
|
||||||
//double distance[nucl_num*nucl_num];
|
|
||||||
//rc = qmckl_get_nucleus_nn_distance(context, distance, nucl_num*nucl_num);
|
|
||||||
//assert(distance[0] == 0.);
|
|
||||||
//assert(distance[1] == distance[nucl_num]);
|
|
||||||
//assert(fabs(distance[1]-2.070304721365169) < 1.e-12);
|
|
||||||
|
|
||||||
#+end_src
|
|
||||||
|
|
||||||
** Nuclear repulsion energy
|
** Nuclear repulsion energy
|
||||||
|
|
||||||
\[
|
\[
|
||||||
|
|||||||
@ -1097,37 +1097,25 @@ qmckl_trexio_read(const qmckl_context context, const char* file_name, const int6
|
|||||||
rc = qmckl_trexio_read_electron_X(context, file);
|
rc = qmckl_trexio_read_electron_X(context, file);
|
||||||
if (rc != QMCKL_SUCCESS) {
|
if (rc != QMCKL_SUCCESS) {
|
||||||
trexio_close(file);
|
trexio_close(file);
|
||||||
return qmckl_failwith( context,
|
return rc;
|
||||||
QMCKL_FAILURE,
|
|
||||||
"qmckl_trexio_read",
|
|
||||||
"Error reading electron");
|
|
||||||
}
|
}
|
||||||
|
|
||||||
rc = qmckl_trexio_read_nucleus_X(context, file);
|
rc = qmckl_trexio_read_nucleus_X(context, file);
|
||||||
if (rc != QMCKL_SUCCESS) {
|
if (rc != QMCKL_SUCCESS) {
|
||||||
trexio_close(file);
|
trexio_close(file);
|
||||||
return qmckl_failwith( context,
|
return rc;
|
||||||
QMCKL_FAILURE,
|
|
||||||
"qmckl_trexio_read",
|
|
||||||
"Error reading nucleus");
|
|
||||||
}
|
}
|
||||||
|
|
||||||
rc = qmckl_trexio_read_ao_X(context, file);
|
rc = qmckl_trexio_read_ao_X(context, file);
|
||||||
if (rc != QMCKL_SUCCESS) {
|
if (rc != QMCKL_SUCCESS) {
|
||||||
trexio_close(file);
|
trexio_close(file);
|
||||||
return qmckl_failwith( context,
|
return rc;
|
||||||
QMCKL_FAILURE,
|
|
||||||
"qmckl_trexio_read",
|
|
||||||
"Error reading AOs");
|
|
||||||
}
|
}
|
||||||
|
|
||||||
rc = qmckl_trexio_read_mo_X(context, file);
|
rc = qmckl_trexio_read_mo_X(context, file);
|
||||||
if (rc != QMCKL_SUCCESS) {
|
if (rc != QMCKL_SUCCESS) {
|
||||||
trexio_close(file);
|
trexio_close(file);
|
||||||
return qmckl_failwith( context,
|
return rc;
|
||||||
QMCKL_FAILURE,
|
|
||||||
"qmckl_trexio_read",
|
|
||||||
"Error reading MOs");
|
|
||||||
}
|
}
|
||||||
|
|
||||||
trexio_close(file);
|
trexio_close(file);
|
||||||
@ -1149,27 +1137,19 @@ qmckl_trexio_read(const qmckl_context context, const char* file_name, const int6
|
|||||||
#ifdef HAVE_TREXIO
|
#ifdef HAVE_TREXIO
|
||||||
|
|
||||||
qmckl_exit_code rc;
|
qmckl_exit_code rc;
|
||||||
char fname[256];
|
char filename[256];
|
||||||
char message[256];
|
|
||||||
|
|
||||||
#ifndef QMCKL_TEST_DIR
|
#ifndef QMCKL_TEST_DIR
|
||||||
#error "QMCKL_TEST_DIR is not defined"
|
#error "QMCKL_TEST_DIR is not defined"
|
||||||
#endif
|
#endif
|
||||||
|
|
||||||
strncpy(fname, QMCKL_TEST_DIR,255);
|
strncpy(filename, QMCKL_TEST_DIR,255);
|
||||||
strncat(fname, "/chbrclf", 255);
|
strncat(filename, "/chbrclf", 255);
|
||||||
printf("Test file: %s\n", fname);
|
printf("Test file: %s\n", filename);
|
||||||
|
|
||||||
rc = qmckl_trexio_read(context, fname, 255);
|
rc = qmckl_trexio_read(context, filename, 255);
|
||||||
|
|
||||||
if (rc != QMCKL_SUCCESS) {
|
qmckl_check(context, rc);
|
||||||
printf("%s\n", qmckl_string_of_error(rc));
|
|
||||||
qmckl_get_error(context, &rc, fname, message);
|
|
||||||
printf("%s\n", fname);
|
|
||||||
printf("%s\n", message);
|
|
||||||
}
|
|
||||||
|
|
||||||
assert ( rc == QMCKL_SUCCESS );
|
|
||||||
|
|
||||||
#+end_src
|
#+end_src
|
||||||
|
|
||||||
@ -1179,11 +1159,11 @@ assert ( rc == QMCKL_SUCCESS );
|
|||||||
printf("Electrons\n");
|
printf("Electrons\n");
|
||||||
int64_t up_num, dn_num;
|
int64_t up_num, dn_num;
|
||||||
rc = qmckl_get_electron_up_num(context, &up_num);
|
rc = qmckl_get_electron_up_num(context, &up_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert (up_num == chbrclf_elec_up_num);
|
assert (up_num == chbrclf_elec_up_num);
|
||||||
|
|
||||||
rc = qmckl_get_electron_down_num(context, &dn_num);
|
rc = qmckl_get_electron_down_num(context, &dn_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert (dn_num == chbrclf_elec_dn_num);
|
assert (dn_num == chbrclf_elec_dn_num);
|
||||||
|
|
||||||
#+end_src
|
#+end_src
|
||||||
@ -1195,13 +1175,13 @@ printf("Nuclei\n");
|
|||||||
|
|
||||||
int64_t nucl_num;
|
int64_t nucl_num;
|
||||||
rc = qmckl_get_nucleus_num(context, &nucl_num);
|
rc = qmckl_get_nucleus_num(context, &nucl_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
assert (nucl_num == chbrclf_nucl_num);
|
assert (nucl_num == chbrclf_nucl_num);
|
||||||
|
|
||||||
printf("Nuclear charges\n");
|
printf("Nuclear charges\n");
|
||||||
double * charge = (double*) malloc (nucl_num * sizeof(double));
|
double * charge = (double*) malloc (nucl_num * sizeof(double));
|
||||||
rc = qmckl_get_nucleus_charge(context, charge, nucl_num);
|
rc = qmckl_get_nucleus_charge(context, charge, nucl_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<nucl_num ; i++) {
|
for (int i=0 ; i<nucl_num ; i++) {
|
||||||
assert (charge[i] == chbrclf_charge[i]);
|
assert (charge[i] == chbrclf_charge[i]);
|
||||||
}
|
}
|
||||||
@ -1211,7 +1191,7 @@ charge = NULL;
|
|||||||
printf("Nuclear coordinates\n");
|
printf("Nuclear coordinates\n");
|
||||||
double * coord = (double*) malloc (nucl_num * 3 * sizeof(double));
|
double * coord = (double*) malloc (nucl_num * 3 * sizeof(double));
|
||||||
rc = qmckl_get_nucleus_coord(context, 'T', coord, 3*nucl_num);
|
rc = qmckl_get_nucleus_coord(context, 'T', coord, 3*nucl_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
int k=0;
|
int k=0;
|
||||||
for (int j=0 ; j<3 ; ++j) {
|
for (int j=0 ; j<3 ; ++j) {
|
||||||
for (int i=0 ; i<nucl_num ; ++i) {
|
for (int i=0 ; i<nucl_num ; ++i) {
|
||||||
@ -1248,7 +1228,7 @@ assert (ao_num == chbrclf_ao_num);
|
|||||||
|
|
||||||
int64_t* nucleus_index = (int64_t*) malloc (nucl_num * sizeof(int64_t));
|
int64_t* nucleus_index = (int64_t*) malloc (nucl_num * sizeof(int64_t));
|
||||||
rc = qmckl_get_ao_basis_nucleus_index(context, nucleus_index, nucl_num);
|
rc = qmckl_get_ao_basis_nucleus_index(context, nucleus_index, nucl_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<nucl_num ; i++) {
|
for (int i=0 ; i<nucl_num ; i++) {
|
||||||
assert (nucleus_index[i] == chbrclf_basis_nucleus_index[i]);
|
assert (nucleus_index[i] == chbrclf_basis_nucleus_index[i]);
|
||||||
}
|
}
|
||||||
@ -1257,7 +1237,7 @@ nucleus_index = NULL;
|
|||||||
|
|
||||||
int64_t* nucleus_shell_num = (int64_t*) malloc (nucl_num * sizeof(int64_t));
|
int64_t* nucleus_shell_num = (int64_t*) malloc (nucl_num * sizeof(int64_t));
|
||||||
rc = qmckl_get_ao_basis_nucleus_shell_num(context, nucleus_shell_num, nucl_num);
|
rc = qmckl_get_ao_basis_nucleus_shell_num(context, nucleus_shell_num, nucl_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<nucl_num ; i++) {
|
for (int i=0 ; i<nucl_num ; i++) {
|
||||||
assert (nucleus_shell_num[i] == chbrclf_basis_nucleus_shell_num[i]);
|
assert (nucleus_shell_num[i] == chbrclf_basis_nucleus_shell_num[i]);
|
||||||
}
|
}
|
||||||
@ -1266,7 +1246,7 @@ nucleus_shell_num = NULL;
|
|||||||
|
|
||||||
int32_t* shell_ang_mom = (int32_t*) malloc (shell_num * sizeof(int32_t));
|
int32_t* shell_ang_mom = (int32_t*) malloc (shell_num * sizeof(int32_t));
|
||||||
rc = qmckl_get_ao_basis_shell_ang_mom(context, shell_ang_mom, shell_num);
|
rc = qmckl_get_ao_basis_shell_ang_mom(context, shell_ang_mom, shell_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<shell_num ; i++) {
|
for (int i=0 ; i<shell_num ; i++) {
|
||||||
assert (shell_ang_mom[i] == chbrclf_basis_shell_ang_mom[i]);
|
assert (shell_ang_mom[i] == chbrclf_basis_shell_ang_mom[i]);
|
||||||
}
|
}
|
||||||
@ -1275,7 +1255,7 @@ shell_ang_mom = NULL;
|
|||||||
|
|
||||||
int64_t* shell_prim_num = (int64_t*) malloc (shell_num * sizeof(int64_t));
|
int64_t* shell_prim_num = (int64_t*) malloc (shell_num * sizeof(int64_t));
|
||||||
rc = qmckl_get_ao_basis_shell_prim_num(context, shell_prim_num, shell_num);
|
rc = qmckl_get_ao_basis_shell_prim_num(context, shell_prim_num, shell_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<shell_num ; i++) {
|
for (int i=0 ; i<shell_num ; i++) {
|
||||||
assert (shell_prim_num[i] == chbrclf_basis_shell_prim_num[i]);
|
assert (shell_prim_num[i] == chbrclf_basis_shell_prim_num[i]);
|
||||||
}
|
}
|
||||||
@ -1284,7 +1264,7 @@ shell_prim_num = NULL;
|
|||||||
|
|
||||||
int64_t* shell_prim_index = (int64_t*) malloc (shell_num * sizeof(int64_t));
|
int64_t* shell_prim_index = (int64_t*) malloc (shell_num * sizeof(int64_t));
|
||||||
rc = qmckl_get_ao_basis_shell_prim_index(context, shell_prim_index, shell_num);
|
rc = qmckl_get_ao_basis_shell_prim_index(context, shell_prim_index, shell_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<shell_num ; i++) {
|
for (int i=0 ; i<shell_num ; i++) {
|
||||||
assert (shell_prim_index[i] == chbrclf_basis_shell_prim_index[i]);
|
assert (shell_prim_index[i] == chbrclf_basis_shell_prim_index[i]);
|
||||||
}
|
}
|
||||||
@ -1293,7 +1273,7 @@ shell_prim_index = NULL;
|
|||||||
|
|
||||||
double* shell_factor = (double*) malloc (shell_num * sizeof(double));
|
double* shell_factor = (double*) malloc (shell_num * sizeof(double));
|
||||||
rc = qmckl_get_ao_basis_shell_factor(context, shell_factor, shell_num);
|
rc = qmckl_get_ao_basis_shell_factor(context, shell_factor, shell_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<shell_num ; i++) {
|
for (int i=0 ; i<shell_num ; i++) {
|
||||||
assert (fabs(shell_factor[i] - chbrclf_basis_shell_factor[i]) < 1.e-6);
|
assert (fabs(shell_factor[i] - chbrclf_basis_shell_factor[i]) < 1.e-6);
|
||||||
}
|
}
|
||||||
@ -1302,7 +1282,7 @@ shell_factor = NULL;
|
|||||||
|
|
||||||
double* exponent = (double*) malloc (prim_num * sizeof(double));
|
double* exponent = (double*) malloc (prim_num * sizeof(double));
|
||||||
rc = qmckl_get_ao_basis_exponent(context, exponent, prim_num);
|
rc = qmckl_get_ao_basis_exponent(context, exponent, prim_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<prim_num ; i++) {
|
for (int i=0 ; i<prim_num ; i++) {
|
||||||
assert (fabs((exponent[i] - chbrclf_basis_exponent[i])/chbrclf_basis_exponent[i]) < 1.e-7);
|
assert (fabs((exponent[i] - chbrclf_basis_exponent[i])/chbrclf_basis_exponent[i]) < 1.e-7);
|
||||||
}
|
}
|
||||||
@ -1311,7 +1291,7 @@ exponent = NULL;
|
|||||||
|
|
||||||
double* coefficient = (double*) malloc (prim_num * sizeof(double));
|
double* coefficient = (double*) malloc (prim_num * sizeof(double));
|
||||||
rc = qmckl_get_ao_basis_coefficient(context, coefficient, prim_num);
|
rc = qmckl_get_ao_basis_coefficient(context, coefficient, prim_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<prim_num ; i++) {
|
for (int i=0 ; i<prim_num ; i++) {
|
||||||
assert (fabs((coefficient[i] - chbrclf_basis_coefficient[i])/chbrclf_basis_coefficient[i]) < 1.e-7);
|
assert (fabs((coefficient[i] - chbrclf_basis_coefficient[i])/chbrclf_basis_coefficient[i]) < 1.e-7);
|
||||||
}
|
}
|
||||||
@ -1320,7 +1300,7 @@ coefficient = NULL;
|
|||||||
|
|
||||||
double* prim_factor = (double*) malloc (prim_num * sizeof(double));
|
double* prim_factor = (double*) malloc (prim_num * sizeof(double));
|
||||||
rc = qmckl_get_ao_basis_prim_factor(context, prim_factor, prim_num);
|
rc = qmckl_get_ao_basis_prim_factor(context, prim_factor, prim_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<prim_num ; i++) {
|
for (int i=0 ; i<prim_num ; i++) {
|
||||||
assert (fabs((prim_factor[i] - chbrclf_basis_prim_factor[i])/chbrclf_basis_prim_factor[i]) < 1.e-7);
|
assert (fabs((prim_factor[i] - chbrclf_basis_prim_factor[i])/chbrclf_basis_prim_factor[i]) < 1.e-7);
|
||||||
}
|
}
|
||||||
@ -1337,13 +1317,14 @@ printf("MOs\n");
|
|||||||
|
|
||||||
int64_t mo_num;
|
int64_t mo_num;
|
||||||
rc = qmckl_get_mo_basis_mo_num(context, &mo_num);
|
rc = qmckl_get_mo_basis_mo_num(context, &mo_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
|
|
||||||
assert (mo_num == chbrclf_mo_num);
|
assert (mo_num == chbrclf_mo_num);
|
||||||
|
|
||||||
printf("MO coefs\n");
|
printf("MO coefs\n");
|
||||||
double * mo_coef = (double*) malloc (ao_num * mo_num * sizeof(double));
|
double * mo_coef = (double*) malloc (ao_num * mo_num * sizeof(double));
|
||||||
rc = qmckl_get_mo_basis_coefficient(context, mo_coef, mo_num*ao_num);
|
rc = qmckl_get_mo_basis_coefficient(context, mo_coef, mo_num*ao_num);
|
||||||
assert (rc == QMCKL_SUCCESS);
|
qmckl_check(context, rc);
|
||||||
for (int i=0 ; i<ao_num * mo_num ; i++) {
|
for (int i=0 ; i<ao_num * mo_num ; i++) {
|
||||||
printf("%d %e %e %e\n", i, mo_coef[i], chbrclf_mo_coef[i],
|
printf("%d %e %e %e\n", i, mo_coef[i], chbrclf_mo_coef[i],
|
||||||
( fabs(mo_coef[i] - chbrclf_mo_coef[i])/fabs(mo_coef[i])) );
|
( fabs(mo_coef[i] - chbrclf_mo_coef[i])/fabs(mo_coef[i])) );
|
||||||
|
|||||||
@ -2,6 +2,6 @@
|
|||||||
requires = [
|
requires = [
|
||||||
"setuptools>=42",
|
"setuptools>=42",
|
||||||
"wheel",
|
"wheel",
|
||||||
"numpy>=1.17.3"
|
"oldest-supported-numpy"
|
||||||
]
|
]
|
||||||
build-backend = "setuptools.build_meta"
|
build-backend = "setuptools.build_meta"
|
||||||
|
|||||||
@ -71,5 +71,7 @@ setup(name = MODULE_NAME,
|
|||||||
"Operating System :: MacOS"
|
"Operating System :: MacOS"
|
||||||
],
|
],
|
||||||
python_requires = ">=3.0",
|
python_requires = ">=3.0",
|
||||||
install_requires = ["numpy>=1.17.3"]
|
# The ABI incompatibility of numpy is a pain, for now set the
|
||||||
|
# min numpy version such that we have wheels for CPython 3.5 & 3.6
|
||||||
|
install_requires = ["numpy>=1.13.3"]
|
||||||
)
|
)
|
||||||
|
|||||||
@ -20,6 +20,7 @@
|
|||||||
%include typemaps.i
|
%include typemaps.i
|
||||||
|
|
||||||
%apply int *OUTPUT { qmckl_exit_code *exit_code};
|
%apply int *OUTPUT { qmckl_exit_code *exit_code};
|
||||||
|
%apply double *OUTPUT { double* det_l};
|
||||||
|
|
||||||
/* Avoid passing file_name length as an additiona argument */
|
/* Avoid passing file_name length as an additiona argument */
|
||||||
%apply (char *STRING, int LENGTH) { (const char* file_name, const int64_t size_max) };
|
%apply (char *STRING, int LENGTH) { (const char* file_name, const int64_t size_max) };
|
||||||
@ -29,6 +30,9 @@
|
|||||||
|
|
||||||
%cstring_bounded_output(char* function_name, 1024);
|
%cstring_bounded_output(char* function_name, 1024);
|
||||||
%cstring_bounded_output(char* message, 1024);
|
%cstring_bounded_output(char* message, 1024);
|
||||||
|
%cstring_bounded_output(char* const basis_type, 2);
|
||||||
|
/* Required for qmckl_last_error function to work */
|
||||||
|
%cstring_bounded_output(char* buffer, 1024);
|
||||||
|
|
||||||
/* This block is needed make SWIG convert NumPy arrays to/from from the C pointer and size_max argument.
|
/* This block is needed make SWIG convert NumPy arrays to/from from the C pointer and size_max argument.
|
||||||
NOTE: `numpy.i` interface file is not part of SWIG but it is included in the numpy distribution (under numpy/tools/swig/numpy.i)
|
NOTE: `numpy.i` interface file is not part of SWIG but it is included in the numpy distribution (under numpy/tools/swig/numpy.i)
|
||||||
@ -48,6 +52,12 @@ import_array();
|
|||||||
/* Include typemaps generated by the process_header.py script */
|
/* Include typemaps generated by the process_header.py script */
|
||||||
%include "qmckl_include.i"
|
%include "qmckl_include.i"
|
||||||
|
|
||||||
|
/* Some custom array typemaps which are not generated by process_header.py */
|
||||||
|
%apply ( double* IN_ARRAY1 , int64_t DIM1 ) { ( const double * A, const int64_t size_max_A) };
|
||||||
|
%apply ( double* IN_ARRAY1 , int64_t DIM1 ) { ( const double * B, const int64_t size_max_B) };
|
||||||
|
%apply ( double* ARGOUT_ARRAY1 , int64_t DIM1 ) { ( double* const C, const int64_t size_max_C) };
|
||||||
|
%apply ( double* ARGOUT_ARRAY1 , int64_t DIM1 ) { ( double* const B, const int64_t size_max_B) };
|
||||||
|
|
||||||
/* Handle properly get_point */
|
/* Handle properly get_point */
|
||||||
|
|
||||||
|
|
||||||
|
|||||||
@ -33,6 +33,9 @@ assert mo_num == 404
|
|||||||
|
|
||||||
pq.set_electron_coord(ctx, 'T', walk_num, coord)
|
pq.set_electron_coord(ctx, 'T', walk_num, coord)
|
||||||
|
|
||||||
|
ao_type = pq.get_ao_basis_type(ctx)
|
||||||
|
assert 'G' in ao_type
|
||||||
|
|
||||||
size_max = 5*walk_num*elec_num*mo_num
|
size_max = 5*walk_num*elec_num*mo_num
|
||||||
|
|
||||||
mo_vgl = pq.get_mo_basis_mo_vgl(ctx, size_max)
|
mo_vgl = pq.get_mo_basis_mo_vgl(ctx, size_max)
|
||||||
|
|||||||
Loading…
Reference in New Issue
Block a user