class_docstr: An algorithm for approximating Gaussian probabilities `P(X \in bounds)`, where
  `X ~ N(0, covariance_matrix)`.
class_docstr_tok: ["An","algorithm","for","approximating","Gaussian","probabilities","`","P","(","X","\\in","bounds",")","`",",","where","`","X","~","N","(","0",",","covariance_matrix",")","`","."]
docstr: Runs the MVNXPB solver instance for a fixed number of steps.
  
  Calculates a bivariate conditional approximation to P(X \in bounds), where
  X ~ N(0, Σ). For details, see [Genz2016numerical] or [Trinh2015bivariate]_.
docstr_tok: ["Runs","the","MVNXPB","solver","instance","for","a","fixed","number","of","steps",".","Calculates","a","bivariate","conditional","approximation","to","P","(","X","\\in","bounds",")",",","where","X","~","N","(","0",",","\u03a3",")",".","For","details",",","see","[","Genz2016numerical","]","or","[","Trinh2015bivariate","]","_","."]
code: def solve(
      self, num_steps: Optional[int] = None, eps: float = 1e-10
  ) -> Tensor:
      r"""Runs the MVNXPB solver instance for a fixed number of steps.
  
      Calculates a bivariate conditional approximation to P(X \in bounds), where
      X ~ N(0, Σ). For details, see [Genz2016numerical] or [Trinh2015bivariate]_.
      """
      if self.step > self.piv_chol.step:
          raise ValueError(
              "Invalid state: solver ran ahead of matrix decomposition."
          )
  
      # Unpack some terms
      start = self.step
      bounds = self.bounds
      piv_chol = self.piv_chol
      L = piv_chol.tril
      y = self.plug_ins
  
      # Subtract marginal log probability of final term from previous result if
      # it did not fit in a block.
      ndim = y.shape[-1]
      if ndim > start and start % 2:
          self.log_prob = self.log_prob - self.log_prob_extra
          self.log_prob_extra = None
  
      # Iteratively compute bivariate conditional approximation
      zero = get_constants_like(
          0, L
      )  # needed when calling `torch.where` below
      num_steps = num_steps or ndim - start
      for i in range(start, start + num_steps):
          should_update_chol = self.step == piv_chol.step
  
          # Determine next pivot element
          if should_update_chol:
              pivot = self.select_pivot()
          else:  # pivot using order specified by precomputed pivoted Cholesky step
              mask = self.perm[..., i:] == piv_chol.perm[..., i : i + 1]
              pivot = i + torch.nonzero(mask, as_tuple=True)[-1]
  
          if pivot is not None and torch.any(pivot > i):
              self.pivot_(pivot=pivot)
  
          # Compute whitened bounds conditional on preceding plug-ins
          Lii = L[..., i, i].clone()
          if should_update_chol:
              Lii = Lii.clip(min=0).sqrt()  # conditional stddev
          inv_Lii = Lii.reciprocal()
          bounds_i = bounds[..., i, :].clone()
          if i != 0:
              bounds_i = bounds_i - torch.sum(
                  L[..., i, :i].clone() * y[..., :i].clone(),
                  dim=-1,
                  keepdim=True,
              )
          lb, ub = (inv_Lii.unsqueeze(-1) * bounds_i).unbind(dim=-1)
  
          # Initialize `i`-th plug-in value as univariate conditional expectation
          Phi_i = Phi(ub) - Phi(lb)
          small = Phi_i <= i * eps
          y[..., i] = case_dispatcher(  # used to select next pivot
              out=(phi(lb) - phi(ub)) / Phi_i,
              cases=(  # fallback cases for enhanced numerical stability
                  (lambda: small & (lb < -9), lambda m: ub[m]),
                  (lambda: small & (lb > 9), lambda m: lb[m]),
                  (lambda: small, lambda m: 0.5 * (lb[m] + ub[m])),
              ),
          )
  
          # Maybe finalize the current block
          if i and i % 2:
              h = i - 1
              blk = slice(h, i + 1)
              Lhh = L[..., h, h].clone()
              Lih = L[..., i, h].clone()
  
              std_i = (Lii.square() + Lih.square()).sqrt()
              istds = 1 / torch.stack([Lhh, std_i], -1)
              blk_bounds = bounds[..., blk, :].clone()
              if i > 1:
                  blk_bounds = blk_bounds - (
                      L[..., blk, : i - 1].clone()
                      @ y[..., : i - 1, None].clone()
                  )
  
              blk_lower, blk_upper = (
                  pair.unbind(-1)  # pair of bounds for `yh` and `yi`
                  for pair in safe_mul(
                      istds.unsqueeze(-1), blk_bounds
                  ).unbind(-1)
              )
              blk_corr = Lhh * Lih * istds.prod(-1)
              blk_prob = bvn(blk_corr, *blk_lower, *blk_upper)
              zh, zi = bvnmom(
                  blk_corr, *blk_lower, *blk_upper, p=blk_prob
              )
  
              # Replace 1D expectations with 2D ones `L[blk, blk]^{-1} y[..., blk]`
              mask = blk_prob > zero
              y[..., h] = torch.where(mask, zh, zero)
              y[..., i] = torch.where(
                  mask, inv_Lii * (std_i * zi - Lih * zh), zero
              )
  
              # Update running approximation to log probability
              self.log_prob = self.log_prob + safe_log(blk_prob)
  
          self.step += 1
          if should_update_chol:
              piv_chol.update_(eps=eps)
  
      # Factor in univariate probability if final term fell outside of a block.
      if self.step % 2:
          self.log_prob_extra = safe_log(Phi_i)
          self.log_prob = self.log_prob + self.log_prob_extra
  
      return self.log_prob
  
imports: from __future__ import annotations
  from typing import Any, Optional, TypedDict
  from warnings import warn
  import torch
  from botorch.utils.probability.bvn import bvn, bvnmom
  from botorch.utils.probability.linalg import augment_cholesky, block_matrix_concat, PivotedCholesky
  from botorch.utils.probability.utils import case_dispatcher, get_constants_like, ndtr, phi, STANDARDIZED_RANGE, swap_along_dim_
  from botorch.utils.safe_math import log, mul
  from linear_operator.utils.cholesky import psd_safe_cholesky
  from linear_operator.utils.errors import NotPSDError
  from torch import LongTensor, Tensor
  from torch.nn.functional import pad
label_desc: Use image node_id 3 for calling the MVNXPB obj's underlying member method code with example usage: obj.solve(num_steps, eps) and returns: self