Tag Info

Hot answers tagged

13

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] := ...


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


8

I think there's a bug in the internal function NDSolve`SPRKDump`CheckSeparability that leads NDSolve to conclude that the system is not separable. I think you should report it and see if WRI can verify it (they would probably appreciate a link to this Q&A). It's a fair amount of work to track it down, and there is a lot of nearly unreadable stuff to ...


6

This problem can be simplified enormously by performing the integral over ϕ symbolically. To do so, first consider the ϕ-dependent exponential in the integrand. Exp[I k (x Sin[θ] Cos[ϕ] + y Sin[θ] Sin[ϕ] + z Cos[θ])] It can be transformed and simplified as follows. expt = Exp[I k ( x Sin[θ] Cos[ϕ] + y Sin[θ] Sin[ϕ] + z Cos[θ])] /. {x -> r ...


5

With two minor changes to the code bm = ToBoundaryMesh["Coordinates" -> coord, "BoundaryElements" -> {LineElement[incidents]}]; nds = NDSolveValue[{Inactive[Div][sigma*Inactive[Grad][u[x, y], {x, y}], {x, y}] == 0, DirichletCondition[u[x, y] == 0, y == 0], DirichletCondition[u[x, y] == 1, y == 1]}, u, {x, y} ∈ mesh]; the ...


5

The the function y is the integral of a logarithmic singularity, so it is relatively easy to understand. The function JacobiCN[s, 7/10] is relatively flat near s == EllipticK[7/10], so that the integral of its reciprocal can be approximated by a logarithm: y[s]/JacobiCN[s, 7/10] ~= y0 / (jp (s - EllipticK[7/10])) where y0 ~= y[s] at s == EllipticK[7/10] ...


5

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


5

Since your integrand does not approach zero but a finite positive number, Limit[Exp[-16.136 (1 - Exp[-0.012*t])], t -> Infinity] (* 9.82255*10^-8 *) the integral over {t, 0, Infinity} does not converge. By the way, the error in the NIntegrate[integrand, {t, 0, 1000}] should be about 10^-7, which seems better than R. In fact, the precision seems ...


4

When you try to NIntegrate your expression, the error messages include: "suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small". Your expression does not have a singularity, its integral is manifestly not zero, it is not oscillatory, so it must be a numerical precision problem. ...


4

To be honest, the code you found isn't a good example of coding in Mathematica. I think the following 3 lines are enough for you: ListStreamPlot[Transpose[{u, v}, {3, 1, 2}], DataRange -> {{0, 1}, {0, 1}}] curl = Most /@ Differences@v - Most@Differences[u, {0, 1}]; ListContourPlot@LowpassFilter[curl, 1] You may also want to try ListDensityPlot: ...


4

Compile f and use a memo-ized version of it Since it seems like NIntegrate decides to symbolically evaluate its argument first, I thought I'd force it not to by compiling the function f. This seems to make a significant difference: Clear[f, f1, g] g[x_] = Nest[f[x] + 1./# &, f[x], 500]; f1 = Compile[{x}, Sum[1/100 Erfc[-(x^2/k)], {k, 100}]]; ...


4

Integrate can be used to handle the DiracDelta, Integrate[Integrand[p, Q2, ν, θ], {p, 0, Infinity}, Assumptions -> p3zero2[Q2, ν, θ] ∈ Reals && p3zero1[Q2, ν, θ] ∈ Reals] (* ((HeavisideTheta[p3zero1[Q2, ν, θ]]*p3zero1[Q2, ν, θ]^2)/ ((M - 2*E3[p3zero1[Q2, ν, θ]])^2*E3[p3zero1[Q2, ν, θ]]* E4[p3zero1[Q2, ν, θ], Q2, ν, θ]*Derivative[1, 0, ...


4

While the unknown g might or might not be amenable to this approach, it's really fast on g == 1: Expectation[ 1, {x1, x2, x3, x4} \[Distributed] MultinormalDistribution[muvec, sigmat]] // RepeatedTiming (23000 times faster than NIntegrate.) Then there's also NExpectation to try, too. Ever since I came across Guess who it is's use of the ...


3

Because the integrand is highly localized at the origin; e.g., Log[10, npdf[{6, 0, 0, 0}]] // N (* -9.16177 *) limiting the range of integration reduces run time by a factor of three and moderately improves accuracy. muvec = ConstantArray[1/10, 4]; sigmat = IdentityMatrix[4]; npdf[x_] := PDF[MultinormalDistribution[muvec, sigmat], x]; ...


3

The error messages come from Integrate/NIntegrate calls on functions with symbolic (non-numerical) parameters a, b, c. These calls come before FindRoot substitutes numerical values for these parameters. (Taking N of Integrate here is effectively the same as calling NIntegrate.) In this case the best approach it seems to me is to do the integral first, ...


3

I do not know what range of $s$ you had in mind, so I arbitrarily chose one below. You should be able to change it at will. This is an example in hopes that you will study this further and get to the point where you can come up with such code yourself. Please take a look at the excellent links provided by Arnoud as well to get you started. Clear[phi] ...


3

Update 3 As ilian correctly states in the comments, without NIntegrate we get the same picture as shown in the "Update 2" section. It proves that memory grow reflects temporary memory allocation by Table, and ClearSystemCache[] does completely eliminate the memory leak. We can also switch off the caching completely: SetSystemOptions[ "CacheOptions" ...


2

The reason for this error is that NIntegrate uses fixed precision when computing the integration ranges, while EllipticK needs to raise the precision internally to obtain a good result. N[EllipticK[7/10], 20] (* 2.0753631352924691439 *) Block[{$MinPrecision = $MaxPrecision = 20}, N[EllipticK[7/10], 20]] (* Divide::infy: Infinite expression ...


2

Clear[integrand] With your definition that uses an unnecessary list and then takes a part integrand[t_?NumericQ] := {1000/Sqrt[t^2 - tMin0^2], 0}; integrand[t][[1]] t This occurs because integrand cannot evaluate with a symbolic input so the first part of the unevaluated expression is the function's argument. Clear[integrand] tMin0 = 1/10; tMax0 ...


2

It is quite tricky! Piecewise[ ] functions work only with "EquationSimplification" -> "Residual"... I'll try to dig up why simplesys = {r'[t] == Piecewise[{{1, 0 <= t <= 10}, {0, 10 <= t <= 20}}, 0], r[0] == 0}; state = First@ NDSolve`ProcessEquations[simplesys, r, t, Method -> {"EquationSimplification" ...


2

Here you have a way without using neither UnitStep nor PieceWise that improves the performance by 75% wrt your code. It computes 200 functions in a very reasonable time for your toy example. The main trick is to use a numeric (black box) function to be able to take Part[... ] inside it. n = 200; af = Array[f, n]; taf[t_] := Through[af[t]] bb[i_?IntegerQ, ...


2

ClearAll[δ, f]; f := 2 Δ^2 Sin[(A + B)/2]^2 Re@NIntegrate[ E^(-I t (2 ϵ + F + G)/Λ) E^(-Z/Λ (t ArcTan@t- Log[1 + t^2]/2)), {t, 0, ∞}] Plot[f, {δ, 0, 100}]


2

The following part is incorrect, and redundant anyway: \[Psi]= \[Integral]u dx= \[Integral]v dy \[Omega]= \[PartialD]v/dx-\[PartialD]u/dy Otherwise, the code you copied has only one major problem, in the definition of n. If you try to execute it, you get Do::iterb: Iterator {j,2,n-1} does not have appropriate bounds. >> which suggests where the ...


2

You can integrate as follows. Integrate[BesselI[-nu, k*x]/x, {x, r1, r}, Assumptions :> {k \[Element] Complexes, r1 \[Element] Reals, r \[Element] Reals, nu \[Element] Reals, r1 > 0, r > r1} The result is a complicated expression in terms of Gamma and HypergeometricPFQRegularized functions. Nevertheless, it can be ...


1

The first thing I did was to rationalizing all calculations, starting with the defintion of σ and minroot. This stops the Solve::ratnz messages. I also made some other improvements to minroot. σ = 6/10; minroot[gg_?NumericQ, bb_?NumericQ] := Module[{b, g, rts, r}, b = Rationalize[bb, 0]; g = Rationalize[gg, 0]; rts = r /. Solve[ ...


1

At the moment the LUM function as written can't be used because it depends upon a number of global variables (h, me, k, Tcmd ...) that you have not supplied. However, in order to produce a list of {x2, LUM[x1,x2]} pairs this can be done as follows. Scalar Assume x2List is the list of values that you want to evaluate and you want to evaluate LUM at a ...


1

About 30X faster integration, no errors... Clear[x, x1, x2, x3, x4] muvec = ConstantArray[1/10, 4]; sigmat = IdentityMatrix[4]; pdf = PDF[MultinormalDistribution[muvec, sigmat], {x1, x2, x3, x4}] // FullSimplify Integrate[ 1*pdf, {x1, -Infinity, Infinity}, {x2, -Infinity, Infinity}, {x3, -Infinity, Infinity}, {x4, -Infinity, Infinity}] // ...


1

This is your function: phi[s_, H_, B_] := NIntegrate[(1 + B Sin[2 H t])/Sqrt[ 1 + B^2 + 2 B Sin[2 H t]], {t, 0, s}] 1/(2 H) Sqrt[ 1 + B^2 + 2 B Sin[2 H s]] But if you want to plot it you have to give the value for s. Plot[{phi[s, 0.5, 1], phi[s, 1, 1], phi[s, 1.5, 1], phi[s, 2, 1]}, {s, 0, 16}]


1

[Background: I added this answer because this question was cited in What are the most common pitfalls awaiting new users? as an example of when NumericQ was required. I took the integrand to be a simple example of a more general problem with NIntegrate, but on reading the comments, it seems to be narrowly the exact integrand of interest to the OP. I ...



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