Simplifying products of DiracDelta

modified

If we look at a simpler example , with only two unknowns (for reasons of visualization)

DiracDelta[z-2]DiracDelta[k-z]

you would conclude

DiracDelta[z-2]DiracDelta[k-z]=?=DiracDelta[z-2]DiracDelta[k-2]

Numerical visualization shows

dirac = Function[{x, eps}, Exp[-(x^2/(2 eps))]/Sqrt[2 Pi eps]] Plot3D[dirac[z - 2, #] dirac [k - 2, #] - (dirac[z - 2, #] dirac [k - z, #]), {z, 1, 3}, {k,1, 3}, PlotRange -> All, PlotPoints -> 100] &[.01] 

difference (distributial difference) between these two expressions!

Final check

NIntegrate[dirac[z - 2, #] dirac[k - 2, #] - (dirac[z - 2, #] dirac[k - z, #]), {z, 1, 3}, {k, 1,3}, Method -> "LocalAdaptive"] &[.01] (*-1.46645*10^-8*) 

evaluates to zero, though conclusion seems to be is correct!

addendum

Back to your original question:

Catch the arguments of DiracDelta's and force them to be zero

Cases[DiracDelta[z - 2] DiracDelta[k - 5] DiracDelta[t - z - k - 9],DiracDelta[p_] :> p, -1] sol = Solve[% == 0, Variables[%]][[1]] Map[DiracDelta[#] &, sol /. Rule -> Subtract] /. List -> Times (*DiracDelta[-5 + k] DiracDelta[-16 + t] DiracDelta[-2 + z]*) 

You probably are supposing, the set of Schwarz distributions is forming an algebra with respect to the set of operation (+,*, evaluate at a point) as we use the algebra of functions with rules of addition and multiplication.

Its an extension of functional algebra with severe limitations for distributions with discontiniuties.

This idea rest on the assumptions, that as for

$$(f⊕g)(x)=f(x)+g(x),(f⊗g)(x)=f(x)∗(g(x))$$

which is going ot work everywhere except for local evalution on the support of discontinuities of any of the distribution. The simples non working examples ar $\delta(x) \Theta[x], \delta^2(x)$, while $\theta^2 = \theta$ is working.

The central point is order of evaluation and the continuity there on different paths to the evaluation point.

$$ \delta_x( f ) = f(x), \delta_x (\delta_y f) = \delta_x (f(y)) = f(y)(x) != f(y,x)$$ $$ \delta_y (\delta_x f) = \delta_y (f(x)) = f(x)(y) != f(x,y)$$ This problem of evaluation order with definite order of the arguments of a function by named symbols is not only conceptual , but fundamental.

Incontrast to vs 6 (the deleted post), at least since vs 13 Mathematica the problem is known. The list of all 6 permutations of evaluations of an integral with a function of three named variables with definite positions: The integration remains unevaluated, if two $\delta$'s with the same argument collide in a step of integration per dimension.

The problem of multidimensional integration is tied to the theorem of Fubini.

In vs 13.3 the evaluation on any of the different orders of integration yields

 ((# -> (int = Integrate[ DiracDelta[z - 2] DiracDelta[k - 5] DiracDelta[t-z-k-9] * f[z, k, t], #1, #2, #3] &) @@ Permute[{{z, -\[Infinity], \[Infinity]}, {k, -\[Infinity], \[Infinity]}, {t, -\[Infinity], \[Infinity]}}, #] &) /@ Permutations[{1, 2, 3}]) /. {_Integrate :> \[Integral]\[Delta][x] f[x] \[DifferentialD]x} 

$$\{1,2,3\}\ \to \ f(2,5,16),\qquad \{1,3,2\}\to \int f(x) \delta (x) \, dx,\qquad \{2,1,3\}\ \to \ f(2,5,16) $$

$$ \{2,3,1\}\ \to \ \int f(x) \delta (x) \, dx,\qquad \{3,1,2\}\ \to \ \int f(x) \delta (x) \, dx, \qquad \{3,2,1\}\ \to \ \int f(x) \delta (x) \, dx$$

History: When Heaviside introduced distributional algebra, for some time, he was the only one able to apply his set of rules producing mathematically correct results. Then he started to study pure mathematics.

Only after Laurent Schwarz introduced distributions as vectors, elements in the dual linear space of evaluation procedures $\mathit S'$ on his Schwarz vector space of smooth functions with compact support $\mathit S_0$, it became clear, what a complex machinery Heaviside had developed in his distribution calculus.

By this theory, product of distributions of different variables are identified as evaluators at a point by their Fourier transformation formulas. Disturbing: The space of distributions contains nearly all integrable, piecewise functions as a core. But limits of series of such functions not converging pointwise everywhere, too.

The multidimensional Schwarz distribution calculus rests on the fact, that for all integrals, order of integration is free. That is Fubini and is working on the space of smooth functions, integrals over all of 1..n-space existing, only.

Here's an extension of @Ulrich's solution, taking into account the Jacobian of the transformation:

diracProductReduce // ClearAll; diracProductReduce::nonlinear = "Some DirectDelta arguments are nonlinear."; diracProductReduce::nonsquare = "The number DirectDelta functions does not equal the number of \ variables `1`."; diracProductReduce[ prod : Verbatim[Times][deltas : DiracDelta[_] ..]] := Module[{args}, args = First /@ {deltas}; diracProductReduce[prod, Variables[args]] ]; diracProductReduce[prod : Verbatim[Times][deltas : DiracDelta[_] ..], vars_] := Module[{args, ca, bvec, jacobian, scale, vars0}, args = First /@ {deltas}; If[Length[args] == Length[vars], ca = CoefficientArrays[args, vars]; If[Length[ca] == 2, {bvec, jacobian} = ca; scale = Det[jacobian]; If[scale =!= 0 && scale =!= 0., vars0 = Simplify@LinearSolve[jacobian, -bvec]; (* Result *) ConditionalExpression[ 1/ Abs@scale (Times @@ MapThread[DiracDelta[#1 - #2] &, {vars, vars0}]), scale != 0 // Simplify], Message[diracProductReduce::sing, jacobian]; (* Failed result *) prod ], Message[diracProductReduce::nonlinear]; (* Failed result *) prod ], Message[diracProductReduce::nonsquare, vars]; (* Failed result *) prod ] ]; 

Examples:

diracProductReduce[ DiracDelta[x - 5] DiracDelta[t - y - 1] * DiracDelta[t + z + y - 3] DiracDelta[t + z + x - y - 3]] (* 2 DiracDelta[-7 + 2 t] DiracDelta[-5 + x] * DiracDelta[-5 + 2 y] DiracDelta[3 + z] *) diracProductReduce[ (* specify variables (t is a parameter) *) DiracDelta[x - 5] DiracDelta[t - y - 1] DiracDelta[t + z + y - 3], {x, y, z}] (* DiracDelta[-5 + x] DiracDelta[1 - t + y] DiracDelta[-4 + 2 t + z] *) diracProductReduce[ (* conditional reduction *) DiracDelta[x - 5] DiracDelta[t - y - 1] DiracDelta[t z + y - 3], {x, y, z}] (* ConditionalExpression[( DiracDelta[-5 + x] DiracDelta[1 - t + y] DiracDelta[1 - 4/t + z])/ Abs[t], t != 0] *)