theochem / fanpy Goto Github PK

Noise is nice to add to the parameters to poke the wavefunction/hamiltonian out of a local minimum near the initial guess. However, this may be something that the user control directly rather than a behaviour triggered by a flag. It doesn't seem useful enough to be a default behaviour and lacks control over the amount of noise added.

There are a couple of options:

1. Remove behaviour for adding noise and the user must manually add noise
2. Remove behaviour for adding noise in the add_params and make method for adding noise
3. Keep original

So right now in APIG and its derivatives, we sort creation operators by their spatial orbitals
i.e. :math:a_1^\dagger a_\bar{1}^\dagger a_2^\dagger a_\bar{2}^\dagger \dots
But we are using indices sorted into alpha and beta blocks
i.e. :math:a_1^\dagger a_{k+1}^\dagger a_2^\dagger a_{k+2}^\dagger \dots
where indices greater than $k$ are beta orbitals.

There may be (or already is) a bug in the future if we continue to assume that
the two are interchangeable.
To convert from one convention to another, we sometimes get a sign change.
(I will use notation :math:142536 to describe
:math:a_1^\dagger a_4^\dagger a_2^\dagger a_5^\dagger a_3^\dagger a_6^\dagger.

For no pairs of electrons, we need an even permutation.
For one pair of electrons, we need an even permutation.
For two pairs of electrons, we need an odd permutation (1324 to 1234).
For three pairs of electrons, we need an odd permutation (142536 to 124536 to 124356 to 123456).
For four pairs of electrons, we need an even permutation (15263748 to 12563748 to 12536748 to 12356748 to 12356478 to 12354678 to 12345678)
and so on.

The signature seems to be :math:(-1)^{(P//2)%2} where :math:P is the number of
electron pairs. We can prove this easily (probably).

I don't this this is a huge problem right now (some of our wavefunctions might need
to be multiplied by -1 to convert to the other convention). As long as we are not
comparing one wavefunction in one convention to another in another, there should
be no problem. (but we probably will do this in the future).

Following changes are proposed for improving usability and ease of development

• a script for running "main" calculations
• a script for generating the code for running a calculation
• remove dtype attribute in BaseWavefunction
• split off cache functions (_olp and _olp_deriv) within load_cache of wavefunction classes to outside (for ease of testing) and wrap these functions within load_cache to remove self from parameters.
• store parameter shape/size as an attribute rather than relying on template_params? add param shape and size as properties that depend on template_params (and overwrite if desired)?
• wrap OneRefApprox around geminal instances rather than using multiple inheritances
• rename wavefunction class names and organize structures

Since cardboardlint checks lints for just the code that has been changed (with master), there's quite a bit of a mess in the master branch. This should be cleaned up.

Intermediate normalization adds a small complexity to the APIG code. We need to divide the first P columns with an appropriate normalization constant. There may be some complications also if the selected determinant is not representative of APIG wavefunction.

However, without intermediate normalization and non pair excitations, we cannot get enough equations to solve the projected schrodinger equation. If we have 1 electron pair and K orbitals, then we can only get K Slater determinants using pair wise excitations. Including the energy, we have K+1 unknowns.

There's a lot of this parent-children business going on (children having children). It would be nice if we can "import" documentation of the parent and update/append to it.

I installed the package using conda and running the nosetests I get the following error:

======================================================================
ERROR: Test RankTwoGeminal.template_params.
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/home/sfias/miniconda3/envs/fanci/lib/python3.6/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/mnt/c/Users/stijn/horton/olsens/wfns/wfn/geminal/test/test_rank2_approx.py", line 80, in test_rank2_geminal_template_params
    template = test.template_params
  File "/mnt/c/Users/stijn/horton/olsens/wfns/wfn/geminal/rank2_approx.py", line 73, in template_params
    return full_to_rank2(template, rmsd=0.01)
  File "/mnt/c/Users/stijn/horton/olsens/wfns/wfn/geminal/rank2_approx.py", line 446, in full_to_rank2
    'matrix'.format(deviation))
ValueError: Rank-2 coefficient matrix has RMSD of 0.20454745240896396 with the full-rank coefficient matrix

----------------------------------------------------------------------

my conda list is the following

# packages in environment at /home/sfias/miniconda3/envs/fanci:
#
# Name                    Version                   Build  Channel
ca-certificates           2017.08.26           h1d4fec5_0
certifi                   2018.1.18                py36_0
cma                       2.5.3                     <pip>
gmp                       6.1.2                h6c8ec71_1
gmpy2                     2.0.8            py36hc8893dd_2
intel-openmp              2018.0.0             hc7b2577_8
libedit                   3.1                  heed3624_0
libffi                    3.2.1                hd88cf55_4
libgcc-ng                 7.2.0                h7cc24e2_2
libgfortran-ng            7.2.0                h9f7466a_2
libstdcxx-ng              7.2.0                h7a57d05_2
mkl                       2018.0.1             h19d6760_4
mpc                       1.0.3                hec55b23_5
mpfr                      3.1.5                h11a74b3_2
ncurses                   6.0                  h9df7e31_2
nose                      1.3.7            py36hcdf7029_2
numpy                     1.14.1           py36h3dfced4_1
openssl                   1.0.2n               hb7f436b_0
pip                       9.0.1                    py36_5
python                    3.6.4                hc3d631a_1
readline                  7.0                  ha6073c6_4
scipy                     1.0.0            py36hbf646e7_0
setuptools                38.5.1                   py36_0
sqlite                    3.22.0               h1bed415_0
tk                        8.6.7                hc745277_3
wfns                      0.0.0                     <pip>
wheel                     0.30.0           py36hfd4bba0_1
xz                        5.2.3                h55aa19d_2
zlib                      1.2.11               ha838bed_2

I'm not sure if I missed something...

The methods _integrate_sd_sd_deriv in UnrestrictedChemicalHamiltonian and RestrictedChemicalHamiltonian have lots of repetitions.

They can be shortened with impure functions that add appropriate values to the variables outside its scope as a side effect. However, this may cause the code to be a little less readable as it deviates from source material (notes on orbital optimization).

I don't know how to proceed with this one for now.

When running AP1roG on H2 with a minimal basis it doesn't converge to the right energy. @kimt33 and I figured out this is because the number of Slater Determinants in the projection space is set to the number of parameters - 1 (because we have one equation from the normalization), whereas AP1roG is already normalised by construction. For AP1roG the projection space needs to be set to the number of parameters - 0.

Needs to be optimized more, while retaining good convergence.

Scripts do not use random initial guess by default, I think. At least the option to add random noise to the parameters should be added.

Some of the module and class names are quite long. e.g. wavefunction, base_wavefunction. These should probably be changed for something shorter.

Doesn't work with zeta and/or etas

Tox can support more linters.

Tests take way to long to finish. The slow tests should be identified and moved labelled (and skipped) appropriately

I came across a few things in the slater module. Please let me know what you think, so I can make the needed changes accordingly:

1. Currently, slater.annihilate(0b00110, 1, 6) returns None because one of the specified indices is not occupied. Isn't it better to have slater.annihilate only return None when none of the specified indices are occupied? In this way slater.annihilate(0b00110, 1, 6) or slater.annihilate(0b00110, 6, 1) will return 0b100.
   The same argument applies to the slater.create function.
2. Function slater.ground only works for the closed-shell systems, right? Isn't it useful to have a more general function taking number of occupied alpha and beta orbitals as arguments.

• excitation operators on wavefunction
• CC operator on wavefunction

Travis yaml was based on the attempt (https://github.com/theochem/cookiecutter-python) to standardize workflow in theochem and quantumelephant repo's. It has not yet been accepted within these organizations and its full set of features have not yet been used. It does, however, add some complexity to the travis yaml.

I think that until cookiecutter is well supported (and maintained), we should use a simpler travis yaml file.

Intermediate normalization adds a small complexity to the APIG code. We need to divide the first P columns with an appropriate normalization constant. There may be some complications also if the selected determinant is not representative of APIG wavefunction.

However, without intermediate normalization and non pair excitations, we cannot get enough equations to solve the projected Schrodinger equation. If we have 1 electron pair and K orbitals, then we can only get K Slater determinants using pair wise excitations. Including the energy, we have K+1 unknowns.

If we want HORTON to be one optional interface to Olsens, this should be fixed.

We should probably copy HORTON's finite difference approximation Jacobian-verifier to Olsens and improve it. It's finicky, and it's difficult even to verify some benchmark functions (SciPy's Rosenbrock implementation). Finicky-ness is normal for finite-difference methods, but I feel like it could be improved/slightly more automated, and it would be useful for other analytical Jacobians we create in the future.

In general, permanents, I think, can be applied to non-square matrices. The determinants too can be applied (except that they are always zero). There is this assert statement that requires the matrix to be square.

I can see why the Geminal.permanent can be a static method and that in the case of APIG, we will always only have permanents of square matrices.

I guess I'm just wondering if this needs to be moved to a "math" location. And if it does, then we should maybe remove the square matrix restriction.

Setting the nuclear-nuclear repulsion in the Hamiltonian does not actually do anything.

It should add it to the computed energies in the objective.

1. Button select_orbitals does not actually do anything. So far we have not implemented the wavefunctions that require you to select orbitals. (e.g. CAS, APsetG, APG)
2. OrbitalSelectionDialog cannot be resized b/c the matplotlib object behaves weirdly when you resize the window. There should also be some way to navigate through the graph.

• Set up linters. See https://github.com/theochem/cardboardlint
• Cleanup code
• Nosetests

The integrals of the Hamiltonian, integrate_sd_sd and integrate_wfn_sd, are returned as tuples of one electron integral, coulomb integral, and the exchange integral. It's a bit of a pain to add together these values.

The returned value can either be the sum or the numpy array of the three integrals. I'm not sure exactly when it would useful to have this information, but it might not be a bad idea to keep them separated. Numpy array should make it easier to manipulate these integrals.

Only the cma solver shows the information on each iteration, whereas the others don't print anything. It will be a nice information to have (especially when deciding to terminate the calculation or restart it with modified solver parameters).

So right now, we have this assumed Slater determinant structure where the even position (starting from zero) from the right describe the alpha orbitals, and the odd position from the right describes the beta orbitals. e.g. 0011 describes the first alpha and beta orbitals being used to create the Slater determinant, 0101 describes the first and second alpha orbitals.

I think we need a Slater determinant class to contain the Slater determinant methods. We need to use the Slater determinants for the projected schrodinger method, and I think it would be easier to maintain the code if we have a separate SD class, rather than hard coding it into the projected schrodinger method. This class does not need to have orbital information (it can just be SDIndecies, or something). But we can easily extend on this to get a "full" SD by having the orbital info elsewhere.

Also, I'm not too sure about the alpha beta indexing. If we decide to use the SD indexing to calculate integrals (overlap, one electron, two electron) then it would be harder to exploit spin symmetry if the alpha and beta orbitals are shuffled. i.e. We need to shuffle the alphas and betas to one side to exploit the block diagonal structure (or block tensor structure). It should be trivial-ish to shuffle them, but current structure will involve an extra step. Also, (much more trivially) we assume that the number of alpha and beta orbitals are always equal. If anyone wants to implement this, it will become quite a bit of a pain in both structures. We can do something "hacky", like have a dummy beta orbital (that's just a zero). I'm not sure when someone might need to do this, but it might happen.

