You could add a dummy dimension in A in dim2:
a = torch.randn(16, 16, 1, 13)
b = torch.randn(16, 16, 13, 49)

res = torch.matmul(a, b).squeeze()

results = torch.zeros(16, 16, 49)
for in_i in range(16):
    for in_j in range(16):
        results[in_i,in_j,:] = torch.matmul(a[in_i, in_j , :], b[in_i, in_j, :, :])

print(torch.allclose(res, results))
> True