feat: general phsp

For chain decay

A -> (R -> D E) B C

generate the phase space as

generate_phsp(
    ma, 
    [ 
        [mr, [md, me]],
         mb,
         mc
    ],
    N = 1000
)