Implement optimized GEMM-like ops for Tensor and Gaussian

[fritzo]

trafficstars

• [x] #281 For Tensor use a GEMM op, fixing the bugs in https://github.com/pyro-ppl/funsor/blob/d18c78e/funsor/torch.py#L495 (which only works for scalar event shape)
• [ ] For Gaussian see pyro.ops.gaussian.gaussian_tensordot() pyro-ppl/pyro#1980 . This is lower priority than Tensor because the naive (x+y).reduce() version is only slightly more expensive.

[fritzo]

@eb8680 can you confirm that I can implement these as something like

@eager.register(Contraction, AddOp, MulOp, frozenset, Variadic[Tensor])
def eager_contract_tensor(red_op, bin_op, reduced_vars, *operands):
    equation = "TODO define an einsum string"
    data = pyro.ops.einsum.contract(equation, *(x.data for x in operands))
    inputs = "TODO"
    output = "TODO"
    return Tensor(data, inputs, output)

@eager.register(Contraction, LogaddexpOp, AddOp, frozenset, Variadic[Tensor])
def eager_contract_tensor(red_op, bin_op, reduced_vars, *operands):
    ...backend='pyro.ops.einsum.torch_log'...

[eb8680]

can you confirm

Yep, this should work as expected.

[fritzo]

@fehiepsi did you already complete this?