Usually one would consult a paper or book to find integration points and weights for unit triangle and tetrahedra. I am  looking for a method to automatically compute such points and weights. The following Mathematica code example computes integration weights and points for unit line (quad/hexahedron) elements:unitGaussianQuadraturePoints[points_] :=   Sort[x /.     Solve[Evaluate[LegendreP[points, x] == 0], {x}], !      OrderedQ[N[{#1, #2}]] &];unitGaussianQuadratureWeights[points_] :=   Module[{gps, f, int, integr, vars, eqns},    gps = unitGaussianQuadraturePoints[points];   f[0, 0] := 1;   f[0., 0] := 1.;   f[x_, n_] := x^n;   int = Integrate[f[x, #], x] & /@ Range[0, points - 1];   integr = Subtract @@@ (int /. x :> {1, -1});   vars = Table[Unique[c], {Length[gps]}];   eqns =     Table[Plus @@ Thread[Times[vars, f[#, i - 1] & /@ gps]] ==       integr[[i]], {i, points}];   Return[(vars /. Solve[eqns, vars])];];unitGaussianQuadratureWeights[2]{{1, 1}}unitGaussianQuadraturePoints[2]{1/Sqrt[3], -(1/Sqrt[3])}I am looking for a paper/book that describes algorithmically how this is done for triangles and/or tetrahedra. Can someone point me to some information about this. Thanks.
Automatic generation of integration points and weights for triangles and tetrahedra
numerics;symbolic computation;reference request