Tag Info

Hot answers tagged

35

The only reason I am attempting to answer this is to perhaps get a Reversal badge. There you go... We will go slowly and this answer is the basis for what comes next. Let's start with two dimensions. You'll see why. We create a rectangular region: Needs["NDSolve`FEM`"] mesh = ToElementMesh[FullRegion[2], {{0, 5}, {0, 1}}, "MeshOrder" -> 1, ...


16

After a lengthy study (I'm using version 8) I conclude that there is a bug in Mathematica in the Integrate function when applied to a Sqrt integrand. Ok. let's go (some patience is required because of the long text) Let us define the functions corresponding to your integrals. Remark: because of the relation $1 + cos(2x) = 2 cos^2(x)$ the two forms of ...


14

The good news is that yes, there is an easy way to put your problem into NDSolve by using the new finite element functionality in v10. The bad news is that it seems the specific problem you're trying to solve is ill-posed. NDSolve can now handle internal boundaries; see e.g. the first figure under "Details" for DirichletCondition. Generating a mesh with ...


13

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


12

In this second answer I give the cause for the mismatch in the integrals, show how to remove it, and make a suggestion to improve the function Integrate[]. Simplified restatement of the problem In order to focus on the core of the problem we consider the simpler integral $\int_0^1 \sqrt{\cos (2 π k r)+1} \, dr$. It has the square root and the cosine ...


12

I'll change your proposed function because it's a constant and for such a function all sums will be equal.So: f[x_] := 1/4 x^2 di = (2 - 1)/50; intervals = Range[1, 2, di]; leftSum = Sum[f[i] di, {i, Most@intervals}]; rightSum = Sum[f[i] di, {i, Rest@intervals}]; middleSum = Sum[f[i + di/2] di, {i, Most@intervals}]; exact = Integrate[f@x, ...


11

You can write your own algorithm and use it from NDSolve. For example, for RK4: CRK4[]["Step"[rhs_, t_, h_, y_, yp_]] := Module[{k0, k1, k2, k3 }, k0 = h yp; k1 = h rhs[t + h/2, y + k0/2]; k2 = h rhs[t + h/2, y + k1/2]; k3 = h rhs[t + h/2, y + k2]; {h, (k0 + 2 k1 + 2 k2 + k3)/6}] CRK4[___]["DifferenceOrder"] := 4 CRK4[___]["StepMode"] := Fixed ...


11

Use NDSolve antiD = NDSolveValue[{f'[x] == Sqrt[1 + x^3], f[0] == 0}, f, {x, 0, 10}] Example usage: Plot[antiD[x], {x, 0, 10}] Alternatively... This works because this function can be antidifferentiated (by Mathematica). antiD = FunctionInterpolation[ Evaluate @ Integrate[Sqrt[1 + x^3], {x, 0, t}, Assumptions -> 0 < t < 10], {t, 0, ...


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


10

Update Almost ten times faster again, or about 90 times faster than the OP's way (0.069 sec v. 5.46 sec): For the second integral, we can find its derivative with respect to x and then integrate with NDSolve. The derivative of the integral has two components, one from differentiating under the integral sign dxdz1 and one from plugging in the limit of ...


9

Update The problem is subtler than my first analysis revealed. There is indeed a problem with the variable et in NIntegrate not being properly blocked. Part of the problem has to do with the extra braces in firstFuncK which has the form {{f -> InterpolatingFunction[<>]}} Somehow that leads to an evaluation of et in the integrand f[et, k] /. ...


9

In Mathematica 10, this computation may be made as follows: Clear @ r volSphere9[r_] = RegionMeasure[Ball[ConstantArray[0, 10], r]] (π^5 r^10)/120 volSphere9[1000.] 2.55016*10^30


9

The interpolation step seems to be unnecessary because the integral into which it enters can be equally well approximated as a Riemann sum. So to get really fast results you could do the following: {r, h} = Transpose[hrdata]; d = Differences[r]; Clear[s]; s[q_] := (4 Pi 0.83 )/q Total[d Rest[r Sin[q r] (h - 1)]] ans2 = Table[{i, s[i]}, {i, 0.05, 11.4, ...


8

It's numeric integration. So it has no means of "knowing" the correct result is zero. In the process error estimates will be formed and if they are larger than the estimated result, this is a problem. But of course they must be larger since the actual result is zero. The way to tame this is to specify an AccuracyGoal that is attainable using the given ...


8

There is a useful attribute, NHoldFirst whose purpose is to protect the function from exactly that. So setting: SetAttributes[a, NHoldFirst]; and then evaluating the integral works the way you want: Integrate[(a[1] + x)^2, {x, 1., 2.}] (*2.33333 + 3. a[1] + 1. a[1]^2*) The relevant example from the documentation cites "indexed" functions that are ...


8

the answer is in the links -- just to demonstrate or validate what it says: First@Last@ Reap@NIntegrate[Sin[x], {x, 0, 1}, EvaluationMonitor :> Sow[x], MaxRecursion -> 0, Method -> "GaussKronrodRule"] xi=First@Last@ Reap@NIntegrate[Sin[x], {x, 0, 1}, EvaluationMonitor :> Sow[x], MaxRecursion -> 0] ...


8

In V10 there has been added some symbolic processing of integrands containing an InterpolatingFunction. In particular if the interpolation grid divides the domain of integration into a number of subintervals, the number being at most the value of the option "MaxSubregions", the integrand will automatically be integrated over each subinterval. In V9, this is ...


8

Using the Method option with the following settings seems to work: NIntegrate[Sin[81 x]/((2^x + 1) (Sin[x])), {x, -Pi/2, Pi/2}, Method -> "LevinRule"] NIntegrate[Sin[81 x]/((2^x + 1) (Sin[x])), {x, -Pi/2, Pi/2}, Method -> "LocalAdaptive"] NIntegrate[Sin[81 x]/((2^x + 1) (Sin[x])), {x, -Pi/2, Pi/2}, Method -> ...


8

Using approximate numbers (e.g. ones with decimal points) can lead to issues with exact solvers such as Integrate. One way around, if the function can be integrated with symbolic parameters, is to use Block to block the numeric values from being substituted until after the integration is complete: Block[{x0, a, b}, Assuming[a > 0 && b > 0 ...


8

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


7

Here is a way to get an approximate symbolic expression for the integral. Some of the coefficients are approximate because at some points in the process we need the value of a definite integral at a = 1, and unfortunately those integrals have to be computed numerically. First, a lemma: Let $f(a) = \int_1^\infty g(x, a) \; dx$. Then $f(a) = f(1) + ...


7

The sum of the squares should be less than or equal to r^2 rather than r. d = 10; r = 1000; F = Piecewise[{{1, Sum[x[i]^2, {i, d}] <= r^2}}, 0]; NIntegrate[F, {x[1], -1000, 1000}, {x[2], -1000, 1000}, {x[3], -1000, 1000}, {x[4], -1000, 1000}, {x[5], -1000, 1000}, {x[6], -1000, 1000}, {x[7], -1000, 1000}, {x[8], -1000, 1000}, {x[9], -1000, 1000}, ...


7

Besides trivial observations that one cannot evaluate numerically integrals involving symbolic constants there are more interesting aspects of the problem at hand. First one should realize that a standard numeric approach is not appropriate for this kind of problems, since the integrand involves singular points (zero in the denominator) thus it is not ...


7

gauMix[x_, means_, vars_] := Total[(E^-(((x - means)^2)/(2*vars)))/Sqrt[2*Pi*vars]]/ Length[means]; means = {-7, 7}; vars = {6, 65/10}; f[x_] := gauMix[x, means, vars]; fxx = Integrate[f[x]*(x^2), {x, -Infinity, Infinity}] (* 221/4 *)


7

The reason why you get a factor of 500 of I have explained in my comment. Let's replace your g with a better behaved function: mass = 100; width = 10^-2; g[x_] := mass/width HeavisideLambda[x/width] This is a triangular peak with area 100 and basewidth 0.01. Now let's impose the desired boundary conditions of $\partial u / \partial x =0$ at the ...


6

There are two ways come to my mind to go. 1. Truncate the upper limits: Since OP has an exponential decay term like what it reads, truncating the integral limits at 2000000 should give a reasonably precise result: NIntegrate[( E^(-(x - 1000000)^2/(2*200000^2)) *E^(-(y - 1000000)^2/(2*200000^2)) *E^(-(z ...


6

Here's a way to get an exact symbolic expression for the integral. Use the notation fi[m,n,a] to denote Integrate[(Exp[-a x](x^2 - 1)^(m/2))/x^n, {x, 1, Infinity}]. Define integration by parts, where the boundary term vanishes. fi[m_, n_, a_] /; m > 0 && n > 1 := (m fi[-2 + m, -2 + n, a] - a fi[m, -1 + n, a])/(-1 + n); See how we're ...


6

This is already bult into Mathematica $Assumptions = {s > 0}; dist = NormalDistribution[m, s]; pdf[x_] = PDF[dist, x]; cf[w_] = CharacteristicFunction[dist, w]; pdf[x] == InverseFourierTransform[cf[w], w, x, FourierParameters -> {1, 1}] // Simplify True cf[w] == FourierTransform[pdf[x], x, w, FourierParameters -> {1, 1}] == ...


6

Just consider the first integral. expr = x (A Ac (m2 + m1 (-1 + x)) + V Vc (m1 + m2 - m1 x))/(8 π^2 (-mh^2 (-1 + x) + (m2^2 + m1^2 (-1 + x)) x)); denominator = Collect[Denominator[expr], x] 8 mh^2 π^2 + 8 (-m1^2 + m2^2 - mh^2) π^2 x + 8 m1^2 π^2 x^2 It has two singular points. sol = Solve[denominator == 0, x] // Simplify If the singular ...


6

In Version 10, once the points have been obtained as per user21's approach, we can tetrahedralize them directly using DelaunayMesh pf = {Cos[u], Sin[u] + Cos[v], Sin[v]}; pp = ParametricPlot3D[pf, {u, 0, 2 Pi}, {v, -Pi, Pi}] data = Reap[ParametricPlot3D[Sow[pf], {u, 0, 2 Pi}, {v, -Pi, Pi}]][[2, 1]]; pts = Cases[data, {_?NumericQ, _?NumericQ, ...



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