Skip to main content

calculus and analysis - Solution to a specific problem caused by generic simplification


I'm trying to get MMA to help me evaluate certain integrals of trig functions. Here is an example: (The actual expressions I want to evaluate are more complicated than this one, but this illustrates the problem.)



Assuming[
n ∈ Integers && m ∈ Integers,
Integrate[Cos[n π x] Cos[m π x], {x, 0, 1}]
]
(* 0 *)

This answer is of course wrong, strictly speaking. The correct answer is $\frac{1}{2}\left(\delta_{n,m}+\delta_{n,-m}\right)$. As discussed here, for instance, MMA aims to produce generically correct results, and the special case $m=\pm n$ ends up being overlooked. I understand all that.


My question is whether anyone can suggest a straightforward workaround to get MMA to produce a more generally correct result for slightly more complex cases such as Integrate[Sin[k π x] Cos[n π x] Sin[m π x], {x,0,1}].


To clarify, I will add that there is no difficulty evaluating the special cases, if you know what they are. For instance:


Assuming[

n ∈ Integers,
Integrate[Cos[n π x] Cos[n π x], {x, 0, 1}]
]
(* 1/2 *)

So this is one of those vexing questions where it's easy to find the answer, once you know what it is. And of course you can use trig identities to get the answer, but getting the signs right is a tedious, fiddly business. I'd like to let MMA do it for me.


In case anyone cares, these integrals arise from PDEs when the solutions are represented as a cosine series.


Thanks.



Answer



Here is a fully automated solution, using ideas from my previous answer and Dr. Hintze's. It is in the form of two functions: trigIntegralToPiecewise and piecewiseToDelta. trigIntegralToPiecewise takes the thing to be integrated (as a pure function of the integration variable) and a list of the integer variables. It returns a Function that evaluates to a Piecewise function of the variables. If trigIntegralToPiecewise returns successfully, this piecewise function should always be equivalent to the integral.



The function can then be used as input to piecewiseToDelta. This is a bit of a kludge but seems to work fairly well. It is called using something like:


{df, check1, check2} = piecewiseToDelta[pwf, {k, m, n}];

df1 is now a pure function of k, m, n that evaluates to a sum of Kronecker deltas. (This is useful not only because the delta form is typically more compact and easier to understand than the Piecewise form, but also because MMA is good at simplifying sums over Kronecker deltas.) It is not guaranteed to be equivalent to pwf. The check1 and check2 returns allow you to reassure yourself. check1 is an Inactive expression using Reduce that, activated, should evaluate to True if the two functions are equivalent. I have not always had good results with Reduce, so check2 is something else. check2[10] evaluates pwf[k,m,n]==df[k,m,n] for $k$, $m$, and $n$ from -10 to 10 and returns True if all equalities hold.


trigIntegralToPiecewise::nonzero = 
"works only for generically zero integrals";

