3
$\begingroup$

I am trying to find the maximum value $t$ such that my function $h$ is bigger than $t$ on (at least) half of the interval $[0,2\pi]$ and smaller than $t$ on the other half, i.e.: $$ \text{sup}\{t\in \mathbb{R} : \lambda(\{ x\in [0, 2\pi]: h(x)\geq t\})>\pi\}.$$

In my case, the function $h$ is the integral of the difference between two other functions, $f$ and $g$. I am trying to calculate the length of the interval where $h(x) \geq t$ by integrating a Boolean. I have tried the following:

Clear[f, g, h, q];
f[x_?NumericQ] = Sum[(1 - Abs[k]/6)*E^(I*k*(x - (Pi/3))), {k, -5, 5}];
g[x_?NumericQ] = Sum[(1 - Abs[k]/6)*E^(I*k*(x - (5*Pi/3))), {k, -5, 5}];
h[t_?NumericQ] = Integrate[f[x] - g[x], {x, 0, t}]/(2*Pi);
q[z_?NumericQ, v_?NumericQ] = Boole[Chop[h[v]] <= Chop[h[z]]];
FindMaxValue[{Chop[h[z]], (Integrate[q[z, v], {v, 0, 2*Pi}]) >= Pi}, z]

Now I am getting the following errors:

NIntegrate::inumr: The integrand q[1.,v] has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,6.28319}}.

General::stop: Further output of NIntegrate::inumr will be suppressed during this calculation.

FindMaxValue::nrlnum: The function value {-[Pi]+!(*SubsuperscriptBox[([Integral]), (0), (2\ [Pi])](q[1.`, v] [DifferentialD]v))} is not a list of real numbers with dimensions {1} at {z} = {1.}.

I guess that the problem is that Mathematica is trying to evaluate q[1,v] which indeed is not a number. How do I get Mathematica to integrate this right?

$\endgroup$
3
  • 1
    $\begingroup$ Concerning your last remark, what the error message is telling you is that it evaluated things like q[1., 2.] and didn't get a number. (Here 2. represents one of the sampling points and may not be one of the actual values for v it used.) If you evaluate q[1., 2.] you will get another clue about what is going wrong. $\endgroup$ – Michael E2 Jan 26 at 15:46
  • 1
    $\begingroup$ May you please clarify what “on half of the interval” mean? The formula expresses “on at least half of the interval” instead. $\endgroup$ – SneezeFor16Min Jan 28 at 7:38
  • $\begingroup$ When the function h is continuous and nowhere constant, then the t I'm looking for will be such that h is larger than t on exactly half the interval. But you're right, it's not the same in general so I edited my question slightly. $\endgroup$ – Evert-Jan Hekkelman Jan 29 at 10:35
5
$\begingroup$

In this case,

  • ?NumericQ condition is not necessary in function definitions, but
  • ?NumericQ is needed for NIntegrate (ref. the documentation for message), so the integral should be split off, here to lambda[z_?NumericQ]:
    • using := for evaluation to happen only when arguments are passed in.
ClearAll[f, g, h, q, lambda];
f[x_] = Sum[(1 - Abs[k]/6)*E^(I*k*(x - (Pi/3))), {k, -5, 5}];
g[x_] = Sum[(1 - Abs[k]/6)*E^(I*k*(x - (5*Pi/3))), {k, -5, 5}];
h[t_] = Integrate[f[x] - g[x], {x, 0, t}]/(2*Pi) // FullSimplify // N;
q[z_, v_] := Boole[Chop[h[v]] <= Chop[h[z]]]
lambda[z_?NumericQ] := NIntegrate[q[z, v], {v, 0, 2 Pi}]
FindMaximum[{Chop[h[z]], 0 <= z <= 2 Pi && lambda[z] > Pi}, {z, Pi}]

{0.8821262327, {z -> 3.141592654}}

$\endgroup$
4
  • 1
    $\begingroup$ Sorry @SneezeFor16Min and @Evert-JanHekkelman, this can not be true. You simply get the total maximium of function h[t]. Plot[{0.8821262327, h[t]}, {t, 0, 2 Pi}] to see this. $\endgroup$ – Akku14 Jan 27 at 20:17
  • 1
    $\begingroup$ @Akku14 On reviewing the question, it seems that $\sup\{t:\lambda(t)>\pi\}$ doesn't correspond to what OP said above, which should be $\sup\{t:\lambda(t)=\pi\}$, exactly half. Otherwise $t_0$ with $\lambda(t_0)=2\pi$ is accepted and thus $h(t_0)$ becomes the simple maximum. $\endgroup$ – SneezeFor16Min Jan 28 at 7:35
  • $\begingroup$ Changing the last line to FindMaximum[{Chop[h[z]], 0 <= z <= 2 Pi && lambda[z] < Pi}, z] works! That returns {0.808616, {z -> 1.5708}}, which is correct. $\endgroup$ – Evert-Jan Hekkelman Jan 29 at 10:36
  • $\begingroup$ Then @Akku14 's interpretation and answer are also correct. :-) $\endgroup$ – SneezeFor16Min Jan 29 at 11:14
2
$\begingroup$

Is this what you want?

I interpreted: "I am trying to find the maximum value h[ts] such that my function h is bigger than h[ts] on half of the interval [0,2π] and smaller than h[ts] on the other half, i.e.:"

f[x_] = Sum[(1 - Abs[k]/6)*E^(I*k*(x - (Pi/3))), {k, -5, 5}];
g[x_] = Sum[(1 - Abs[k]/6)*E^(I*k*(x - (5*Pi/3))), {k, -5, 5}];
h[t_] = Integrate[f[x] - g[x], {x, 0, t}]/(2*Pi) // FullSimplify

(*   (4*(35 + 40*Cos[t] + 13*Cos[2*t] + 2*Cos[3*t])*
  Sin[t/2]^4)/(5*Sqrt[3]*Pi)   *)

{p1, p2} = {t1, t2} /. 
  First@Solve[{h[t2] == h[t1], t2 - t1 == Pi, 0 < t1 < Pi}, {t1, t2},
  Reals] // FullSimplify

(*   {Pi/2, (3*Pi)/2}   *)

Plot[{h[Pi/2], h[t]}, {t, 0, 2 Pi}, 
  Epilog -> Point@{{p1, h[p1]}, {p2, h[p2]}}
]

enter image description here

Edit Extension for functions with more humps.

Use UnitStep and provide the exact zeros of h[t]- hs (hs searched for) to get better integration results.

nsol[hs_?NumericQ] := 
  t /. NSolve[
  0 < t < 2 Pi &&   TrigExpand[h[t] - hs] == 0, t
]

nint[hs_?NumericQ] := 
 NIntegrate[UnitStep[h[t] - hs], 
  Evaluate[Flatten[{t, 0, nsol[hs], 2 Pi}]]]

FindRoot[nint[hs] == Pi, {hs, .7}]

(*   {hs -> 0.808616}   *)

h[Pi/2] // N[#, 15] &

(*   0.808615713285295   *)
$\endgroup$
2
  • $\begingroup$ That's a clever workaround, thanks for the idea! I am new to Mathematica and trying to learn, and I am wondering if you know a way to do such a thing if the function has more than one hump? Or would that become a complicated undertaking? $\endgroup$ – Evert-Jan Hekkelman Jan 26 at 12:27
  • 1
    $\begingroup$ @Evert-JanHekkelman please see my edit. $\endgroup$ – Akku14 Jan 26 at 14:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.