I would like to have an implementation of CI or VB code eventually for some quasibasis stuff, so it will be useful there also.

There should be a method for calculating the integral as a function of the transformation matrix or of the jacobi rotation indices and angle. It will make orbital rotation a bit easier.

Since the type of each parameter is thoroughly checked, the check for the parameter type raises an error when numpy array is used. Turns out np.int64 is not a child of int, so we cannot just check that the input is an integer. I've hacked in a fix (also check if the parameter is an instance of np.int64), but there should be a better way to handle these inputs.

This problem also occurs with float and np.float64.

Fix the path to test files in ${HORTONDATA}/test folder.

Like horton tests, if $HORTONDATA is set, then we can access the test files like this:

hf_dict = hartreefock(fn="test/h2.xyz", basis="6-31g**", nelec=nelec)

And it should work for all of us no matter where you are running your tests from.

All of the init's have some from ... import ... statements, which can all be removed I think. If they are removed, then modules must be imported by absolute name (e.g. import wfns.wavefunction.base_wavefunction ) and cannot be imported from the directory/supermodule (e.g. from wfns.wavefunction import base_wavefunction). Since explicit absolute imports are recommended anyways, we can remove all of the contents of init. There are no dependence on specific import order (as far as I know), so the necessary modules can be imported explicitly.

