Thanks @tom I think you’re right this solution will just do it.
Here is what it looks like by the way:
def p_mean_pytorch(tensor,power):
          mean_tensor = tensor.pow(power).mean(0)
          return (  (mean_tensor * mean_tensor.sign() ).pow(1/power) * mean_pow.sign() )
     

i.e.
p_mean_pytorch(torch.tensor([[-100,33,99],[39,9,-10000],[1,3,4],[0,0,0]]).float(), 3)
tensor([  -61.7249,    20.9335, -6299.6050])

I would have liked to get the same root resolution as their implementation of the paper but I guess this will do  . I will check the precomputing part later, the only challenging part is telling torch.ger to be conditional but I think I can do something with a mask.
 def power_mean_precompute_3(tensor):
     magical_number=torch.tensor([0.5,0.866]) # np.power(-1+0j,1/3)
     mean_tensor = tensor.pow(3).mean(0)
     mask = torch.le(mean_tensor,0.)
     le_result = torch.ger((mean_tensor * mask.float() * mean_tensor.sign()).pow(1/3),magical_number)
     ge_result = torch.ger((mean_tensor * (~mask).float()).pow(1/3),torch.tensor([1.,0]))
     return le_result+ge_result

power_mean_precompute_3(torch.tensor([[-100,33,99],[39,9,-10000],[1,3,4],[0,0,0]]).float())
tensor([[  30.8625,   53.4538],
        [  20.9335,    0.0000],
        [3149.8025, 5455.4580]])

gen_mean(torch.tensor([[-100,33,99],[39,9,-10000],[1,3,4],[0,0,0]]).numpy(),3)
array([  30.86246739  +53.45536157j,   20.93346287   +0.        j,
       3149.80160592+5455.61641521j])

It seems to be doing what intended but it should get a bit optimized.
Also, this won’t work for a batch as torch.ger won’t allow it but to adapt it shouldn’t be too hard.