# Why can't Mathematica evaluate this integral?

I want to work with the rectangle function, which I define by

f[x_, m_] := Limit[1/((2*(x - m))^(2*k) + 1), k -> Infinity];


(I know that in theory I can use HeavisidePi[x-m] to obtain the same function, but using my definition evaluates to 1/2 at the discontinuities, which I want to retain.)

Clearly, the integral of f with respect to x over {x, m - 1/2, m + 1/2} should evaluate to 1. However, I can't persuade Mathematica to define the integral at all: the code

Clear["Global*"];
f[x_, m_] := Limit[1/((2*(x - m))^(2*k) + 1), k -> Infinity];
Module[{m = 0}, Plot[f[x, m], {x, -1, 1}, GridLines -> Automatic]]
Module[{m = 0}, Integrate[f[x, m], {x, m - 1/2, m + 1/2}]]


yields the desired plot, but returns Undefined for the integral.

I have tried adding Assumptions but clearly not the right ones, because it doesn't work. I have tried using NIntegrate, but it doesn't work. I have tried changing the definition of f to HeavisidePi[x-m] (which I don't want to do anyway), but it doesn't work. I have tried using a Piecewise definition - same non-result...

How do I obtain the correct result?

• Does f[x, 0] evaluate correctly? Oct 16 '20 at 16:13
• That limit evaluates to ConditionalExpression[0, Log[2 x] > 0]. It's not clear to me what Integrate will be able to do with it. I will suggest that if you want a Heaviside function but with a value at the origin, create it with Piecewise. Oct 16 '20 at 16:28
• Why not using HeavisidePi? It keeps HeavisidePi[1/2] and HeavisidePi[-1/2] unevaluated. If some of your outputs contain HeavisidePi[1/2] you can simply replace it with desired value - like HeavisidePi[1/2]->1/2 or HeavisidePi[-1/2]->1/2. Oct 16 '20 at 16:36
• Hi all. In answer to your questions: (1) f[0,0] evaluates correctly to 1; (2) @Daniel, I'm not sure I understand the reservation: at least until the limit ->Infinity is reached, f is an entire (and therefore integrable) function - and in the limit the definite integral described clearly produces 1 (having said this, I did try Piecewise but that didn't work either, any least with the half-values at the discontinuities); (3) As mentioned in the OP, HeavisidePi also produces an Undefined result under integration. Oct 16 '20 at 16:44
• When did you get Undefined? Integrate[HeavisidePi[x - m], {x, m - 1/2, m + 1/2}] /. m -> 0 nicely evaluates into 1. Oct 16 '20 at 16:53

Do integration first and then take the limit k-> inf.

Use indefinite integration. Integrate  only finds a solution for m==0.

g[x_, m_, k_] = 1/((2*(x - m))^(2*k) + 1)

mint[x_, k_] = Integrate[g[x, 0, k], x]

(*   x Hypergeometric2F1[1, 1/(2 k), 1 + 1/(2 k), -4^k x^(2 k)]   *)


Use a trick. Tell Limit  that x^(2 k) is always positive. (Valid for Integer k) That is only neccessary for x== -1/2. For all other positive and negative x Limit[mint[x, k], k -> \[Infinity]]  is x. Maybe this is a removable singularity at x== -1/2 ?

mint2[x_, k_] = mint[x, k] /. x^(2 k) -> Abs[x]^(2 k)

Limit[mint2[1/2, k] - mint2[-1/2, k], k -> \[Infinity]]

(*   1   *)


A graphic shows the same.

Manipulate[Plot[mint[x, k], {x, -.5, .5}], {k, 1, 1000}]


Rubi (https://rulebasedintegration.org/) does the integral with arbitrary m.

rint[x_, m_, k_] = Int[g[x, m, k], x]

(*   (-m + x) Hypergeometric2F1[1, 1/(2 k),
1/2 (2 + 1/k), -4^k (-m + x)^(2 k)]   *)

rint2[x_, m_, k_] =
rint[x, m, k] /. (-m + x)^(2 k) -> Abs[(-m + x)]^(2 k)

Limit[rint2[m + 1/2, m, k] - rint2[m - 1/2, m, k], k -> \[Infinity]]

(*   1   *)


Does this help?

f2[x_?NumericQ, m_] := Limit[1/((2*(x - m))^(2*k) + 1), k -> Infinity];

Module[{m = 0}, NIntegrate[f2[x, m], {x, m - 1/2, m + 1/2}]]

(* 1. *)
`
• Many thanks to you both. Oct 16 '20 at 22:13