Feature suggestion: Support for adding custom subspaces

[obackhouse]

When attempting to form a Tensor which has two axes in different MoSpaces, an libtensor exception is raised if the axes are both either 'v1' or 'f' subspaces. This isn't exactly when both are in the same subspaces, because it seems to work fine if they are both 'o1'.

import numpy as np
from adcc import Tensor, ReferenceState
from adcc.functions import direct_sum

mol_a = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='sto3g', verbose=0)
rhf_a = scf.RHF(mol_a).run(conv_tol=1e-12)
rs_a = ReferenceState(rhf_a)

mol_b = gto.M(atom='O 0 0 0; H 0 0 1; H 0 1 0', basis='cc-pvdz', verbose=0)
rhf_b = scf.RHF(mol_b).run(conv_tol=1e-12)
rs_b = ReferenceState(rhf_b)

small_space = rs_a.mospaces
large_space = rs_b.mospaces

# f, o1 - works
try:
    c = direct_sum('s,l->sl', Tensor(small_space, 'f'), Tensor(large_space, 'o1'))
    c.set_from_ndarray(np.eye(small_space.n_orbs('f'), large_space.n_orbs('o1')))
except Exception as e:
    print(e)
    print('f, o1 failed')

# f, v1 - works
try:
    c = direct_sum('s,l->sl', Tensor(small_space, 'f'), Tensor(large_space, 'v1'))
    c.set_from_ndarray(np.eye(small_space.n_orbs('f'), large_space.n_orbs('v1')))
except Exception as e:
    print(e)
    print('f, v1 failed')

# o1, o1 - works
try:
    c = direct_sum('s,l->sl', Tensor(small_space, 'o1'), Tensor(large_space, 'o1'))
    c.set_from_ndarray(np.eye(small_space.n_orbs('o1'), large_space.n_orbs('o1')))
except Exception as e:
    print(e)
    print('o1, o1 failed')

# v1, v1 - fails
try:
    c = direct_sum('s,l->sl', Tensor(small_space, 'v1'), Tensor(large_space, 'v1'))
    c.set_from_ndarray(np.eye(small_space.n_orbs('v1'), large_space.n_orbs('v1')))
except Exception as e:
    print(e)
    print('v1, v1 failed')

# f, f - fails
try:
    c = direct_sum('s,l->sl', Tensor(small_space, 'f'), Tensor(large_space, 'f'))
    c.set_from_ndarray(np.eye(small_space.n_orbs('f'), large_space.n_orbs('f')))
except Exception as e:
    print(e)
    print('f, f failed')

[mfherbst]

Cheers for the bug report @obackhouse ! Yes that definitely is not intended behaviour.

Update: See below.

[obackhouse]

On second thought, this is probably going to happen for 'o1' as well, but since both of my spaces here have the same occupancy it didn't come up.

[mfherbst]

Ah wait ... I was reading to quickly. You are combining in one direct_sum tensor axes from different Reference states (and thus different mospaces). I would not have expected this to work for any combination of axes actually ... I'm surprised you get no segfault or anything nasty like that. I don't think libtensor has been written with that even in mind. With that said ... I doubt we will be able to support such kind of workflow. The only thing I can realistically see being supported is to effectively split up the v1 space into the orbital indices also present in sto-3g and the indices only in cc-pvtz. That does, however, not yet work and will be a bit of effort to get working.

For the foreseeable time, the only way I think this sort of stuff can be achieved is to zero-pad the sto-3g results and place them into the full cc-pvtz spaces using numpy arrays. This is not super efficient ... and of course it means that your sto-3g tensors blow up into tensors with quite a number of zero blocks.

Hope this helps.

[mfherbst]

Also: Can you give a bit of context what you actually want to do with this sort of workflow?

[obackhouse]

Understood - if this is just something that libtensor doesn't support then fair enough. The padding would indeed be a bit of a pain.

As for workflows where this would be needed, for me that is in AGF2 where one iterates a quasiparticle basis and requires to partially transform ERIs i.e. (pq|rs) -> (pi|jk), where pqrs are MOs and ijkl are a basis defined by a subsequent iteration. The iterated basis is not the same size as the MO space, hence why I just used two different basis set MO spaces in my example - it wouldn't be exactly this but would be functionally similar. This is of course only tangential to ADC(2) and likely won't be a problem for any standard ADC method.

Possibly a more common place this may be needed is in density fitting. In that case you have i.e. (Q|ij) which is a tensor in two different spaces (auxiliary, MO, MO). I suppose that this is different to my above case, since there isn't an occupancy in that space - does libtensor provide support for DF integrals? If so, how does it handle the space of the auxiliary basis?

[mfherbst]

partially transform ERIs i.e. (pq|rs) -> (pi|jk),

Yes indeed that is related to density-fitting, which is supported in libtensor, but we have not yet made effort to get that into adcc. I should say that this is on the medium-term list of things we would like to support. Whether there is occupancy or not is not quite the issue by the way. It's more that right now the "o1", "v1" etc spaces are hard-coded and for now it is not possible to "define your own subspace". This is not a limitation in libtensor, but just of the interface.

Actually support for defining custom subspaces is quite a bit simpler to achieve than what I originally thought you wanted. I first thought you wanted to split an existing orbital subspace into two. But your alternative suggestion (add a new subspace and than transform existing tensors to that subspace using some form of projection / unitary rotation) is much simpler and something I could well see supported.

I'll keep this in mind when I think about density-fitting.