Fock matrix containers.
ebcc.ham.fock.RFock(cc, array=None, space=None, mo_coeff=None, g=None)
Bases: BaseFock
, BaseRHamiltonian
Restricted Fock matrix container class.
Initialise the Fock matrix.
Parameters: |
|
---|
Source code in ebcc/ham/base.py
def __init__(
self,
cc: BaseEBCC,
array: Optional[Any] = None,
space: Optional[tuple[Any, ...]] = None,
mo_coeff: Optional[tuple[Any, ...]] = None,
g: Optional[Namespace[Any]] = None,
) -> None:
"""Initialise the Fock matrix.
Args:
cc: Coupled cluster object.
array: Fock matrix in the MO basis.
space: Space object for each index.
mo_coeff: Molecular orbital coefficients for each index.
g: Namespace containing blocks of the electron-boson coupling matrix.
"""
Namespace.__init__(self)
# Parameters:
self.__dict__["cc"] = cc
self.__dict__["space"] = space if space is not None else (cc.space,) * 2
self.__dict__["mo_coeff"] = mo_coeff if mo_coeff is not None else (cc.mo_coeff,) * 2
self.__dict__["array"] = array if array is not None else self._get_fock()
self.__dict__["g"] = g if g is not None else cc.g
# Boson parameters:
self.__dict__["shift"] = cc.options.shift if g is not None else None
self.__dict__["xi"] = cc.xi if g is not None else None
ebcc.ham.fock.RFock.__getitem__(key)
Just-in-time getter.
Parameters: |
|
---|
Returns: |
|
---|
Source code in ebcc/ham/fock.py
def __getitem__(self, key: str) -> NDArray[T]:
"""Just-in-time getter.
Args:
key: Key to get.
Returns:
Fock matrix for the given spaces.
"""
if key not in self._members:
i = self.space[0].mask(key[0])
j = self.space[1].mask(key[1])
self._members[key] = self.array[i][:, j].copy()
if self.shift:
xi = self.xi
g = self.g.__getattr__(f"b{key}").copy()
g += self.g.__getattr__(f"b{key[::-1]}").transpose(0, 2, 1)
self._members[key] -= util.einsum("I,Ipq->pq", xi, g)
return self._members[key]
ebcc.ham.fock.UFock(cc, array=None, space=None, mo_coeff=None, g=None)
Bases: BaseFock
, BaseUHamiltonian
Unrestricted Fock matrix container class.
Initialise the Fock matrix.
Parameters: |
|
---|
Source code in ebcc/ham/base.py
def __init__(
self,
cc: BaseEBCC,
array: Optional[Any] = None,
space: Optional[tuple[Any, ...]] = None,
mo_coeff: Optional[tuple[Any, ...]] = None,
g: Optional[Namespace[Any]] = None,
) -> None:
"""Initialise the Fock matrix.
Args:
cc: Coupled cluster object.
array: Fock matrix in the MO basis.
space: Space object for each index.
mo_coeff: Molecular orbital coefficients for each index.
g: Namespace containing blocks of the electron-boson coupling matrix.
"""
Namespace.__init__(self)
# Parameters:
self.__dict__["cc"] = cc
self.__dict__["space"] = space if space is not None else (cc.space,) * 2
self.__dict__["mo_coeff"] = mo_coeff if mo_coeff is not None else (cc.mo_coeff,) * 2
self.__dict__["array"] = array if array is not None else self._get_fock()
self.__dict__["g"] = g if g is not None else cc.g
# Boson parameters:
self.__dict__["shift"] = cc.options.shift if g is not None else None
self.__dict__["xi"] = cc.xi if g is not None else None
ebcc.ham.fock.UFock.__getitem__(key)
Just-in-time getter.
Parameters: |
|
---|
Returns: |
|
---|
Source code in ebcc/ham/fock.py
def __getitem__(self, key: str) -> RFock:
"""Just-in-time getter.
Args:
key: Key to get.
Returns:
Fock matrix for the given spin.
"""
if key not in ("aa", "bb"):
raise KeyError(f"Invalid key: {key}")
if key not in self._members:
i = "ab".index(key[0])
self._members[key] = RFock(
self.cc,
array=self.array[i],
space=(self.space[0][i], self.space[1][i]),
mo_coeff=(self.mo_coeff[0][i], self.mo_coeff[1][i]),
g=self.g[key] if self.g is not None else None,
)
return self._members[key]
ebcc.ham.fock.GFock(cc, array=None, space=None, mo_coeff=None, g=None)
Bases: BaseFock
, BaseGHamiltonian
Generalised Fock matrix container class.
Initialise the Fock matrix.
Parameters: |
|
---|
Source code in ebcc/ham/base.py
def __init__(
self,
cc: BaseEBCC,
array: Optional[Any] = None,
space: Optional[tuple[Any, ...]] = None,
mo_coeff: Optional[tuple[Any, ...]] = None,
g: Optional[Namespace[Any]] = None,
) -> None:
"""Initialise the Fock matrix.
Args:
cc: Coupled cluster object.
array: Fock matrix in the MO basis.
space: Space object for each index.
mo_coeff: Molecular orbital coefficients for each index.
g: Namespace containing blocks of the electron-boson coupling matrix.
"""
Namespace.__init__(self)
# Parameters:
self.__dict__["cc"] = cc
self.__dict__["space"] = space if space is not None else (cc.space,) * 2
self.__dict__["mo_coeff"] = mo_coeff if mo_coeff is not None else (cc.mo_coeff,) * 2
self.__dict__["array"] = array if array is not None else self._get_fock()
self.__dict__["g"] = g if g is not None else cc.g
# Boson parameters:
self.__dict__["shift"] = cc.options.shift if g is not None else None
self.__dict__["xi"] = cc.xi if g is not None else None
ebcc.ham.fock.GFock.__getitem__(key)
Just-in-time getter.
Parameters: |
|
---|
Returns: |
|
---|
Source code in ebcc/ham/fock.py
def __getitem__(self, key: str) -> NDArray[T]:
"""Just-in-time getter.
Args:
key: Key to get.
Returns:
Fock matrix for the given spin.
"""
if key not in self._members:
i = self.space[0].mask(key[0])
j = self.space[1].mask(key[1])
self._members[key] = self.array[i][:, j].copy()
if self.shift:
xi = self.xi
g = self.g.__getattr__(f"b{key}").copy()
g += self.g.__getattr__(f"b{key[::-1]}").transpose(0, 2, 1)
self._members[key] -= util.einsum("I,Ipq->pq", xi, g)
return self._members[key]