support time dependent system

hiro_models/model_base.py

@@ -32,8 +32,8 @@ class Model(ABC):
     description: str = ""
     """A free-form description of the model instance."""
 
-    # t: NDArray[np.float64] = np.linspace(0, 10, 1000)
-    # """The simulation time points."""
+    t: NDArray[np.float64] = np.linspace(0, 10, 1000)
+    """The simulation time points."""
 
     __base_version__: int = 1
     """The version of the model base."""
@@ -146,7 +146,7 @@ class Model(ABC):
 
     @property
     @abstractmethod
-    def system(self) -> qt.Qobj:
+    def system(self) -> DynamicMatrix:
         """The system hamiltonian."""
 
         pass
@@ -230,7 +230,7 @@ class Model(ABC):
     ###########################################################################
 
     def system_expectation(
-        self, data: HIData, operator: qt.Qobj, **kwargs
+        self, data: HIData, operator: DynamicMatrix, **kwargs
     ) -> EnsembleValue:
         """Calculates the expectation value of ``operator`` from the
         hierarchy data ``data``.
@@ -246,7 +246,8 @@ class Model(ABC):
 
         return hopsflow.util.operator_expectation_ensemble(
             data.valid_sample_iterator(data.stoc_traj),  # type: ignore
-            operator.full(),
+            operator,
+            self.t,
             normalize=True,  # always nonlinear
             save=f"{operator_hash}_{self.hexhash}",
             N=N,

hiro_models/one_qubit_model.py

@@ -129,6 +129,13 @@ class QubitModel(Model):
     processes.
     """
 
+    H: DynamicMatrix = field(
+        default_factory=lambda: ConstantMatrix(1 / 2 * qt.sigmaz().full())
+    )  # type: ignore
+    """
+    The system hamiltonian :math:`H` with shape ``(2, 2)``.
+    """
+
     @property
     def coupling_operators(self) -> list[DynamicMatrix]:
         """The bath coupling operators :math:`L`."""
@@ -136,10 +143,10 @@ class QubitModel(Model):
         return [self.L]
 
     @property
-    def system(self) -> qt.Qobj:
+    def system(self) -> DynamicMatrix:
         """The system hamiltonian."""
 
-        return 1 / 2 * qt.sigmaz()  # type: ignore
+        return self.H  # type: ignore
 
     @property
     def bcf_norm(self) -> float:
@@ -318,7 +325,7 @@ class QubitModel(Model):
         g, w = self.bcf_coefficients(self.bcf_terms)
 
         system = params.SysP(
-            H_sys=self.system.full(),
+            H_sys=self.system,
             L=self.coupling_operators,
             g=g,
             w=w,