1
0
mirror of https://github.com/TREX-CoE/qmckl.git synced 2024-07-18 08:53:47 +02:00

Added kappa (distance rescale factor) to electron and nucleus structs along with tests. #15

This commit is contained in:
vijay gopal chilkuri 2021-05-26 11:59:03 +05:30
parent f416799623
commit ef51934c7c
3 changed files with 247 additions and 48 deletions

View File

@ -117,7 +117,6 @@ typedef struct qmckl_context_struct {
qmckl_ao_basis_struct ao_basis;
/* To be implemented:
qmckl_nucleus_struct nucleus;
qmckl_mo_struct mo;
qmckl_determinant_struct det;
,*/
@ -174,6 +173,7 @@ qmckl_context qmckl_context_check(const qmckl_context context) {
To create a new context, ~qmckl_context_create()~ should be used.
- Upon success, it returns a pointer to a new context with the ~qmckl_context~ type
- It returns ~QMCKL_NULL_CONTEXT~ upon failure to allocate the internal data structure
- A new context always has all its members initialized with a NULL value
# Header
#+begin_src c :comments org :tangle (eval h_func) :exports none
@ -222,8 +222,8 @@ qmckl_context qmckl_context_create() {
ctx->numprec.range = QMCKL_DEFAULT_RANGE;
ctx->ao_basis.uninitialized = (1 << 10) - 1;
ctx->nucleus.uninitialized = (1 << 3) - 1;
ctx->electron.uninitialized = (1 << 2) - 1;
ctx->nucleus.uninitialized = (1 << 4) - 1;
ctx->electron.uninitialized = (1 << 3) - 1;
/* Allocate qmckl_memory_struct */
{

View File

@ -63,19 +63,25 @@ int main() {
The following data stored in the context:
| ~uninitialized~ | int32_t | Keeps bit set for uninitialized data |
| ~num~ | int64_t | Total number of electrons |
| ~up_num~ | int64_t | Number of up-spin electrons |
| ~down_num~ | int64_t | Number of down-spin electrons |
| ~walk_num~ | int64_t | Number of walkers |
| ~provided~ | bool | If true, ~electron~ is valid |
| ~coord_new~ | double[walk_num][3][num] | New set of electron coordinates |
| ~coord_old~ | double[walk_num][3][num] | Old set of electron coordinates |
| ~coord_new_date~ | uint64_t | Last modification date of the coordinates |
| ~ee_distance~ | double[walk_num][num][num] | Electron-electron distances |
| ~ee_distance_date~ | uint64_t | Last modification date of the electron-electron distances |
| ~en_distance~ | double[walk_num][nucl_num][num] | Electron-nucleus distances |
| ~en_distance_date~ | uint64_t | Last modification date of the electron-electron distances |
| ~uninitialized~ | int32_t | Keeps bit set for uninitialized data |
| ~num~ | int64_t | Total number of electrons |
| ~up_num~ | int64_t | Number of up-spin electrons |
| ~down_num~ | int64_t | Number of down-spin electrons |
| ~walk_num~ | int64_t | Number of walkers |
| ~kappa_ee~ | double | The distance scaling factor |
| ~kappa_en~ | double | The distance scaling factor |
| ~provided~ | bool | If true, ~electron~ is valid |
| ~coord_new~ | double[walk_num][3][num] | New set of electron coordinates |
| ~coord_old~ | double[walk_num][3][num] | Old set of electron coordinates |
| ~coord_new_date~ | uint64_t | Last modification date of the coordinates |
| ~ee_distance~ | double[walk_num][num][num] | Electron-electron distances |
| ~ee_distance_date~ | uint64_t | Last modification date of the electron-electron distances |
| ~en_distance~ | double[walk_num][nucl_num][num] | Electron-nucleus distances |
| ~en_distance_date~ | uint64_t | Last modification date of the electron-electron distances |
| ~ee_distance_rescaled~ | double[walk_num][num][num] | Electron-electron distances |
| ~ee_distance_rescaled_date~ | uint64_t | Last modification date of the electron-electron distances |
| ~en_distance_rescaled~ | double[walk_num][nucl_num][num] | Electron-nucleus distances |
| ~en_distance_rescaled_date~ | uint64_t | Last modification date of the electron-electron distances |
** Data structure
@ -85,13 +91,19 @@ typedef struct qmckl_electron_struct {
int64_t up_num;
int64_t down_num;
int64_t walk_num;
double kappa_ee;
double kappa_en;
int64_t coord_new_date;
int64_t ee_distance_date;
int64_t en_distance_date;
int64_t ee_distance_rescaled_date;
int64_t en_distance_rescaled_date;
double* coord_new;
double* coord_old;
double* ee_distance;
double* en_distance;
double* ee_distance_rescaled;
double* en_distance_rescaled;
int32_t uninitialized;
bool provided;
} qmckl_electron_struct;
@ -165,7 +177,7 @@ qmckl_get_electron_num (const qmckl_context context, int64_t* const num) {
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
int32_t mask = 1;
int32_t mask = 1 << 0;
if ( (ctx->electron.uninitialized & mask) != 0) {
return QMCKL_NOT_PROVIDED;
@ -193,7 +205,7 @@ qmckl_get_electron_up_num (const qmckl_context context, int64_t* const up_num) {
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
int32_t mask = 1;
int32_t mask = 1 << 0;
if ( (ctx->electron.uninitialized & mask) != 0) {
return QMCKL_NOT_PROVIDED;
@ -221,7 +233,7 @@ qmckl_get_electron_down_num (const qmckl_context context, int64_t* const down_nu
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
int32_t mask = 1;
int32_t mask = 1 << 0;
if ( (ctx->electron.uninitialized & mask) != 0) {
return QMCKL_NOT_PROVIDED;
@ -260,7 +272,7 @@ qmckl_get_electron_walk_num (const qmckl_context context, int64_t* const walk_nu
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
int32_t mask = 2;
int32_t mask = 1 << 1;
if ( (ctx->electron.uninitialized & mask) != 0) {
return QMCKL_NOT_PROVIDED;
@ -272,6 +284,71 @@ qmckl_get_electron_walk_num (const qmckl_context context, int64_t* const walk_nu
}
#+end_src
*** Scaling factors Kappa
#+begin_src c :comments org :tangle (eval h_func) :exports none
qmckl_exit_code qmckl_get_kappa_ee (const qmckl_context context, double* const kappa_ee);
qmckl_exit_code qmckl_get_kappa_en (const qmckl_context context, double* const kappa_en);
#+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code
qmckl_get_kappa_ee (const qmckl_context context, double* const kappa_ee) {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_INVALID_CONTEXT;
}
if (kappa_ee == NULL) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_2,
"qmckl_get_kappa_ee",
"kappa_ee is a null pointer");
}
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
int32_t mask = 1 << 2;
if ( (ctx->electron.uninitialized & mask) != 0) {
return QMCKL_NOT_PROVIDED;
}
// TODO: assert (ctx->electron.kappa_ee > (double) 0);
,*kappa_ee = ctx->electron.kappa_ee;
return QMCKL_SUCCESS;
}
qmckl_exit_code
qmckl_get_kappa_en (const qmckl_context context, double* const kappa_en) {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_INVALID_CONTEXT;
}
if (kappa_en == NULL) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_2,
"qmckl_get_kappa_en",
"kappa_en is a null pointer");
}
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
int32_t mask = 1 << 2;
if ( (ctx->electron.uninitialized & mask) != 0) {
return QMCKL_NOT_PROVIDED;
}
// TODO: assert (ctx->electron.kappa_en > (double) 0);
,*kappa_en = ctx->electron.kappa_en;
return QMCKL_SUCCESS;
}
#+end_src
*** Electron coordinates
Returns the current electron coordinates. The pointer is assumed
@ -359,6 +436,7 @@ qmckl_get_electron_coord (const qmckl_context context, const char transp, double
#+begin_src c :comments org :tangle (eval h_func)
qmckl_exit_code qmckl_set_electron_num (qmckl_context context, const int64_t up_num, const int64_t down_num);
qmckl_exit_code qmckl_set_kappa (qmckl_context context, const double kappa_ee, const double kappa_en);
qmckl_exit_code qmckl_set_electron_walk_num (qmckl_context context, const int64_t walk_num);
qmckl_exit_code qmckl_set_electron_coord (qmckl_context context, const char transp, const double* coord);
#+end_src
@ -437,7 +515,7 @@ qmckl_set_electron_num(qmckl_context context,
"down_num <= 0");
}
int32_t mask = 1;
int32_t mask = 1 << 0;
ctx->electron.up_num = up_num;
ctx->electron.down_num = down_num;
@ -447,6 +525,7 @@ qmckl_set_electron_num(qmckl_context context,
}
#+end_src
The following function sets the number of walkers.
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code
@ -461,13 +540,46 @@ qmckl_set_electron_walk_num(qmckl_context context, const int64_t walk_num) {
"walk_num <= 0");
}
int32_t mask = 2;
int32_t mask = 1 << 1;
ctx->electron.walk_num = walk_num;
<<post2>>
}
#+end_src
Next we set the rescale parameter for the rescaled distance metric.
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code
qmckl_set_kappa(qmckl_context context,
const double kappa_ee,
const double kappa_en) {
<<pre2>>
// TODO: Check for 0 values
//if (kappa_ee != 0) {
// return qmckl_failwith( context,
// QMCKL_INVALID_ARG_2,
// "qmckl_set_kappa",
// "kappa_ee == 0");
//}
//if (kappa_en <= 0) {
// return qmckl_failwith( context,
// QMCKL_INVALID_ARG_3,
// "qmckl_set_kappa",
// "kappa_en == 0");
//}
int32_t mask = 1 << 2;
ctx->electron.kappa_ee = kappa_ee;
ctx->electron.kappa_en = kappa_en;
<<post2>>
}
#+end_src
The following function sets the electron coordinates of all the
walkers. When this is done, the pointers to the old and new sets
of coordinates are swapped, and the new coordinates are
@ -568,11 +680,14 @@ int64_t walk_num = chbrclf_walk_num;
int64_t elec_num = chbrclf_elec_num;
int64_t elec_up_num = chbrclf_elec_up_num;
int64_t elec_dn_num = chbrclf_elec_dn_num;
double kappa_ee = 1.0; // TODO Get kappa_ee from chbrclf
double kappa_en = 1.0; // TODO Get kappa_en from chbrclf
double* elec_coord = &(chbrclf_elec_coord[0][0][0]);
int64_t nucl_num = chbrclf_nucl_num;
double* charge = chbrclf_charge;
double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
double nucl_kappa = 1.0; // TODO Change get kappa from chbrclf example
/* --- */
@ -607,6 +722,26 @@ rc = qmckl_get_electron_num (context, &n);
assert(rc == QMCKL_SUCCESS);
assert(n == elec_num);
double k_ee;
double k_en;
rc = qmckl_get_kappa_ee (context, &k_ee);
assert(rc == QMCKL_NOT_PROVIDED);
rc = qmckl_get_kappa_en (context, &k_en);
assert(rc == QMCKL_NOT_PROVIDED);
rc = qmckl_set_kappa (context, kappa_ee, kappa_en);
assert(rc == QMCKL_SUCCESS);
assert(!qmckl_electron_provided(context));
rc = qmckl_get_kappa_ee (context, &k_ee);
assert(rc == QMCKL_SUCCESS);
assert(k_ee == kappa_ee);
rc = qmckl_get_kappa_en (context, &k_en);
assert(rc == QMCKL_SUCCESS);
assert(k_en == kappa_en);
int64_t w;
rc = qmckl_get_electron_walk_num (context, &w);
@ -1104,6 +1239,9 @@ assert(qmckl_electron_provided(context));
rc = qmckl_set_nucleus_num (context, nucl_num);
assert(rc == QMCKL_SUCCESS);
rc = qmckl_set_nucleus_kappa (context, nucl_kappa);
assert(rc == QMCKL_SUCCESS);
rc = qmckl_set_nucleus_charge (context, charge);
assert (rc == QMCKL_SUCCESS);
@ -1152,7 +1290,7 @@ assert(fabs(en_distance[1][0][1] - 3.1804527583077356) < 1.e-12);
}
#+end_src
**✸ Compute file names
*** Compute file names
#+begin_src emacs-lisp
; The following is required to compute the file names

View File

@ -62,15 +62,18 @@ int main() {
The following data stored in the context:
| ~uninitialized~ | int32_t | Keeps bit set for uninitialized data |
| ~num~ | int64_t | Total number of nuclei |
| ~provided~ | bool | If true, ~nucleus~ is valid |
| ~charge~ | double[num] | Nuclear charges |
| ~coord~ | double[3][num] | Nuclear coordinates, in transposed format |
| ~nn_distance~ | double[num][num] | Nucleus-nucleus distances |
| ~nn_distance_date~ | int64_t | Date when Nucleus-nucleus distances were computed |
| ~repulsion~ | double | Nuclear repulsion energy |
| ~repulsion_date~ | int64_t | Date when the nuclear repulsion energy was computed |
| ~uninitialized~ | int32_t | Keeps bit set for uninitialized data |
| ~num~ | int64_t | Total number of nuclei |
| ~provided~ | bool | If true, ~nucleus~ is valid |
| ~charge~ | double[num] | Nuclear charges |
| ~coord~ | double[3][num] | Nuclear coordinates, in transposed format |
| ~nn_distance~ | double[num][num] | Nucleus-nucleus distances |
| ~nn_distance_date~ | int64_t | Date when Nucleus-nucleus distances were computed |
| ~nn_distance_rescaled~ | double[num][num] | Nucleus-nucleus rescaled distances |
| ~nn_distance_rescaled_date~ | int64_t | Date when Nucleus-nucleus rescaled distances were computed |
| ~repulsion~ | double | Nuclear repulsion energy |
| ~repulsion_date~ | int64_t | Date when the nuclear repulsion energy was computed |
| ~kappa~ | double | The distance scaling factor |
** Data structure
@ -85,21 +88,23 @@ typedef struct qmckl_nucleus_struct {
double* nn_distance;
double* nn_distance_rescaled;
double repulsion;
double kappa;
int32_t uninitialized;
bool provided;
} qmckl_nucleus_struct;
#+end_src
The ~uninitialized~ integer contains one bit set to one for each
initialization function which has not bee called. It becomes equal
initialization function which has not been called. It becomes equal
to zero after all initialization functions have been called. The
struct is then initialized and ~provided == true~.
** Access functions
#+begin_src c :comments org :tangle (eval h_func) :exports none
qmckl_exit_code qmckl_get_nucleus_num (const qmckl_context context, int64_t* const num);
qmckl_exit_code qmckl_get_nucleus_charge (const qmckl_context context, double* const charge);
qmckl_exit_code qmckl_get_nucleus_kappa (const qmckl_context context, double* const kappa);
qmckl_exit_code qmckl_get_nucleus_coord (const qmckl_context context, const char transp, double* const coord);
#+end_src
@ -181,6 +186,31 @@ qmckl_get_nucleus_charge (const qmckl_context context, double* const charge) {
}
qmckl_exit_code
qmckl_get_nucleus_kappa (const qmckl_context context, double* const kappa) {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_INVALID_CONTEXT;
}
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
int32_t mask = 1 << 2;
if ( (ctx->nucleus.uninitialized & mask) != 0) {
return qmckl_failwith( context,
QMCKL_NOT_PROVIDED,
"qmckl_get_nucleus_kappa",
"nucleus data is not provided");
}
(*kappa) = ctx->nucleus.kappa;
return QMCKL_SUCCESS;
}
qmckl_exit_code
qmckl_get_nucleus_coord (const qmckl_context context, const char transp, double* const coord) {
@ -205,7 +235,7 @@ qmckl_get_nucleus_coord (const qmckl_context context, const char transp, double*
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
int32_t mask = 1 << 2;
int32_t mask = 1 << 3;
if ( (ctx->nucleus.uninitialized & mask) != 0) {
return qmckl_failwith( context,
@ -265,6 +295,7 @@ bool qmckl_nucleus_provided(const qmckl_context context) {
#+begin_src c :comments org :tangle (eval h_func)
qmckl_exit_code qmckl_set_nucleus_num (qmckl_context context, const int64_t num);
qmckl_exit_code qmckl_set_nucleus_charge (qmckl_context context, const double* charge);
qmckl_exit_code qmckl_set_nucleus_kappa (qmckl_context context, const double kappa);
qmckl_exit_code qmckl_set_nucleus_coord (qmckl_context context, const char transp, const double* coord);
#+end_src
@ -349,6 +380,32 @@ qmckl_set_nucleus_charge(qmckl_context context, const double* charge) {
ctx->nucleus.charge= memcpy(ctx->nucleus.charge, charge, num*sizeof(double));
assert (ctx->nucleus.charge != NULL);
<<post2>>
}
#+end_src
The following function 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_kappa(qmckl_context context, const double kappa) {
<<pre2>>
//TODO: Check for small values of kappa
//if (kappa == 0) {
// return qmckl_failwith( context,
// QMCKL_INVALID_ARG_2,
// "qmckl_set_nucleus_kappa",
// "kappa cannot be 0");
//}
int32_t mask = 1 << 2;
int64_t num;
qmckl_exit_code rc;
ctx->nucleus.kappa = kappa;
<<post2>>
}
#+end_src
@ -364,7 +421,7 @@ qmckl_set_nucleus_coord(qmckl_context context, const char transp, const double*
int64_t nucl_num = (int64_t) 0;
qmckl_exit_code rc;
int32_t mask = 1 << 2;
int32_t mask = 1 << 3;
rc = qmckl_get_nucleus_num(context, &nucl_num);
if (rc != QMCKL_SUCCESS) return rc;
@ -405,6 +462,7 @@ qmckl_set_nucleus_coord(qmckl_context context, const char transp, const double*
const int64_t nucl_num = chbrclf_nucl_num;
const double* nucl_charge = chbrclf_charge;
const double* nucl_coord = &(chbrclf_nucl_coord[0][0]);
const double nucl_kappa = 1.0; // TODO Change get kappa from chbrclf example
/* --- */
@ -425,6 +483,19 @@ rc = qmckl_get_nucleus_num (context, &n);
assert(rc == QMCKL_SUCCESS);
assert(n == nucl_num);
double k;
rc = qmckl_get_nucleus_kappa (context, &k);
assert(rc == QMCKL_NOT_PROVIDED);
rc = qmckl_set_nucleus_kappa (context, nucl_kappa);
assert(rc == QMCKL_SUCCESS);
assert(!qmckl_nucleus_provided(context));
rc = qmckl_get_nucleus_kappa (context, &k);
assert(rc == QMCKL_SUCCESS);
assert(k == nucl_kappa);
double nucl_coord2[3*nucl_num];
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2);
@ -476,7 +547,7 @@ assert(qmckl_nucleus_provided(context));
current date is stored.
** Nucleus-nucleus distances
*** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes
@ -816,7 +887,7 @@ qmckl_exit_code qmckl_compute_nn_distance_rescaled (
\]
*** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_nucleus_repulsion(qmckl_context context, double* energy);
#+end_src
@ -960,16 +1031,6 @@ qmckl_exit_code qmckl_compute_nucleus_repulsion (
*** Test
#+begin_src c :tangle (eval c_test)
/* Reference input data */
assert(qmckl_nucleus_provided(context));
double rep;
rc = qmckl_get_nucleus_repulsion(context, &rep);
assert(rep - 318.2309879436158 < 1.e-10);
#+end_src
* End of files :noexport: