TREX/trexio
1
0
Fork 0
mirror of https://github.com/TREX-CoE/trexio.git synced 2024-10-02 06:21:05 +02:00
Code Issues Projects Releases Wiki Activity

Merge branch 'configure' of github.com:TREX-CoE/trexio into configure

Browse Source
This commit is contained in:
Anthony Scemama 2021-11-01 10:16:14 +01:00
commit 3154e01021
10 changed files with 202 additions and 95 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

7
ChangeLog
View File

@ -1,6 +1,13 @@
CHANGES
=======
2.0
---
- Fixed data structure of the <ecp> group in trex.org in order to simplify the
access to both local and non-local quantities
- Adapted data structure of the <basis> group to be consistent with the <ecp>
1.1
---

10
README.md
View File

@ -50,9 +50,15 @@ TREX library for efficient I/O.
The `make install` command takes care of installing the TREXIO shared library on the user machine.
Once installed, add `-ltrexio` to the list of compiler options.
In some cases (e.g. when using custom `prefix` during configuration), the TREXIO library might end up installed in a directory, which is absent in the default `$LIBRARY_PATH`.
In order to link the program against TREXIO, the search paths in the current shell can be modified as follows: `export LIBRARY_PATH=$LIBRARY_PATH:<path_to_trexio>/lib` (same holds for `$LD_LIBRARY_PATH`). Do not forget to change `<path_to_trexio>`.
If your compilation relies on some build tools (like Autotools or CMake), feel free to use the built-in functions to locate and link external dependencies automatically.
In order to link the program against TREXIO, the search paths can be modified as follows:
`export LIBRARY_PATH=$LIBRARY_PATH:<path_to_trexio>/lib`
(same holds for `$LD_LIBRARY_PATH`). The `<path_to_trexio>` has to be replaced by the prefix used during the installation.
If your compilation relies on some build tools (like Autotools or CMake), feel free to use the built-in solutions to locate and link external dependencies automatically.
In Fortran applications, make sure that the `trexio_f.f90` module file is included in the source tree.
You might have to manually copy it into your program source directory.

2
configure.ac
View File

