docstr: Compute the expected (feasible) hypervolume improvement given MC samples.
  
  Args:
      samples: A `sample_shape x batch_shape x q' x m`-dim tensor of samples.
      X: A `batch_shape x q x d`-dim tensor of inputs.
  
  Returns:
      A `batch_shape x (model_batch_shape)`-dim tensor of expected hypervolume
      improvement for each batch.
docstr_tok: ["Compute","the","expected","(","feasible",")","hypervolume","improvement","given","MC","samples",".","Args",":","samples",":","A","`","sample_shape","x","batch_shape","x","q","'","x","m","`","-dim","tensor","of","samples",".","X",":","A","`","batch_shape","x","q","x","d","`","-dim","tensor","of","inputs",".","Returns",":","A","`","batch_shape","x","(","model_batch_shape",")","`","-dim","tensor","of","expected","hypervolume","improvement","for","each","batch","."]
code: def _compute_log_qehvi(
      self, samples: Tensor, X: Optional[Tensor] = None
  ) -> Tensor:
      r"""Compute the expected (feasible) hypervolume improvement given MC samples.
  
      Args:
          samples: A `sample_shape x batch_shape x q' x m`-dim tensor of samples.
          X: A `batch_shape x q x d`-dim tensor of inputs.
  
      Returns:
          A `batch_shape x (model_batch_shape)`-dim tensor of expected hypervolume
          improvement for each batch.
      """
      # Note that the objective may subset the outcomes (e.g. this will usually happen
      # if there are constraints present).
      obj = self.objective(
          samples, X=X
      )  # mc_samples x batch_shape x q x m
      q = obj.shape[-2]
      if self.constraints is not None:
          log_feas_weights = compute_smoothed_feasibility_indicator(
              constraints=self.constraints,
              samples=samples,
              eta=self.eta,
              log=True,
              fat=self.fat,
          )
      device = self.ref_point.device
      q_subset_indices = self.compute_q_subset_indices(
          q_out=q, device=device
      )
      batch_shape = obj.shape[:-2]  # mc_samples x batch_shape
      # areas tensor is `mc_samples x batch_shape x num_cells x 2`-dim
      log_areas_per_segment = torch.full(
          size=(
              *batch_shape,
              self.cell_lower_bounds.shape[-2],  # num_cells
              2,  # for even and odd terms
          ),
          fill_value=-torch.inf,
          dtype=obj.dtype,
          device=device,
      )
  
      cell_batch_ndim = self.cell_lower_bounds.ndim - 2
      # conditionally adding mc_samples dim if cell_batch_ndim > 0
      # adding ones to shape equal in number to to batch_shape_ndim - cell_batch_ndim
      # adding cell_bounds batch shape w/o 1st dimension
      sample_batch_view_shape = torch.Size(
          [
              batch_shape[0] if cell_batch_ndim > 0 else 1,
              *[
                  1
                  for _ in range(
                      len(batch_shape) - max(cell_batch_ndim, 1)
                  )
              ],
              *self.cell_lower_bounds.shape[1:-2],
          ]
      )
      view_shape = (
          *sample_batch_view_shape,
          self.cell_upper_bounds.shape[-2],  # num_cells
          1,  # adding for q_choose_i dimension
          self.cell_upper_bounds.shape[-1],  # num_objectives
      )
  
      for i in range(1, self.q_out + 1):
          # TODO: we could use batches to compute (q choose i) and (q choose q-i)
          # simultaneously since subsets of size i and q-i have the same number of
          # elements. This would decrease the number of iterations, but increase
          # memory usage.
          q_choose_i = q_subset_indices[
              f"q_choose_{i}"
          ]  # q_choose_i x i
          # this tensor is mc_samples x batch_shape x i x q_choose_i x m
          obj_subsets = obj.index_select(
              dim=-2, index=q_choose_i.view(-1)
          )
          obj_subsets = obj_subsets.view(
              obj.shape[:-2] + q_choose_i.shape + obj.shape[-1:]
          )  # mc_samples x batch_shape x q_choose_i x i x m
  
          # NOTE: the order of operations in non-log _compute_qehvi is 3), 1), 2).
          # since 3) moved above 1), _log_improvement adds another Tensor dimension
          # that keeps track of num_cells.
  
          # 1) computes log smoothed improvement over the cell lower bounds.
          # mc_samples x batch_shape x num_cells x q_choose_i x i x m
          log_improvement_i = self._log_improvement(
              obj_subsets, view_shape
          )
  
          # 2) take the minimum log improvement over all i subsets.
          # since all hyperrectangles share one vertex, the opposite vertex of the
          # overlap is given by the component-wise minimum.
          # negative of maximum of negative log_improvement is approximation to min.
          log_improvement_i = self._smooth_min(
              log_improvement_i,
              dim=-2,
          )  # mc_samples x batch_shape x num_cells x q_choose_i x m
  
          # 3) compute the log lengths of the cells' sides.
          # mc_samples x batch_shape x num_cells x q_choose_i x m
          log_lengths_i = self._log_cell_lengths(
              log_improvement_i, view_shape
          )
  
          # 4) take product over hyperrectangle side lengths to compute area (m-dim).
          # after, log_areas_i is mc_samples x batch_shape x num_cells x q_choose_i
          log_areas_i = log_lengths_i.sum(
              dim=-1
          )  # areas_i = lengths_i.prod(dim=-1)
  
          # 5) if constraints are present, apply a differentiable approximation of
          # the indicator function.
          if self.constraints is not None:
              log_feas_subsets = log_feas_weights.index_select(
                  dim=-1, index=q_choose_i.view(-1)
              ).view(log_feas_weights.shape[:-1] + q_choose_i.shape)
              log_areas_i = log_areas_i + log_feas_subsets.unsqueeze(
                  -3
              ).sum(dim=-1)
  
          # 6) sum over all subsets of size i, i.e. reduce over q_choose_i-dim
          # after, log_areas_i is mc_samples x batch_shape x num_cells
          log_areas_i = logsumexp(
              log_areas_i, dim=-1
          )  # areas_i.sum(dim=-1)
  
          # 7) Using the inclusion-exclusion principle, set the sign to be positive
          # for subsets of odd sizes and negative for subsets of even size
          # in non-log space: areas_per_segment += (-1) ** (i + 1) * areas_i,
          # but here in log space, we need to keep track of sign:
          log_areas_per_segment[..., i % 2] = logplusexp(
              log_areas_per_segment[..., i % 2],
              log_areas_i,
          )
  
      # 8) subtract even from odd log area terms
      log_areas_per_segment = logdiffexp(
          log_a=log_areas_per_segment[..., 0],
          log_b=log_areas_per_segment[..., 1],
      )
  
      # 9) sum over segments (n_cells-dim) and average over MC samples
      return logmeanexp(logsumexp(log_areas_per_segment, dim=-1), dim=0)
  
imports: from __future__ import annotations
  from typing import Callable, List, Optional, Tuple, Union
  import torch
  from botorch.acquisition.logei import TAU_MAX, TAU_RELU
  from botorch.acquisition.multi_objective import MultiObjectiveMCAcquisitionFunction
  from botorch.acquisition.multi_objective.objective import MCMultiOutputObjective
  from botorch.models.model import Model
  from botorch.sampling.base import MCSampler
  from botorch.utils.multi_objective.box_decompositions.non_dominated import NondominatedPartitioning
  from botorch.utils.multi_objective.hypervolume import NoisyExpectedHypervolumeMixin, SubsetIndexCachingMixin
  from botorch.utils.objective import compute_smoothed_feasibility_indicator
  from botorch.utils.safe_math import fatmin, log_fatplus, log_softplus, logdiffexp, logmeanexp, logplusexp, logsumexp, smooth_amin
  from botorch.utils.transforms import concatenate_pending_points, is_ensemble, match_batch_shape, t_batch_mode_transform
  from torch import Tensor
label_desc: Use image node_id 2 for calling the qLogExpectedHypervolumeImprovement obj's underlying member method code with example usage: obj._compute_log_qehvi(samples, X) and returns: logmeanexp