# Tag Info

35

General comments First, if you plan to use multi-dimensional integrals it is better to test with multi-dimensional integrals not with one dimensional ones. One might think that the test in the question is an appropriate one if multi-dimensional integration is done by the integrator in a recursive manner. This seems to be case for scipy.integrate.nquad (see ...

19

I'll preface this answer first with a complaint: NExpectation[] and NProbability[] are not sufficiently resilient obviously adjustable. Ideally, these two functions are an "interface" to NIntegrate[], allowing the user to formulate his expression purely in distributional terms. Unfortunately, when one hits cases like this, the things one might usually ...

17

This question comes up often enough. See this discussion at community.wolfram.com : Integration method used in NIntegrate , and the notebook Finding the applied NIntegrate methods attached to my second response in the discussion. That notebook contains examples of usage of the undocumented function NIntegrateSamplingPoints and NIntegrate's option ...

15

The correct syntax is NIntegrate[1, {x} ∈ ImplicitRegion[(x > 5 && x < 9) || (x > 11 && x < 13), x], Method -> "MonteCarlo"] The {x} has moved out in front. Alternatively you can do: NIntegrate[Boole[(x > 5 && x < 9) || (x > 11 && x < 13)], {x, 5, 13}, Method -> "MonteCarlo"] Also, if you ...

15

It turns out ListSurfacePlot3D does a terribly poor job of approximating the surface in the OP, otherwise one will just apply DiscretizeGraphics to the output obtained from ListSurfacePlot3D and be done with it. But since that's not applicable here, we present an approach that uses alpha shapes to approximate the shape of the given point set by tuning a ...

15

The memory leak in NIntegrate is a bug and has been fixed as of version 10.2.0. Earlier versions would lose around 720 bytes per evaluation for this example, which could not be recovered without restarting the kernel. ClearSystemCache[] should be used to make sure the memory is released. Using version 10.2: NI[z_?NumericQ, b0_?NumericQ] := ...

14

Introduction I think there are several questions on this site about ODEs of the form $$(x-a)^2 u''(x) = F(x,u,u')$$ with an initial condition at $x=a$. There is no general guarantee that solutions exist over an interval $(a,b]$, but sometimes it is possible as in this case. Outline We transform the equation $u''(x) = F(x,u,u')$ over the infinite interval ...

14

In my experience FindRoot works best for such problems: In[1]:= fun[a_?NumericQ] := NIntegrate[(x/(Exp[x] - 1)), {x, 0, a}] In[2]:= FindRoot[1 - (2/a)*((1/a)*fun[a] + (a/6) - 1) == 0.7, {a, 0.1}] Out[2]= {a -> 58.3073}

13

NDSolve currently can't handle this kind of differential equation, LaplaceTransform is your friend. Since in this case inverse Laplace transform can't be done analytically by InverseLaplaceTransform, you need the help of numerical Laplace inversion package in addition: eq = {0.01 - 6.25 x[t] + (1.2 Integrate[x'[t - τ]/Sqrt[τ], {τ, 0, t}])/10^7 == 16 ...

13

As a minor addition to J.M.'s excellent answer, If something breaks while evaluating NExpectation[] or NProbability[], a number of the things that otherwise can be adjusted in NIntegrate[] aren't there. Options can be passed to NIntegrate, for example try something like Table[NExpectation[X, X \[Distributed] JohnsonDistribution["SB", γ, δ, 0, 1], ...

13

The main problem is that your pos is not seen as a 3D vector. The cross product is therefore interpreted as a scalar: q*Cross[D[pos[t], t], b] when adding this to the vector q.e this 'scalar' term is added to each of the vector components: q*e + q*Cross[D[pos[t], t], b] This won't work, instead do: b = {1, 0, 0}; e = {0, 0, 1}; q = 1; m = 1; ...

12

$Version (* Out[228]= "8.0 for Microsoft Windows (64-bit) (October 7, 2011)" *) The result provided by Mathematica is correct: Integrate[a/(a^2 + Sin[t]^2), {t, 0, 2 π}] (* Out[213]= (2 π)/(Sqrt[1 + 1/a^2] a) *) Now the same procedure as "always" which "explains" the zero result. The indefinite integral is Integrate[a/(a^2 + Sin[t]^2), t] (* ... 12 For what it's worth, you can monitor the progress of Integrate and NIntegrate. I'm not sure how helpful the tools below are, but I feel it is probably worth mentioning them. Integrate InternalIntegratedebugSwitch If you set InternalIntegratedebugSwitch to the magic number 10, it will print its (major) steps in searching for the answer. For instance: ... 12 If you're running Mathematica 10 or above, there's a dead simple method using the new Area function: newimg = img /. Ball[x___, 0.6] -> Sphere[x, 0.5] Total[Map[Area[RegionIntersection[#, Cuboid[{-1, -1, -1}, {1, 1, 1}]]] &, newimg]] (* 25.1327 *) Note a few things here: Mathematica considers Ball to be a 3D object; if we want its surface area ... 12 WhenEvent is working. Try WhenEvent[x[t] < 0, Print[t]; x[t] -> 0] to see that every crossing is detected. The problem is that it only detects crossings. So changing x[t] -> 0 does not reset the event. At the next step x[t] becomes negative and no event is detected. (This is how it is supposed to behave.) The way to deal with this is to use ... 12 My question is: how to set that NDSolve will not save whole InterpolationFunction for the result? There is actually a very simple way to do this: instead of specifying a list of functions in the second argument, specify an empty list instead. This now begs the question of how one can obtain results. The solution is to use the event location ... 12 Some explanations first The substitution in the question introduces the reduced wave function$u(r)$by solving the original radial equation in polar coordinates, $$-\frac{1}{2}\left(R''(r)+\frac {1}{r}R'(r)\right) - \frac{1}{r}R(r) + \frac {m^2}{2r^2}R(r) = E R(r)$$ using the ansatz $$R(r)\equiv \frac{1}{\sqrt{r}}u(r)$$ The apparently divergent ... 11 I've got a theory.... Let's try to get the antiderivative: Integrate[(1 + 16*Tan[2*x - y]^2)/(1 + 4*Tan[2*x - y]^2), x, Assumptions -> y \[Element] Reals] (*returns 5 x - (5 y)/2 - ArcTan[2 Tan[2 x - y]] *) Seems legit. You can test with D[5 x - (5 y)/2 - ArcTan[2 Tan[2 x - y]], x] and plot or rearrange. This antiderivative is technically correct, ... 11 Working solution One can manually implement the shooting method with ParametricNDSolveValue and FindRoot: psol = ParametricNDSolveValue[{k f[z] + f'[z] == 0, f[-1] == Exp[2], WhenEvent[z == 0, f[z] -> r0 f[z]]}, f, {z, -1, 1}, {k, r0}]; k0 = k /. FindRoot[psol[k, 0.5][1] == 1, {k, 1}] (* 0.653426 *) Plot[psol[k0, 0.5][z], {z, -1, 1}] ... 11 s = ParametricNDSolveValue[{x'[t] == -y[t] + x[t]*Log[x[t]], y'[t] == x[t] + y[t]*Log[x[t]], x[0] == x0, y[0] == 0}, {x, y}, {t, 1}, x0] f[x0_, t_] := Through[Through[s@x0]@t] pts = Table[f[x0, t], {x0, 1, 2, .2}, {t, 0, 1, .1}]; Show[Graphics[{Green, Arrow /@ pts, Black, Point /@ pts}, ... 10 It has to do with the default behavior of GenerateConditions in multivariate integrals. Setting it explicitly to True will help in this case. Some explanation may be found here or here. The gist is that, for multiple integration, automatically checking conditions and issuing provisos for all but the final integration is typically both too costly (in speed) ... 10 (Update: I forgot to copy some of the code) For values of z = x beyond the critical value (just below z == 7500), there are three real equilibria, two stable spirals and one (unstable) saddle. eq[x_] := Module[{q = x}, w0 = 7000; G0 = 50; Q = 14000; hU = 0.6*G0; w00 = w0 - q + Q*Abs[A0[t]]^2; At = 10^-6; G1 = 1.0 G0; tmax = 4; e1 ... 10 EDIT #2 My error was useful. It brought me to the conclusion that the difficulties in solving the PDE of the OP are due to the drift term $$\frac{\partial (x u(x,t))}{\partial x}$$ If the drift term is included, many boundary problems are ill defined. It turns out that there are cases where mathematically there is only a trivial solution u = 0 but ... 9 Mathematica seems to split the integrand component, E^(-((-m + Log[x])^2/(2 s^2))) into E^(-((m^2 + Log[x]^2)/(2 s^2))) times the sort-of "coefficient" E^((m Log[x])/s^2) (* == x^(m/s^2) *) in order to calculate the integral in terms of Meijer$G$. For reasons that are obscure to me, it seems to want the coefficient of m in the exponent to be ... 9 First, in general, I would advise you not to trust numerical algorithms. If there are doubts about the outcomes then solve the same problem with different (numerical or not) methods and see do their results agree. For the integral in the question I assume you can evaluate it with several different invocations of the Monte Carlo method and compare the ... 9 This is the best I can come up with, I'm very interested to see if anyone else has a better solution. The idea here is to just run through values of$t$, and do a DFT on $$E(t+\frac{\tau}{2}) E ^*(t-\frac{\tau}{2})$$ So I set up the time/frequency resolution for my DFT, using a dt value I know gives a broad enough spectrum, dt = 0.025; num = 2^14; df = ... 9 Using the undocumented IntegrationMonitor: {val, {vals}} = Reap@NIntegrate[x y^2, {x, 0, 1}, {y, 0, 1}, PrecisionGoal -> 2, Method -> "MonteCarlo", IntegrationMonitor :> ((Sow[Total@Through[#["Integral"]]]) &)] (* {0.165268, {{0.186623, 0.172189, 0.168129, 0.166339, 0.165429, 0.16988, 0.173145, 0.171675, 0.173355, 0.177199, ... 8 One can use Picard-type iteration to get the solution: Using an approximation to x'[t] (in the integral), we can integrate the ODE to obtain a new approximation. Remarkably, it converges in just two steps. My original thought was to step through the integration using the tools from tutorial/NDSolveStateData to build an interpolation of x'[t] at each step ... 8 It looks like classic catastrophic round-off error. (Look at those exponents on$e$!). {sol} = DSolve[{x'[t] == y[t]/100, y'[t] == -100 x[t] - 100 y[t] + 2020, x[0] == 0, y[0] == 20}, {x, y}, t]; Now consider y[9.] vs. y[9] and y[9.20]: y[9.] /. sol N::meprec: Internal precision limit$MaxExtraPrecision = 50. reached while evaluating ...

8

You can make a change of variable to solve the problem. Here I'll use dChange for this task: r0 = 0.5; eqn = {k[z] f[z] + f'[z] == 0, k'[z] == 0, f[-1] == Exp[2], f[1] == 1}; c = Piecewise[{{r0, z > 0}}, 1]; neweqn = dChange[eqn, f[z] == c g[z]]; {solg, solk} = NDSolveValue[neweqn, {g, k}, {z, -1, 1}]; Plot[c solg[z], {z, -1, 1}]

Only top voted, non community-wiki answers of a minimum length are eligible