qc2.data ======== .. py:module:: qc2.data .. autoapi-nested-parse:: qc2 DATA package. Submodules ---------- .. toctree:: :maxdepth: 1 /autoapi/qc2/data/data/index Classes ------- .. autoapisummary:: qc2.data.qc2Data Package Contents ---------------- .. py:class:: qc2Data(filename: str = 'qchem_data.hdf5', molecule: ase.Atoms = Atoms(), algorithm: qc2.algorithms.base.base_algorithm.BaseAlgorithm = BaseAlgorithm(), *, schema: str = 'qcschema') Main qc2 class. This class orchestrates classical qchem programs and python libraries for quantum computing. .. attribute:: _schema Format in which to save qchem data. Options are ``qcschema`` or ``fcidump``. Defaults to ``qcschema``. :type: str .. attribute:: filename The path to the HDF5 or fcidump file used to save/read qchem data. :type: str .. attribute:: molecule Attribute representing the molecular structure as an ASE :class:`ase.atoms.Atoms` instance. :type: Atoms .. attribute:: algorithm Instance of the algorithm to be run. Examples are :class:`~qc2.algorithm.qiskit.vqe.VQE` and :class:`~qc2.algorithm.pennylane.oo_vqe.oo_VQE`. :type: BaseAlgorithm .. py:attribute:: _schema :value: 'qcschema' .. py:attribute:: _filename :value: 'qchem_data.hdf5' .. py:attribute:: _molecule :value: None .. py:property:: molecule :type: ase.Atoms Returns the molecule attribute. :returns: Molecule as an ASE :class:`ase.atoms.Atoms` object. .. py:attribute:: _algorithm :value: None .. py:property:: algorithm :type: qc2.algorithms.base.base_algorithm.BaseAlgorithm Returns the chosen algorithm. :returns: Instance of an algorithm class, *e.g.*, VQE. .. py:method:: _check_filename_extension() -> None Ensures that files have proper extensions. .. py:method:: run() -> None Runs ASE qchem calculator and saves the data into a formated file. :returns: None **Example** >>> from ase.build import molecule >>> from qc2.ase import DIRAC >>> from qc2.data import qc2Data >>> >>> mol = molecule('H2') >>> >>> hdf5_file = 'h2.hdf5' >>> qc2data = qc2Data(hdf5_file, mol, schema='qcschema') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() >>> >>> fcidump_file = 'h2.fcidump' >>> qc2data = qc2Data(fcidump_file, mol, schema='fcidump') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() .. py:method:: read_schema() -> Union[qiskit_nature.second_q.formats.qcschema.QCSchema, qiskit_nature.second_q.formats.fcidump.FCIDump] Reads and stores data in :class:`QCSchema` or :class:`FCIDump`. Reads and stores the required data from an HDF5 or FCIDump file as either a :class:`QCSchema` or :class:`FCIDump` dataclass instance. :returns: Instance of :class:`QCSchema` or :class:`FCIDump` dataclass. :rtype: Union[QCSchema, FCIDump] .. rubric:: Notes See qiskit_nature/second_q/formats for more information on the supported data formats. **Example** >>> from ase.build import molecule >>> from qc2.ase import DIRAC >>> from qc2.data import qc2Data >>> >>> mol = molecule('H2') >>> >>> hdf5_file = 'h2.hdf5' >>> qc2data = qc2Data(hdf5_file, mol, schema='qcschema') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() >>> qcschema = qc2data.read_schema() >>> >>> fcidump_file = 'h2.fcidump' >>> qc2data = qc2Data(fcidump_file, mol, schema='fcidump') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() >>> fcidump = qc2data.read_schema() .. py:method:: process_schema(*, basis: str = 'molecular') -> qiskit_nature.second_q.problems.ElectronicStructureProblem Creates an instance of :class:`ElectronicStructureProblem`. Reads data using the :meth:`~.read_schema` method and converts it into an instance of :class:`ElectronicStructureProblem` based on the specified schema format (``fcidump`` or ``qcschema``) and electronic basis as defined by :class:`qiskit.ElectronicBasis`. :param basis: The basis in which to construct the :class:`ElectronicStructureProblem`. Options are ``atomic`` or ``molecular``. Defaults to ``molecular``. :type basis: str, optional :returns: An instance representing the :class:`ElectronicStructureProblem`. :rtype: ElectronicStructureProblem .. rubric:: Notes - For ``fcidump`` schema, the conversion is done using the `fcidump_to_problem` function from qiskit_nature/second_q/formats/fcidump_translator.py. - For ``qcschema`` schema, the conversion is done using the `qcschema_to_problem` function from qiskit_nature/second_q/formats/qcschema_translator.py. - Dipoles are excluded when converting `QCSchema` data. **Example** >>> from ase.build import molecule >>> from qc2.ase import DIRAC >>> from qc2.data import qc2Data >>> >>> mol = molecule('H2') >>> >>> hdf5_file = 'h2.hdf5' >>> qc2data = qc2Data(hdf5_file, mol, schema='qcschema') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() >>> es_problem = qc2data.process_schema(basis='atomic') >>> >>> fcidump_file = 'h2.fcidump' >>> qc2data = qc2Data(fcidump_file, mol, schema='fcidump') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() >>> es_problem = qc2data.process_schema() .. py:method:: get_active_space_hamiltonian(num_electrons: Union[int, Tuple[int, int]], num_spatial_orbitals: int, *, initial_es_problem: Optional[qiskit_nature.second_q.problems.ElectronicStructureProblem] = None) -> Tuple[qiskit_nature.second_q.problems.ElectronicStructureProblem, float, qiskit_nature.second_q.hamiltonians.ElectronicEnergy] Builds the active-space reduced Hamiltonian. :param num_electrons: The number of active electrons. If a tuple is provided, it represents alpha and beta active electrons. :type num_electrons: Union[int, Tuple[int, int]] :param num_spatial_orbitals: The number of spatial orbitals. :type num_spatial_orbitals: int :param initial_es_problem: Initial instance of :class:`ElectronicStructureProblem`. If None, it is instantiated internally. Defaults to None. :type initial_es_problem: Optional[ElectronicStructureProblem] :returns: - active_space_es_problem (ElectronicStructureProblem): final active space transformed :class:`ElectronicStructureProblem`. - core_energy (float): The core energy, which includes the nuclear repulsion energy and the energy of inactive orbitals. - active_space_hamiltonian (ElectronicEnergy): Instance of :class:`ElectronicEnergy`, the active-space reduced Hamiltonian. :rtype: Tuple[ElectronicStructureProblem, float, ElectronicEnergy] .. rubric:: Notes - The active-space reduced Hamiltonian is obtained by transforming the original electronic structure problem's Hamiltonian using an ActiveSpaceTransformer. - The core energy is computed as the sum of the nuclear repulsion energy and the energy of inactive orbitals. **Example** >>> from ase.build import molecule >>> from qc2.ase import DIRAC >>> from qc2.data import qc2Data >>> >>> mol = molecule('H2') >>> hdf5_file = 'h2.hdf5' >>> qc2data = qc2Data(hdf5_file, mol, schema='qcschema') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() >>> n_electrons = (1, 1) >>> n_spatial_orbitals = 2 >>> (es_problem, e_core, ham) = qc2data.get_active_space_hamiltonian( ... n_electrons, n_spatial_orbitals ... ) .. py:method:: get_transformed_hamiltonian(*, initial_es_problem: qiskit_nature.second_q.problems.ElectronicStructureProblem, matrix_transform_a: numpy.ndarray, matrix_transform_b: Optional[numpy.ndarray] = None, initial_basis: str = 'atomic', final_basis: str = 'molecular') -> Tuple[qiskit_nature.second_q.problems.ElectronicStructureProblem, qiskit_nature.second_q.hamiltonians.ElectronicEnergy] Transforms the Hamiltonian from one basis set to another. :param initial_es_problem: The original electronic structure problem. :type initial_es_problem: ElectronicStructureProblem :param matrix_transform_a: The transformation matrix for alpha spin orbitals. :type matrix_transform_a: np.ndarray :param matrix_transform_b: The transformation matrix for beta spin orbitals. :type matrix_transform_b: np.ndarray, optional :param initial_basis: The initial basis set. Defaults to ``atomic``. :type initial_basis: str, optional :param final_basis: The final basis set to transform to. Defaults to ``molecular``. :type final_basis: str, optional :returns: - transformed_es_problem (ElectronicStructureProblem): An instance of the transformed :class:`ElectronicStructureProblem`. - transformed_hamiltonian (ElectronicEnergy): An instance of :class:`ElectronicEnergy`, the transformed Hamiltonian. :rtype: Tuple[ElectronicStructureProblem, ElectronicEnergy] **Example** >>> from ase.build import molecule >>> from qc2.ase import PySCF >>> from qc2.data import qc2Data >>> >>> mol = molecule('H2') >>> hdf5_file = 'h2.hdf5' >>> qc2data = qc2Data(hdf5_file, mol, schema='qcschema') >>> qc2data.molecule.calc = PySCF(...) # => specify qchem calculator >>> qc2data.run() >>> ao_es_problem = qc2data.process_schema(basis='atomic') >>> mo_coeff_a = np.array( ... [[0.54884228, 1.21245192], ... [0.54884228, -1.21245192]] ... ) >>> mo_es_problem, hamiltonian = qc2data.get_transformed_hamiltonian( ... initial_es_problem=ao_es_problem, ... matrix_transform_a=mo_coeff_a, ... initial_basis="atomic", ... final_basis="molecular" ... ) .. py:method:: get_fermionic_hamiltonian(num_electrons: Union[int, Tuple[int, int]], num_spatial_orbitals: int, *, transform: bool = False, initial_es_problem: Optional[qiskit_nature.second_q.problems.ElectronicStructureProblem] = None, matrix_transform_a: Optional[numpy.ndarray] = None, matrix_transform_b: Optional[numpy.ndarray] = None, initial_basis: str = 'atomic', final_basis: str = 'molecular') -> Tuple[qiskit_nature.second_q.problems.ElectronicStructureProblem, float, qiskit_nature.second_q.operators.FermionicOp] Builds the fermionic Hamiltonian of a target molecule. This method constructs the electronic Hamiltonian in 2nd-quantization based on the provided parameters. It can optionally perform a basis set transformation if the `transform` flag is set to True. :param num_electrons: The number of active electrons. If this is a tuple, it represents the number of alpha- and beta-spin electrons, respectively. If this is a number, it is interpreted as the total number of active electrons, should be even, and implies that the number of alpha and beta electrons equals half of this value, respectively. :type num_electrons: Union[int, Tuple[int, int]] :param num_spatial_orbitals: The number of active orbitals. :type num_spatial_orbitals: int :param transform: If True, performs a basis transformation. Defaults to False. :type transform: bool, optional :param initial_es_problem: The initial electronic structure problem. Required if `transform` is True. Defaults to None. :type initial_es_problem: ElectronicStructureProblem, optional :param matrix_transform_a: Transformation matrix for alpha spin orbitals. Required if `transform` is True. Defaults to None. :type matrix_transform_a: np.ndarray, optional :param matrix_transform_b: Transformation matrix for beta spin orbitals. Required if `transform` is True. Defaults to None. :type matrix_transform_b: np.ndarray, optional :param initial_basis: The initial basis set. Defaults to "atomic". :type initial_basis: str, optional :param final_basis: The final basis set to transform to. Defaults to "molecular". :type final_basis: str, optional :returns: - core_energy (float): The core energy after active space transformation. - es_problem (ElectronicStructureProblem): An instance of the :class:`ElectronicStructureProblem`. - second_q_op (FermionicOp): An instance of :class:`FermionicOp` representing the fermionic Hamiltonian in 2nd quantization. :rtype: Tuple[float, ElectronicStructureProblem, FermionicOp] :raises ValueError: If :attr:`num_electrons` or :attr:`num_spatial_orbitals` is None. Or If :attr:`initial_es_problem` is None for :attr:`transform` equal True. .. rubric:: Notes Based on the qiskit-nature modules: qiskit_nature/second_q/problems/electronic_structure_problem.py qiskit_nature/second_q/transformers/active_space_transformer.py **Example** >>> from ase.build import molecule >>> from qc2.ase import DIRAC >>> from qc2.data import qc2Data >>> >>> mol = molecule('H2') >>> hdf5_file = 'h2.hdf5' >>> qc2data = qc2Data(hdf5_file, mol, schema='qcschema') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() >>> n_electrons = (1, 1) >>> n_spatial_orbitals = 2 >>> (es_prob, e_core, op) = qc2data.get_fermionic_hamiltonian( ... n_electrons, n_spatial_orbitals ... ) .. py:method:: get_qubit_hamiltonian(num_electrons: Union[int, Tuple[int, int]], num_spatial_orbitals: int, mapper: qiskit_nature.second_q.mappers.QubitMapper = JordanWignerMapper(), *, format: str = 'qiskit', transform: bool = False, initial_es_problem: Optional[qiskit_nature.second_q.problems.ElectronicStructureProblem] = None, matrix_transform_a: Optional[numpy.ndarray] = None, matrix_transform_b: Optional[numpy.ndarray] = None, initial_basis: str = 'atomic', final_basis: str = 'molecular') -> Tuple[float, Union[qiskit.quantum_info.SparsePauliOp, PennyLaneOperatorType]] Generates the qubit Hamiltonian of a target molecule. This method generates the qubit Hamiltonian representation of a target molecule. It can optionally perform a basis set transformation if the ``transform`` flag is True. :param num_electrons: The number of active electrons. If this is a tuple, it represents the number of alpha- and beta-spin electrons, respectively. If this is a number, it is interpreted as the total number of active electrons, should be even, and implies that the number of alpha and beta electrons equals half of this value, respectively. :type num_electrons: Union[int, Tuple[int, int]] :param num_spatial_orbitals: The number of active orbitals. :type num_spatial_orbitals: int :param mapper: The qubit mapping strategy to convert fermionic operators to qubit operators. Defaults to ``JordanWignerMapper()``. :type mapper: QubitMapper, optional :param format: The format in which to return the qubit Hamiltonian. Supported formats are ``qiskit`` and ``pennylane``. Defaults to ``qiskit``. :type format: str, optional :param transform: If True, performs a basis transformation. Defaults to False. :type transform: bool, optional :param initial_es_problem: The initial electronic structure problem. Required if `transform` is True. Defaults to None. :type initial_es_problem: ElectronicStructureProblem, optional :param matrix_transform_a: Transformation matrix for alpha spin orbitals. Required if `transform` is True. Defaults to None. :type matrix_transform_a: np.ndarray, optional :param matrix_transform_b: Transformation matrix for beta spin orbitals. Required if `transform` is True. Defaults to None. :type matrix_transform_b: np.ndarray, optional :param initial_basis: The initial basis set. Defaults to ``atomic``. :type initial_basis: str, optional :param final_basis: The final basis set to transform to. Defaults to ``molecular``. :type final_basis: str, optional :returns: - core_energy (float): The core energy after active space transformation. - qubit_op (Union[SparsePauliOp, Operator]): If the format is ``qiskit``, it returns a :class:`SparsePauliOp` representing the qubit Hamiltonian in the qiskit format. If the format is ``pennylane``, it returns a :class:`Operator` instance representing the qubit Hamiltonian in the PennyLane format. :rtype: Tuple[float, Union[SparsePauliOp, Operator]] :raises TypeError: If the provided `format` is not supported (not ``qiskit`` nor ``pennylane``). **Example** >>> from ase.build import molecule >>> from qc2.ase import DIRAC >>> from qc2.data import qc2Data >>> >>> mol = molecule('H2') >>> hdf5_file = 'h2.hdf5' >>> qc2data = qc2Data(hdf5_file, mol, schema='qcschema') >>> qc2data.molecule.calc = DIRAC(...) # => specify qchem calculator >>> qc2data.run() >>> n_electrons = (1, 1) >>> n_spatial_orbitals = 2 >>> mapper = BravyiKitaevMapper() >>> (e_core, qubit_op) = qc2data.get_qubit_hamiltonian( ... n_electrons, n_spatial_orbitals, mapper, format='qiskit' ... )