trigIntegralToPiecewise[
func_,
vars_

] := Module[{ri, num, den, dvs, cases, res, pw, pwfunc},
ri = Integrate[func[x], {x, 0, 1}];
num = Numerator[Together[ri]];
den = Denominator[Together[ri]];
If[0 =!= Simplify[num, Assumptions -> vars \[Element] Integers],
Message[trigIntegralToPiecewise::nonzero];
Return[$Failed]
];
cases = BooleanConvert[Reduce[0 == den, vars]];
cases = If[Head[cases] === Or,

List @@ cases,
{cases}
];
cases = Reverse[Subsets[cases]];
cases = And @@@ cases;
cases = DeleteDuplicates[Reduce[#, vars] & /@ cases];
pw = Table[
res = Simplify[
Integrate[func[x], {x, 0, 1},
Assumptions -> case],

vars \[Element] Integers && case
];
{res, case},
{case, cases}
];
pw = Append[pw, {0, True}]; (* To be safe *)
pwfunc = Function @@ {
vars,
Piecewise[
pw

]
};
pwfunc
]

piecewiseToDelta[
pwf_,
vars_,
max_: 3,
extrabasis_: {}

] := Module[{nv, basis, kds, iterators, lb, ub, cs, eqs, cvals,
dfunc, check1, check2},
nv = Length[vars];
iterators =
Transpose[{vars, ConstantArray[lb, nv], ConstantArray[ub, nv]}];
basis = Join[{1}, vars, extrabasis];
kds = Plus @@@ Tuples[Outer[Times, Drop[vars, 1], {1, -1}]];
kds = Sum[
(Sum[C[i, j] basis[[j]], {j, Length[basis]}]) KroneckerDelta[
vars[[1]], kds[[i]]],

{i, Length[kds]}
];
cs = Flatten@Table[
C[i, j],
{j, Length[basis]}, {i, Length[kds]}
];
eqs = Table @@
Prepend[iterators /. {lb -> -max, ub -> max},
kds == pwf @@ vars];
eqs = DeleteCases[Flatten[eqs], True];

cvals = Solve[eqs, cs][[1]];
dfunc = Function @@ {
vars,
kds /. cvals /. C[_, _] -> 0
};
check1 = Inactivate[
Assuming[
vars \[Element] Integers,
Simplify[
Reduce[pwf @@ vars == dfunc @@ vars]

]
]
];
check2[cmax_, it_: iterators] := Module[{iters},
iters = it /. {lb -> -cmax, ub -> cmax};
And @@ Flatten[
Table @@ Prepend[
iters,
pwf @@ vars == dfunc @@ vars
]

]
];
{dfunc, check1, check2}
]

Example:


pwf1 = trigIntegralToPiecewise[x ↦ Sin[k 𝝿 x] Cos[n 𝝿 x] Sin[m 𝝿 x], {k, m, n}];

{df1, check11, check21} = piecewiseToDelta[pwf1, {k, m, n}];
df1 // TraditionalForm


$$\{k,m,n\}↦-\frac{1}{4} \delta _{k,n-m}-\frac{1}{4} \delta _{k,-m-n}+\frac{\delta _{k,m-n}}{4}+\frac{\delta _{k,m+n}}{4}$$


check11[20]
(* True *)

Activate@check21
(* True *)

Comments

Popular posts from this blog

front end - keyboard shortcut to invoke Insert new matrix

I frequently need to type in some matrices, and the menu command Insert > Table/Matrix > New... allows matrices with lines drawn between columns and rows, which is very helpful. I would like to make a keyboard shortcut for it, but cannot find the relevant frontend token command (4209405) for it. Since the FullForm[] and InputForm[] of matrices with lines drawn between rows and columns is the same as those without lines, it's hard to do this via 3rd party system-wide text expanders (e.g. autohotkey or atext on mac). How does one assign a keyboard shortcut for the menu item Insert > Table/Matrix > New... , preferably using only mathematica? Thanks! Answer In the MenuSetup.tr (for linux located in the $InstallationDirectory/SystemFiles/FrontEnd/TextResources/X/ directory), I changed the line MenuItem["&New...", "CreateGridBoxDialog"] to read MenuItem["&New...", "CreateGridBoxDialog", MenuKey["m", Modifiers-...

How to thread a list

I have data in format data = {{a1, a2}, {b1, b2}, {c1, c2}, {d1, d2}} Tableform: I want to thread it to : tdata = {{{a1, b1}, {a2, b2}}, {{a1, c1}, {a2, c2}}, {{a1, d1}, {a2, d2}}} Tableform: And I would like to do better then pseudofunction[n_] := Transpose[{data2[[1]], data2[[n]]}]; SetAttributes[pseudofunction, Listable]; Range[2, 4] // pseudofunction Here is my benchmark data, where data3 is normal sample of real data. data3 = Drop[ExcelWorkBook[[Column1 ;; Column4]], None, 1]; data2 = {a #, b #, c #, d #} & /@ Range[1, 10^5]; data = RandomReal[{0, 1}, {10^6, 4}]; Here is my benchmark code kptnw[list_] := Transpose[{Table[First@#, {Length@# - 1}], Rest@#}, {3, 1, 2}] &@list kptnw2[list_] := Transpose[{ConstantArray[First@#, Length@# - 1], Rest@#}, {3, 1, 2}] &@list OleksandrR[list_] := Flatten[Outer[List, List@First[list], Rest[list], 1], {{2}, {1, 4}}] paradox2[list_] := Partition[Riffle[list[[1]], #], 2] & /@ Drop[list, 1] RM[list_] := FoldList[Transpose[{First@li...

functions - Get leading series expansion term?

Given a function f[x] , I would like to have a function leadingSeries that returns just the leading term in the series around x=0 . For example: leadingSeries[(1/x + 2)/(4 + 1/x^2 + x)] x and leadingSeries[(1/x + 2 + (1 - 1/x^3)/4)/(4 + x)] -(1/(16 x^3)) Is there such a function in Mathematica? Or maybe one can implement it efficiently? EDIT I finally went with the following implementation, based on Carl Woll 's answer: lds[ex_,x_]:=( (ex/.x->(x+O[x]^2))/.SeriesData[U_,Z_,L_List,Mi_,Ma_,De_]:>SeriesData[U,Z,{L[[1]]},Mi,Mi+1,De]//Quiet//Normal) The advantage is, that this one also properly works with functions whose leading term is a constant: lds[Exp[x],x] 1 Answer Update 1 Updated to eliminate SeriesData and to not return additional terms Perhaps you could use: leadingSeries[expr_, x_] := Normal[expr /. x->(x+O[x]^2) /. a_List :> Take[a, 1]] Then for your examples: leadingSeries[(1/x + 2)/(4 + 1/x^2 + x), x] leadingSeries[Exp[x], x] leadingSeries[(1/x + 2 + (1 - 1/x...