code: def test_gamma_gaussian_tensordot(
      dot_dims,
      x_batch_shape,
      x_dim,
      x_rank,
      y_batch_shape,
      y_dim,
      y_rank,
  ):
      x_rank = min(x_rank, x_dim)
      y_rank = min(y_rank, y_dim)
      x = random_gamma_gaussian(x_batch_shape, x_dim, x_rank)
      y = random_gamma_gaussian(y_batch_shape, y_dim, y_rank)
      na = x_dim - dot_dims
      nb = dot_dims
      nc = y_dim - dot_dims
      try:
          torch.linalg.cholesky(
              x.precision[..., na:, na:] + y.precision[..., :nb, :nb]
          )
      except RuntimeError:
          pytest.skip(
              "Cannot marginalize the common variables of two Gaussians."
          )
  
      z = gamma_gaussian_tensordot(x, y, dot_dims)
      assert z.dim() == x_dim + y_dim - 2 * dot_dims
  
      # We make these precision matrices positive definite to test the math
      x.precision = x.precision + 3 * torch.eye(x.dim())
      y.precision = y.precision + 3 * torch.eye(y.dim())
      z = gamma_gaussian_tensordot(x, y, dot_dims)
      # compare against broadcasting, adding, and marginalizing
      precision = pad(x.precision, (0, nc, 0, nc)) + pad(
          y.precision, (na, 0, na, 0)
      )
      info_vec = pad(x.info_vec, (0, nc)) + pad(y.info_vec, (na, 0))
      covariance = torch.inverse(precision)
      loc = (
          covariance.matmul(info_vec.unsqueeze(-1)).squeeze(-1)
          if info_vec.size(-1) > 0
          else info_vec
      )
      z_covariance = torch.inverse(z.precision)
      z_loc = z_covariance.matmul(
          z.info_vec.view(z.info_vec.shape + (int(z.dim() > 0),))
      ).sum(-1)
      assert_close(loc[..., :na], z_loc[..., :na])
      assert_close(loc[..., x_dim:], z_loc[..., na:])
      assert_close(
          covariance[..., :na, :na], z_covariance[..., :na, :na]
      )
      assert_close(
          covariance[..., :na, x_dim:], z_covariance[..., :na, na:]
      )
      assert_close(
          covariance[..., x_dim:, :na], z_covariance[..., na:, :na]
      )
      assert_close(
          covariance[..., x_dim:, x_dim:], z_covariance[..., na:, na:]
      )
  
      s = torch.randn(z.batch_shape).exp()
      # Assume a = c = 0, integrate out b
      num_samples = 200000
      scale = 10
      # generate samples in [-10, 10]
      value_b = (
          torch.rand((num_samples,) + z.batch_shape + (nb,)) * scale
          - scale / 2
      )
      value_x = pad(value_b, (na, 0))
      value_y = pad(value_b, (0, nc))
      expect = torch.logsumexp(
          x.log_density(value_x, s) + y.log_density(value_y, s), dim=0
      )
      expect += math.log(scale**nb / num_samples)
      actual = z.log_density(torch.zeros(z.batch_shape + (z.dim(),)), s)
      assert_close(
          actual.clamp(max=10.0),
          expect.clamp(max=10.0),
          atol=0.1,
          rtol=0.1,
      )
  
imports: import math
  import pytest
  import torch
  from torch.nn.functional import pad
  import pyro.distributions
  from pyro.distributions.util import broadcast_shape
  from pyro.ops.gamma_gaussian import GammaGaussian, gamma_and_mvn_to_gamma_gaussian, gamma_gaussian_tensordot, matrix_and_mvn_to_gamma_gaussian
  from tests.common import assert_close
  from tests.ops.gamma_gaussian import assert_close_gamma_gaussian, random_gamma, random_gamma_gaussian
  from tests.ops.gaussian import random_mvn
label_desc: Use image node_id 13 for calling a global function with example usage: test_gamma_gaussian_tensordot(dot_dims, x_batch_shape, x_dim, x_rank, y_batch_shape, y_dim, y_rank) without return types