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 12(δn,m+δn,m). As discussed here, for instance, MMA aims to produce generically correct results, and the special case m=±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}14δk,nm14δk,mn+δk,mn4+δk,m+n4


check11[20]
(* True *)

Activate@check21
(* True *)

Comments

Popular posts from this blog

mathematical optimization - Minimizing using indices, error: Part::pkspec1: The expression cannot be used as a part specification

I want to use Minimize where the variables to minimize are indices pointing into an array. Here a MWE that hopefully shows what my problem is. vars = u@# & /@ Range[3]; cons = Flatten@ { Table[(u[j] != #) & /@ vars[[j + 1 ;; -1]], {j, 1, 3 - 1}], 1 vec1 = {1, 2, 3}; vec2 = {1, 2, 3}; Minimize[{Total@((vec1[[#]] - vec2[[u[#]]])^2 & /@ Range[1, 3]), cons}, vars, Integers] The error I get: Part::pkspec1: The expression u[1] cannot be used as a part specification. >> Answer Ok, it seems that one can get around Mathematica trying to evaluate vec2[[u[1]]] too early by using the function Indexed[vec2,u[1]] . The working MWE would then look like the following: vars = u@# & /@ Range[3]; cons = Flatten@{ Table[(u[j] != #) & /@ vars[[j + 1 ;; -1]], {j, 1, 3 - 1}], 1 vec1 = {1, 2, 3}; vec2 = {1, 2, 3}; NMinimize[ {Total@((vec1[[#]] - Indexed[vec2, u[#]])^2 & /@ R...

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...

What is and isn't a valid variable specification for Manipulate?

I have an expression whose terms have arguments (representing subscripts), like this: myExpr = A[0] + V[1,T] I would like to put it inside a Manipulate to see its value as I move around the parameters. (The goal is eventually to plot it wrt one of the variables inside.) However, Mathematica complains when I set V[1,T] as a manipulated variable: Manipulate[Evaluate[myExpr], {A[0], 0, 1}, {V[1, T], 0, 1}] (*Manipulate::vsform: Manipulate argument {V[1,T],0,1} does not have the correct form for a variable specification. >> *) As a workaround, if I get rid of the symbol T inside the argument, it works fine: Manipulate[ Evaluate[myExpr /. T -> 15], {A[0], 0, 1}, {V[1, 15], 0, 1}] Why this behavior? Can anyone point me to the documentation that says what counts as a valid variable? And is there a way to get Manpiulate to accept an expression with a symbolic argument as a variable? Investigations I've done so far: I tried using variableQ from this answer , but it says V[1...