1
$\begingroup$

I encounter a problem when I try to integrate the piecewise function defined in the following code:

 ClearAll["Global`*"];

(* Parameters of the problem*)

M = 1; (* Mass scale*)
A = M/4 ;(* Background Amplitude*)
a = A; (* Perturbation Amplitude*)
L = 0.1/M; (* Background Length scale*)
Lb = 30*L; (* Boundary distance*)
l = L/20; (* Perturbation Length scale*)
d = L/3; (* Spatial translation of the perturbation *)
tmax = 10*L; (* Maximal time*)
\[Kappa] = 3*M; (*Velociy field stiffness*)
\[Delta]\[Chi]0[x_] = 
 a/(1 + ((x - d)/l)^2) + 
  a/(1 + ((x + d)/l)^2); (* Initial perturbation profile*)

(* Construct an initial characteristic with the correct asymptotical \
behaviour *)   

f1[x_] = -\[Kappa]*
  x; (* Desired profile close to x=0 to create the caustic *)
f2[x_] = 1; (* Asymptotical profile *)
v[x_] = Piecewise[{
   {f1[x], -L < x < L}, 
   {f2[x], x < -Lb/2 ||  x > Lb/2},
   {Interpolation[{{{-Lb/2}, f2[-Lb/2], f2'[-Lb/2], f2''[-Lb/2], 
        f2'''[-Lb/2], f2''''[-Lb/2]}, {{-L}, f1[-L], f1'[-L]}}][
     x], -Lb/2 < x < -L},
   {Interpolation[{{{Lb/2}, f2[Lb/2], f2'[Lb/2]}, {{L}, f1[L], 
        f1'[L]}}][x], L < x < Lb/2}}];


(* Define the initial condition for the field *)

s[x_] = M^2*Sqrt[(2 - 2*(v[x])^2)/(3 - (v[x])^2)];
q[x_] = Integrate[s[x], x, 
   Assumptions -> -Lb < x < Lb]; (* Integrate the spatial derivative*)

The piecewise function is called v[x] and it is defined explicitly on the intervals: [-L,L] ; ]-inf, Lb/2] and [Lb/2, +inf[ and defined as an interpolation on the remaining intervals. Then I define s[x] by using this piecewise function v[x]. When I try to define a function q[x] as the integral of s[x] , I have the following error:

Integrate::ilim: Invalid integration variable or limit(s) in -2.99988.

I have read the documentation on this error but I don't manage to relate this to my problem. Can anybody help me please?

Thanks a lot.

$\endgroup$
7
  • $\begingroup$ Do you think you could reduce this into a minimal-working-example? $\endgroup$ May 12, 2018 at 18:29
  • $\begingroup$ I think you mean Integrate[s[x],{ x, -Lb , Lb}] rather than Integrate[s[x], x, Assumptions -> -Lb < x < Lb] $\endgroup$
    – andre314
    May 12, 2018 at 18:49
  • $\begingroup$ Integrate[s[x],{x,-.1,.1}] and Integrate[s[x],{x,1.5,3}] and Integrate[s[x],{x,-3,-1.5}] appear to work with the rest of your code unchanged, but anything which strays into your Interpolation zones fails. Does that give you any ideas what to start looking at? Defining q[x_] and then using x as a variable of integration worries me. $\endgroup$
    – Bill
    May 12, 2018 at 18:50
  • $\begingroup$ @AccidentalFourierTransform : if you ignore the definition of the parameters, the example is really short. Please read from (* Construct an initial characteristic with the correct asymptotical \ behaviour *) $\endgroup$
    – Targamas
    May 12, 2018 at 20:42
  • $\begingroup$ @andre Actually I want to define the function q[x] as an undefined integral but with [-Lb,Lb] as domains $\endgroup$
    – Targamas
    May 12, 2018 at 20:44

2 Answers 2

2
$\begingroup$

You can use NDSolve to compute a numerical antiderivative:

q1 = NDSolveValue[{qq'[x] == s[x], qq[-Lb] == 0}, qq, {x, -Lb, Lb}]
(* Some extrapolation warnings that are due to how NDSolve rewrites Piecewise[].
   They are unimportant. *)

ListLinePlot@q0

Mathematica graphics

Since the interpolations are between just two points, one can use InterpolatingPolynomial for slightly better results. (The OP's Interpolation leads to some small imaginary parts, due to a square-root of slightly negative numbers. For some reason, we don't get that with the InterpolatingPolynomial. The negative numbers are due to an overshoot around x == 1.5 (try 1 - (v[-1.5 + 10^-6])^2).

v1[x_] = Piecewise[{{f1[x], -L < x < L},
      {f2[x], x < -Lb/2 || x > Lb/2},
     {InterpolatingPolynomial[
      {{{-Lb/2}, f2[-Lb/2], f2'[-Lb/2], f2''[-Lb/2], f2'''[-Lb/2],
       f2''''[-Lb/2]}, {{-L}, f1[-L], f1'[-L]}}, x],
      -Lb/2 < x < -L},
    {InterpolatingPolynomial[{{{Lb/2}, f2[Lb/2], f2'[Lb/2]}, {{L}, f1[L], f1'[L]}}, x],
     L < x < Lb/2}}];

Block[{v = v1},
 q2 = NDSolveValue[{qq'[x] == s[x], qq[-Lb] == 0}, qq, {x, -Lb, Lb}]
 ]
$\endgroup$
0
$\begingroup$

It seems, MMA can't do Integration of squared interpolation functions.

Generate a new interpolationfunction of the squared function values.

Plot[(v[x])^2, {x, -5, 5}, PlotRange -> All]

enter image description here

s2 = Interpolation[
       Table[{x, M^2*Sqrt[(2 - 2*(v[x])^2)/(3 - (v[x])^2)]}, {x, -3, 3, 
              1/10}]];

Plot[s2[x], {x, -3, 3}, PlotRange -> All]

enter image description here

is2[x_] = Integrate[s2[x], x]

Plot[is2[x], {x, -3, 3}]

enter image description here

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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