Right now, we are assuming that there are always pairs of electrons (even number). But the Slater determinants in projection space is expressed with respect to spin orbitals. For APIG derived things, we will, in most cases, work with spatial orbitals: all excitations will be either redundant (excite to same orbital) or pair excitations. So in these cases, we do not need to have spin orbitals in the Slater determinants.

However, if we have any odd electron system, the existing code breaks. We would need to build coefficient matrices that have approx. twice as many rows and select twice as many orbitals and so on. It would be necessary at one point to be able to do unrestricted systems.

So what do we do? Do we have a restricted projecedt sch solver and unrestricted sch solver, or a general sch solver and we manually double the number of pairs and spatial orbitals, or a general sch solver with a flag that specifies which of the two to use?

In the case of restricted projected sch solver, do we need to have slater determinants wrt spin orbitals? We will certainly need it sometimes (due to single excitation in over projections). It doesn't seem like we will lose much by having it expressed wrt spin orbitals, though.

At the moment, the cached functions _olp and _olp_deriv are stored in the _cache_fns attribute. This attribute is created in the __init__ of the class of wavefunctions that will be caching the overlaps (and its derivative). This approach seems to be quite ugly and I think can be cleaned up a bit.

Maybe it'll be better to give the option of caching the overlaps to all wavefunctions and let the user decide accordingly. This would suggest that the _cache_fns be moved to the load_cache method.