@ -2,7 +2,7 @@
# Process this file with autoconf to produce a configure script.
AC_PREREQ([2.69])
AC_INIT([trexio], [1.1.0], [https://github.com/TREX-CoE/trexio/issues])
AC_INIT([trexio], [2.0.0], [https://github.com/TREX-CoE/trexio/issues])
AM_INIT_AUTOMAKE([subdir-objects color-tests parallel-tests silent-rules 1.11])
AM_MAINTAINER_MODE()

2
python/pytrexio/_version.py
View File

@ -1 +1 @@
__version__ = "0.1.0"
__version__ = "1.0.0"

18
python/test/test_api.py
View File

@ -71,11 +71,15 @@ charges_np = np.array(charges, dtype=np.int32)
# from the size of the list/array by SWIG using typemaps from numpy.i
trexio.write_nucleus_charge(test_file, charges_np)
basis_shell_num = 24
# initialize arrays of nuclear indices as a list and convert it to numpy array
indices = [i for i in range(nucleus_num)]
indices = [i for i in range(basis_shell_num)]
# type cast is important here because by default numpy transforms a list of integers into int64 array
indices_np = np.array(indices, dtype=np.int64)
# first write basis_shell_num because it is needed to check dimensions of basis_nucleus_index in TREXIO >= 2.0.0
trexio.write_basis_shell_num(test_file, basis_shell_num)
# function call below works with both lists and numpy arrays, dimension needed for memory-safety is derived
# from the size of the list/array by SWIG using typemacs from numpy.i
trexio.write_basis_nucleus_index(test_file, indices_np)
@ -160,20 +164,20 @@ except trexio.Error:
print("Unsafe call to safe API: checked")
# safe call to read array of int values (nuclear indices)
rindices_np_16 = trexio.read_basis_nucleus_index(test_file2, dim=nucleus_num, dtype=np.int16)
rindices_np_16 = trexio.read_basis_nucleus_index(test_file2, dim=basis_shell_num, dtype=np.int16)
assert rindices_np_16.dtype is np.dtype(np.int16)
for i in range(nucleus_num):
for i in range(basis_shell_num):
assert rindices_np_16[i]==indices_np[i]
rindices_np_32 = trexio.read_basis_nucleus_index(test_file2, dim=nucleus_num, dtype=np.int32)
rindices_np_32 = trexio.read_basis_nucleus_index(test_file2, dim=basis_shell_num, dtype=np.int32)
assert rindices_np_32.dtype is np.dtype(np.int32)
for i in range(nucleus_num):
for i in range(basis_shell_num):
assert rindices_np_32[i]==indices_np[i]
rindices_np_64 = trexio.read_basis_nucleus_index(test_file2)
assert rindices_np_64.dtype is np.dtype(np.int64)
assert rindices_np_64.size==nucleus_num
for i in range(nucleus_num):
assert rindices_np_64.size==basis_shell_num
for i in range(basis_shell_num):
assert rindices_np_64[i]==indices_np[i]
# read nuclear coordinates without providing optional argument dim

13
python/test/test_pytrexio.py
View File

@ -59,11 +59,16 @@ assert rc==0
# charges[i] = 1.
#rc = trexio_write_nucleus_charge(test_file, charges)
basis_num = 24
# initialize arrays of nuclear indices as a list and convert it to numpy array
indices = [i for i in range(nucleus_num)]
indices = [i for i in range(basis_num)]
# type cast is important here because by default numpy transforms a list of integers into int64 array
indices_np = np.array(indices, dtype=np.int32)
# first write basis_num because it is needed to check dimensions of basis_nucleus_index in TREXIO >= 2.0.0
rc = trexio_write_basis_shell_num(test_file, basis_num)
assert rc==0
# function call below works with both lists and numpy arrays, dimension needed for memory-safety is derived
# from the size of the list/array by SWIG using typemacs from numpy.i
rc = trexio_write_safe_basis_nucleus_index(test_file, indices_np)
@ -124,8 +129,12 @@ assert rc==23
#for i in range(nucleus_num):
# assert charges2[i]==charges[i]
result_basis = trexio_read_basis_shell_num(test_file2)
assert result[0]==0
assert result[1]==basis_num
# safe call to read_safe array of int values
rc, rindices_np = trexio_read_safe_basis_nucleus_index(test_file2, nucleus_num)
rc, rindices_np = trexio_read_safe_basis_nucleus_index(test_file2, basis_num)
assert rc==0
assert rindices_np.dtype is np.dtype(np.int32)
for i in range(nucleus_num):

4
tests/io_dset_int_hdf5.c
View File

@ -25,7 +25,7 @@ static int test_write_dset (const char* file_name, const back_end_t backend) {
assert (file != NULL);
// write numerical attribute in an empty file
rc = trexio_write_nucleus_num(file, num);
rc = trexio_write_basis_shell_num(file, num);
assert (rc == TREXIO_SUCCESS);
// write numerical (integer) dataset in a file
@ -91,7 +91,7 @@ static int test_read_dset (const char* file_name, const back_end_t backend) {
assert (file != NULL);
// read numerical attribute from the file
rc = trexio_read_nucleus_num(file, &num);
rc = trexio_read_basis_shell_num(file, &num);
assert (rc == TREXIO_SUCCESS);
assert (num == 12);

4
tests/io_dset_int_text.c
View File

@ -25,7 +25,7 @@ static int test_write_dset (const char* file_name, const back_end_t backend) {
assert (file != NULL);
// write numerical attribute in an empty file
rc = trexio_write_nucleus_num(file, num);
rc = trexio_write_basis_shell_num(file, num);
assert (rc == TREXIO_SUCCESS);
// write numerical (integer) dataset in a file
@ -91,7 +91,7 @@ static int test_read_dset (const char* file_name, const back_end_t backend) {
assert (file != NULL);
// read numerical attribute from the file
rc = trexio_read_nucleus_num(file, &num);
rc = trexio_read_basis_shell_num(file, &num);
assert (rc == TREXIO_SUCCESS);
assert (num == 12);

14
tests/test_f.f90
View File

@ -44,9 +44,9 @@ subroutine test_write(file_name, back_end)
integer :: rc = 1
integer :: num
integer :: num, basis_shell_num
integer :: basis_nucleus_index(12)
integer :: basis_nucleus_index(24)
double precision :: charge(12)
double precision :: coord(3,12)
@ -70,7 +70,8 @@ subroutine test_write(file_name, back_end)
0.00000000d0, 2.47304151d0 , 0.00000000d0 /), &
shape(coord) )
basis_nucleus_index = (/ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12 /)
basis_shell_num = 24
basis_nucleus_index = (/ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24 /)
label = [character(len=8) :: 'C', 'Na','C', 'C 66', 'C','C', 'H 99', 'Ru', 'H', 'H', 'H', 'H' ]
@ -104,6 +105,8 @@ subroutine test_write(file_name, back_end)
deallocate(sym_str)
call trexio_assert(rc, TREXIO_SUCCESS, 'SUCCESS WRITE POINT GROUP')
rc = trexio_write_basis_shell_num(trex_file, basis_shell_num)
call trexio_assert(rc, TREXIO_SUCCESS, 'SUCCESS WRITE BASIS NUM')
rc = trexio_write_basis_nucleus_index(trex_file, basis_nucleus_index)
call trexio_assert(rc, TREXIO_SUCCESS, 'SUCCESS WRITE INDEX')
@ -137,9 +140,9 @@ subroutine test_read(file_name, back_end)
integer :: i, j, k, ind, offset, flag
integer :: rc = 1
integer :: num, num_read
integer :: num, num_read, basis_shell_num
integer :: basis_nucleus_index(12)
integer :: basis_nucleus_index(24)
double precision :: charge(12)
double precision :: coord(3,12)
@ -152,6 +155,7 @@ subroutine test_read(file_name, back_end)
character*(128) :: str
num = 12
basis_shell_num = 24
! ================= START OF TEST ===================== !

195
trex.org
View File

@ -118,15 +118,18 @@ arrays are 0-based. Hence, we introduce the ~index~ type which is an
* Effective core potentials (ecp group)
An effective core potential (ECP) $V_A^{\text{ECP}}$ replacing the
core electrons of atom $A$ is expressed as
core electrons of atom $A$ can be expressed as
\[
V_A^{\text{ECP}} =
V_{A \ell_{\max}} +
\sum_{\ell=0}^{\ell_{\max} -1}
V_{A \ell_{\max}+1} +
\sum_{\ell=0}^{\ell_{\max}}
\sum_{m=-\ell}^{\ell} | Y_{\ell m} \rangle \left[
V_{A \ell} - V_{A \ell_{\max}} \right] \langle Y_{\ell m} |
V_{A \ell} - V_{A \ell_{\max}+1} \right] \langle Y_{\ell m} |
\]
The first term in the equation above is sometimes attributed to the local channel,
while the remaining terms correspond to the non-local channel projections.
The functions $V_{A\ell}$ are parameterized as:
\[
V_{A \ell}(\mathbf{r}) =
@ -135,23 +138,34 @@ arrays are 0-based. Hence, we introduce the ~index~ type which is an
e^{-\alpha_{A q \ell} |\mathbf{r}-\mathbf{R}_{A}|^2 }
\]
See http://dx.doi.org/10.1063/1.4984046 for more info.
See http://dx.doi.org/10.1063/1.4984046 or https://doi.org/10.1063/1.5121006 for more info.
#+NAME: ecp
| Variable | Type | Dimensions | Description |
|-----------------------+---------+------------------------------------------+----------------------------------------------------------------------------------------------|
| ~lmax_plus_1~ | ~int~ | ~(nucleus.num)~ | $\ell_{\max} + 1$, one higher than the maximum angular momentum in the removed core orbitals |
| ~z_core~ | ~float~ | ~(nucleus.num)~ | Charges to remove |
| ~local_n~ | ~int~ | ~(nucleus.num)~ | Number of local functions $N_{q \ell}$ |
| ~local_num_n_max~ | ~dim~ | | Maximum value of ~local_n~, used for dimensioning arrays |
| ~local_exponent~ | ~float~ | ~(ecp.local_num_n_max, nucleus.num)~ | $\alpha_{A q \ell_{\max}}$ |
| ~local_coef~ | ~float~ | ~(ecp.local_num_n_max, nucleus.num)~ | $\beta_{A q \ell_{\max}}$ |
| ~local_power~ | ~int~ | ~(ecp.local_num_n_max, nucleus.num)~ | $n_{A q \ell_{\max}}$ |
| ~non_local_n~ | ~int~ | ~(nucleus.num)~ | $N_{q \ell_{\max}}$ |
| ~non_local_num_n_max~ | ~dim~ | | Maximum value of ~non_local_n~, used for dimensioning arrays |
| ~non_local_exponent~ | ~float~ | ~(ecp.non_local_num_n_max, nucleus.num)~ | $\alpha_{A q \ell}$ |
| ~non_local_coef~ | ~float~ | ~(ecp.non_local_num_n_max, nucleus.num)~ | $\beta_{A q \ell}$ |
| ~non_local_power~ | ~int~ | ~(ecp.non_local_num_n_max, nucleus.num)~ | $n_{A q \ell}$ |
|----------------------+---------+-----------------+----------------------------------------------------------------------------------------|
| ~max_ang_mom_plus_1~ | ~int~ | ~(nucleus.num)~ | $\ell_{\max}+1$, one higher than the max angular momentum in the removed core orbitals |
| ~z_core~ | ~int~ | ~(nucleus.num)~ | Number of core electrons to remove per atom |
| ~num~ | ~dim~ | | Total number of ECP functions for all atoms and all values of $\ell$ |
| ~ang_mom~ | ~int~ | ~(ecp.num)~ | One-to-one correspondence between ECP items and the angular momentum $\ell$ |
| ~nucleus_index~ | ~index~ | ~(ecp.num)~ | One-to-one correspondence between ECP items and the atom index |
| ~exponent~ | ~float~ | ~(ecp.num)~ | $\alpha_{A q \ell}$ all ECP exponents |
| ~coefficient~ | ~float~ | ~(ecp.num)~ | $\beta_{A q \ell}$ all ECP coefficients |
| ~power~ | ~int~ | ~(ecp.num)~ | $n_{A q \ell}$ all ECP powers |
There might be some confusion in the meaning of the $\ell_{\max}$.
It can be attributed to the maximum angular momentum occupied
in the core orbitals, which are removed by the ECP.
On the other hand, it can be attributed to the maximum angular momentum of the
ECP that replaces the core electrons.
*Note*, that the latter $\ell_{\max}$ is always higher by 1 than the former.
*Note for developers*: avoid having variables with similar prefix in their name.
HDF5 back end might cause issues due to the way ~find_dataset~ function works.
For example, in the ECP group we use ~max_ang_mom~ and not ~ang_mom_max~.
The latter causes issues when written before ~ang_mom~ in the TREXIO file.
#+CALL: json(data=ecp, title="ecp")
@ -159,22 +173,74 @@ arrays are 0-based. Hence, we introduce the ~index~ type which is an
:RESULTS:
#+begin_src python :tangle trex.json
"ecp": {
"lmax_plus_1" : [ "int" , [ "nucleus.num" ] ]
, "z_core" : [ "float", [ "nucleus.num" ] ]
, "local_n" : [ "int" , [ "nucleus.num" ] ]
, "local_num_n_max" : [ "dim" , [] ]
, "local_exponent" : [ "float", [ "nucleus.num", "ecp.local_num_n_max" ] ]
, "local_coef" : [ "float", [ "nucleus.num", "ecp.local_num_n_max" ] ]
, "local_power" : [ "int" , [ "nucleus.num", "ecp.local_num_n_max" ] ]
, "non_local_n" : [ "int" , [ "nucleus.num" ] ]
, "non_local_num_n_max" : [ "dim" , [] ]
, "non_local_exponent" : [ "float", [ "nucleus.num", "ecp.non_local_num_n_max" ] ]
, "non_local_coef" : [ "float", [ "nucleus.num", "ecp.non_local_num_n_max" ] ]
, "non_local_power" : [ "int" , [ "nucleus.num", "ecp.non_local_num_n_max" ] ]
"max_ang_mom_plus_1" : [ "int" , [ "nucleus.num" ] ]
, "z_core" : [ "int" , [ "nucleus.num" ] ]
, "num" : [ "dim" , [] ]
, "ang_mom" : [ "int" , [ "ecp.num" ] ]
, "nucleus_index" : [ "index", [ "ecp.num" ] ]
, "exponent" : [ "float", [ "ecp.num" ] ]
, "coefficient" : [ "float", [ "ecp.num" ] ]
, "power" : [ "int" , [ "ecp.num" ] ]
} ,
#+end_src
:END:
** Example
For example, consider H_2 molecule with the following
[[https://pseudopotentiallibrary.org/recipes/H/ccECP/H.ccECP.gamess][effective core potential]]
(in GAMESS input format for the H atom):
#+BEGIN_EXAMPLE
H-ccECP GEN 0 1
3
1.00000000000000 1 21.24359508259891
21.24359508259891 3 21.24359508259891
-10.85192405303825 2 21.77696655044365
1
0.00000000000000 2 1.000000000000000
#+END_EXAMPLE
In TREXIO representation this would be:
#+BEGIN_EXAMPLE
num = 8
# lmax+1 per atom
max_ang_mom_plus_1 = [ 1, 1 ]
# number of core electrons to remove per atom
zcore = [ 0, 0 ]
# first 4 ECP elements correspond to the first H atom ; the remaining 4 elements are for the second H atom
nucleus_index = [
0, 0, 0, 0,
1, 1, 1, 1
]
# 3 first ECP elements correspond to potential of the P orbital (l=1), then 1 element for the S orbital (l=0) ; similar for the second H atom
ang_mom = [
1, 1, 1, 0,
1, 1, 1, 0
]
# ECP quantities that can be attributed to atoms and/or angular momenta based on the aforementioned ecp_nucleus and ecp_ang_mom arrays
coefficient = [
1.00000000000000, 21.24359508259891, -10.85192405303825, 0.00000000000000,
1.00000000000000, 21.24359508259891, -10.85192405303825, 0.00000000000000
]
exponent = [
21.24359508259891, 21.24359508259891, 21.77696655044365, 1.000000000000000,
21.24359508259891, 21.24359508259891, 21.77696655044365, 1.000000000000000
]
power = [
-1, 1, 0, 0,
-1, 1, 0, 0
]
#+END_EXAMPLE
* Basis set (basis group)
We consider here basis functions centered on nuclei. Hence, we enable
@ -216,16 +282,14 @@ arrays are 0-based. Hence, we introduce the ~index~ type which is an
#+NAME: basis
| Variable | Type | Dimensions | Description |
|---------------------+---------+--------------------+----------------------------------------------------------|
|-----------------+---------+---------------------+--------------------------------------------------------------|
| ~type~ | ~str~ | | Type of basis set: "Gaussian" or "Slater" |
| ~num~ | ~dim~ | | Total Number of shells |
| ~prim_num~ | ~dim~ | | Total number of primitives |
| ~nucleus_index~ | ~index~ | ~(nucleus.num)~ | Index of the first shell of each nucleus ($A$) |
| ~nucleus_shell_num~ | ~int~ | ~(nucleus.num)~ | Number of shells for each nucleus |
| ~shell_ang_mom~ | ~int~ | ~(basis.num)~ | Angular momentum ~0:S, 1:P, 2:D, ...~ |
| ~shell_prim_num~ | ~int~ | ~(basis.num)~ | Number of primitives in the shell ($N_{\text{prim}}$) |
| ~shell_factor~ | ~float~ | ~(basis.num)~ | Normalization factor of the shell ($\mathcal{N}_s$) |
| ~shell_prim_index~ | ~index~ | ~(basis.num)~ | Index of the first primitive in the complete list |
| ~shell_num~ | ~dim~ | | Total number of shells |
| ~nucleus_index~ | ~index~ | ~(basis.shell_num)~ | One-to-one correspondence between shells and atomic indices |
| ~shell_ang_mom~ | ~int~ | ~(basis.shell_num)~ | One-to-one correspondence between shells and angular momenta |
| ~shell_factor~ | ~float~ | ~(basis.shell_num)~ | Normalization factor of each shell ($\mathcal{N}_s$) |
| ~shell_index~ | ~index~ | ~(basis.prim_num)~ | One-to-one correspondence between primitives and shell index |
| ~exponent~ | ~float~ | ~(basis.prim_num)~ | Exponents of the primitives ($\gamma_{ks}$) |
| ~coefficient~ | ~float~ | ~(basis.prim_num)~ | Coefficients of the primitives ($a_{ks}$) |
| ~prim_factor~ | ~float~ | ~(basis.prim_num)~ | Normalization coefficients for the primitives ($f_{ks}$) |
@ -237,14 +301,12 @@ arrays are 0-based. Hence, we introduce the ~index~ type which is an
#+begin_src python :tangle trex.json
"basis": {
"type" : [ "str" , [] ]
, "num" : [ "dim" , [] ]
, "prim_num" : [ "dim" , [] ]
, "nucleus_index" : [ "index", [ "nucleus.num" ] ]
, "nucleus_shell_num" : [ "int" , [ "nucleus.num" ] ]
, "shell_ang_mom" : [ "int" , [ "basis.num" ] ]
, "shell_prim_num" : [ "int" , [ "basis.num" ] ]
, "shell_factor" : [ "float", [ "basis.num" ] ]
, "shell_prim_index" : [ "index", [ "basis.num" ] ]
, "shell_num" : [ "dim" , [] ]
, "nucleus_index" : [ "index", [ "basis.shell_num" ] ]
, "shell_ang_mom" : [ "int" , [ "basis.shell_num" ] ]
, "shell_factor" : [ "float", [ "basis.shell_num" ] ]
, "shell_index" : [ "index", [ "basis.prim_num" ] ]
, "exponent" : [ "float", [ "basis.prim_num" ] ]
, "coefficient" : [ "float", [ "basis.prim_num" ] ]
, "prim_factor" : [ "float", [ "basis.prim_num" ] ]
@ -252,6 +314,8 @@ arrays are 0-based. Hence, we introduce the ~index~ type which is an
#+end_src
:END:
** Example
For example, consider H_2 with the following basis set (in GAMESS
format), where both the AOs and primitives are considered normalized:
@ -272,31 +336,44 @@ P 1
P 1
1 3.880000E-01 1.000000E+00
D 1
1 1.057000E+00 1.0000000
1 1.057000E+00 1.000000E+00
#+END_EXAMPLE
we have:
In TREXIO representaion we have:
#+BEGIN_EXAMPLE
type = "Gaussian"
num = 12
prim_num = 20
shell_num = 12
nucleus_index = [0 , 6]
shell_ang_mom = [0 , 0 , 0 , 1 , 1 , 2 , 0 , 0 , 0 , 1 , 1 , 2 ]
shell_prim_num = [5 , 1 , 1 , 1 , 1 , 1 , 5 , 1 , 1 , 1 , 1 , 1 ]
shell_prim_index = [0 , 5 , 6 , 7 , 8 , 9 , 10, 15, 16, 17, 18, 19]
shell_factor = [1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.]
# 6 shells per H atom
nucleus_index =
[ 0, 0, 0, 0, 0, 0,
1, 1, 1, 1, 1, 1 ]
# 3 shells in S (l=0), 2 in P (l=1), 1 in D (l=2)
shell_ang_mom =
[ 0, 0, 0, 1, 1, 2,
0, 0, 0, 1, 1, 2 ]
# no need to renormalize shells
shell_factor =
[ 1., 1., 1., 1., 1., 1.,
1., 1., 1., 1., 1., 1. ]
# 5 primitives for the first S shell and then 1 primitive per remaining shells in each H atom
shell_index =
[ 0, 0, 0, 0, 0, 1, 2, 3, 4, 5,
6, 6, 6, 6, 6, 7, 8, 9, 10, 11 ]
# parameters of the primitives (10 per H atom)
exponent =
[ 33.87, 5.095, 1.159, 0.3258, 0.1027, 0.3258, 0.1027, 1.407,
0.388, 1.057, 33.87, 5.095, 1.159, 0.3258, 0.1027, 0.3258, 0.1027, 1.407,
0.388, 1.057]
[ 33.87, 5.095, 1.159, 0.3258, 0.1027, 0.3258, 0.1027, 1.407, 0.388, 1.057,
33.87, 5.095, 1.159, 0.3258, 0.1027, 0.3258, 0.1027, 1.407, 0.388, 1.057 ]
coefficient =
[ 0.006068, 0.045308, 0.202822, 0.503903, 0.383421, 1.0, 1.0,
1.0, 1.0, 1.0, 0.006068, 0.045308, 0.202822, 0.503903, 0.383421, 1.0, 1.0,
1.0, 1.0, 1.0]
[ 0.006068, 0.045308, 0.202822, 0.503903, 0.383421, 1.0, 1.0, 1.0, 1.0, 1.0,
0.006068, 0.045308, 0.202822, 0.503903, 0.383421, 1.0, 1.0, 1.0, 1.0, 1.0 ]
prim_factor =
[ 1.0006253235944540e+01, 2.4169531573445120e+00, 7.9610924849766440e-01