New answers tagged

3

Here is a method for estimating the error of integration. By sampling at two different rates, one can compare the integrals and estimate the error from the degree of precision of the integral (equals the InterpolationOrder in the cases at hand). If we interpolate $f(x)$ by a polynomial $p(x)$, then the error $\int_a^b f(x)\,dx - \int_a^b p(x) \, dx$ is ...


2

In this code we can check GaussianQuadratureWeights and FindRoot for potential errors. Let us evaluate GaussianQuadratureError[2 M, (1 + z) Sin[u[x/2 (z + 1)]], -1, 1] and we have answer for $u(x)=x^2$ -6.5402263142525195*^-105* Derivative[64][(1 + z)*Sin[(1/4)*x^2*(1 + z)^2]] Since $-1\le z\le 1, 0\le x\le 1$ we can conclude that Gauss quadrature not ...


3

Using the substitution u==a x the identity is transformed to Integrate[u/((u^2 + ab^2) Sin[u]), {u, 0, Infinity}] == Pi/(2 Sinh[ab]) The integral is singular at u==k Pi, k=1,2,... and the integration range is splitted into subintervals (similar to @yarchik answer) containing only one singularity (thanks to @flinty and @ ChipHurst comments): int[ab_?...


0

I rewrote this without the subscripts, I added NumericQ in places where it was needed, I used SetDelayed (:=) for some of your functions, and I assumed by i you meant I the imaginary number: Z = 6; n = 6; K1 = 1.55; fpn1[q_?NumericQ] := (18.25974896615874*(1 - 0.8748275119106319*I)*I)/E^(0.2115*q^2) Fpn1[q_?NumericQ] := ((4*Pi)/q)*NIntegrate[ρn[r]*Sin[q*r]*...


4

I want to focus on the numerics here. There is a way to do this integral taking infinite number of singularities into account. One can split the integral into the domains containing only one singularity at $a x_n=\pi n$, i.e. $x\in[x_n-\frac{\pi}{2a},x_n+\frac{\pi}{2a}]$ and an integral in the interval $x\in[0,\frac{\pi}{2a}]$. On the last step we use NSum ...


3

Following the @flinty's comment, we obtain Residue[x/(x^2 + b^2)/Sin[a*x], {x, π/a*n}, Assumptions -> n ∈ Integers] (*((-1)^-n n π)/(a^2 b^2 + n^2 π^2)*) Sum[%, {n, -∞, ∞}] (*0*) and 2*π*I*Residue[x/(x^2 + b^2)/Sin[a*x], {x, I*b}] (*π Csch[a b]*) Now, making use of the Jordan's lemma, we conclude that $$PV\int_{-\infty}^\infty\frac x {(x^2+b^2)\sin(...


9

We construct an NDSolve method which can pass an NIntegrate method to NIntegrate to set up an integration rule. We define a method nintegrate implements such a method. The requirements are the ODE is of the form y'[x] == f[x], and the NIntegrate method returns an interpolatory rule. Example: foo = NDSolveValue[{y'[x] == Sin[x^2], y[0] == 0}, y, {x, 0, 15}...


2

Too long for a comment. First, your F[w_,s_]:= (1/(s + (I*a*w)/mx + ((\[Sigma]^2 + a^2)*w^2)/(2*mx) + balpha/mx^2*s^(alpha - 1)*(I*w*a + 1/2*(\[Sigma]^2 + a^2)*w^2))) is a complex valued function. I don't see any reason why its inverse Fourier transform should be a real valued function. Second, up to the definition used by Mathematica, its inverse Fourier ...


2

I now realize that the question can be solved analytically, for the most part, although not with ImplicitRegion. The constraints embodied in R1 and R2 can be solved to obtain θ in terms of v and parameters. r1c1 = Reduce[θ (v + r) > r && v > 0, θ] // Last (* θ > 1/(1 + 2 v) *) r1c2 = Reduce[2 θ v s + (1 - θ) r > λ1 && v > 0 &...


1

The integration NIntegrate[g[NIntegrate[f[x],{x,0,y}]],{y,0,1}]fails because the inner NIntegrate doesn't know about y to be numeric! Try int[y_?NumericQ]:=NIntegrate[f[x],{x,0,y}] NIntegrate[int[y],{y,0,1}]


4

This interesting problem can be solved numerically by computing InterpolationFunctions for the two sums of integrals in the last two lines of code in the Question. λ2 = -(1/2) + r + 1/(1 + r) - r^2 Log[1 + 1/r]; k = 0.5; r = 0.5; t = Flatten[Table[ R1 = ImplicitRegion[θ (v + r) > r && 2 θ v s + (1 - θ) r > λ1, {{θ, 0, 1}, {v, 0, 1}}]; R2 = ...


0

You need to estimate the value of $\rho_0$ using NIntegrate: Subscript[a, p] = 0.428; Subscript[α, p] = 1.04; R = 2.45; Subscript[ρ, 0] = 1/NIntegrate[1/( 1 + (E^(2.336448598130841` (-2.45` + r)) + E^(-2.336448598130841` (2.45` + r))) (0.5` + 0.08329862557267803` r^2)^1.04`), {r, 0, ∞}] Subscript[ρ, p1][r_] = Subscript[ρ, 0]/((E^((r - R)/...


0

Here's a variation of bbgodfrey's answer that doesn't require quieting: g[a_,b_]:=Activate @ Block[{NIntegrate=Inactive[NIntegrate]},Derivative[1,0][f][a,b]] Check: g[1, 1] 1.80525 The trick is that Mathematica knows how to take derivatives of inactive integrals, so temporarily inactivating NIntegrate avoids having NIntegrate trying to integrate a ...


3

J = ImplicitRegion[ Abs[xa] < 1 && Abs[xa - xb] < 1 && Abs[-L + xb] < 1 && Abs[xa - xc] < 1 && Abs[-L + xc] < 1 && L - 2 <= xa <= 1 && L - 1 <= xb <= 2 && L - 1 <= xc <= 2, {xa, xb, xc}]; Assuming[2 < L < 3, RegionMeasure[J]] (* 1/3 (27 - 27 L + 9 L^2 - L^...


3

reg = ImplicitRegion[{Abs[xa] < 1, Abs[xa - xb] < 1, Abs[-L + xb] < 1, Abs[xa - xc] < 1, Abs[-L + xc] < 1, L - 2 <= xa <= 1, L - 1 <= xb <= 2, L - 1 <= xc <= 2}, {xa, xb, xc}]; Assuming[2 < L < 3, Integrate[1, Element[{xa, xb, xc}, reg]]] // Simplify // AbsoluteTiming (*{0.142462, -(1/3) (-3 + L)^3}*)


1

If you you only need a numerical result for given L try int[L_?NumericQ] :=NIntegrate[ Boole[Abs[xa] < 1] Boole[Abs[xa - xb] < 1] Boole[ Abs[-L + xb] < 1] Boole[Abs[xa - xc] < 1] Boole[ Abs[-L + xc] < 1], {xa, L - 2, 1}, {xb, L - 1, 2}, {xc, L - 1,2} ] int[2.5] // AbsoluteTiming (*{0.0222419, 0.0416667}*)


1

This is not a complete answer but an extended comment. ImplicitRegion does not like the usage of function Tmaxc12, so we can construct it inline: region[U_?NumericQ] := ImplicitRegion[ (0.5*MC12 (vx^2. + vy^2. + vz^2.) + (1 - Cos[theta])*(Sqrt[Te[U]*(Te[U] + 2 m*c^2)/c^2] + MC12*vz)* Sqrt[Te[U]*(Te[U] + 2 m*c^2)/c^2]/MC12 - Sqrt[...


7

First of all, there do not seem to be any singularities in the integration region: $Assumptions = And[-Pi<=x<=Pi,Pi<=y<=Pi,Pi<=z<=Pi,-5<=r<=2]; den = 478 + (96*I)*br + 32*br^2 + 96*r - (64*I)*br*r - 32*r^2 + 11*Cos[2*x] - 264*Cos[y] - (48*I)*br*Cos[y] - 48*r*Cos[y] + 11*Cos[2*y] - 336*Cos[z] + 144*Cos[y]*Cos[z] + 12*Cos[...


3

The root cause of the difficulties encountered by the OP is that MeijerG, as implemented in Mathematica, often is not accurate when evaluated at machine precision. For instance, with the parameters, params = {p -> 1, q1 -> 1, q2 -> 1, q13 -> 1, qi1 -> 1, mu -> 3/2, L -> 3, gammak1 -> 1, BI -> 9, C1 -> 4, Bsd -> 4, M -&...


1

On such a high dimensional integral, the default rule is the Monte Carlo rule. You can increase the number of points. I also increased the PrecisionGoal, so that the error estimate will be reported. NIntegrate[..., Method -> {"MonteCarloRule", "Points" -> 10^6}, PrecisionGoal -> 6] NIntegrate::maxp: The integral failed to ...


3

The symbolic processing misses the obvious way to compute these integrals, namely, Method -> "Trapezoidal", about 100 times faster than "SymbolicProcessing" -> 0: NIntegrate[ Exp[2 I s] Exp[2 I t] ((Cos[s] - Cos[t])^2 + (Sin[s] - Sin[t])^2) + 1, {s, 0, 2 π}, {t, 0, 2 π}, Method -> "Trapezoidal"] // RepeatedTiming ...


4

The first NIntegrate spends some time doing symbolic pre-processing of the integrand. You can turn that off and the integrals complete in about the same time: NIntegrate[ Exp[2 I s] Exp[2 I t] ((Cos[s] - Cos[t])^2 + (Sin[s] - Sin[t])^2) + 1, {s, 0, 2 π}, {t, 0, 2 π}, Method -> {Automatic, "SymbolicProcessing" -> 0}] // Timing (* result: {...


1

The trick seems to be to FunctionExpand the MeijerG function before integration: a[θ_] = MeijerG[{{-2, -17/2, -7/2}, {}}, {{0}, {-3/2}}, θ] // FunctionExpand // FullSimplify (* (Sqrt[π] Sqrt[θ] (128 + θ (4416 + θ (50048 + θ (213120 - 7 θ (-36608 + 3 θ (3168 + 5 θ (-432 + θ (400 + 9 θ (-40 + θ (18 + 11 θ)))))))))) - 32 E^(1/θ) π (4 + 7 θ (20 + θ (233 ...


6

We can integrate in 3 steps: Integrate[(yp/(b + yp^2)^(3/2)), {yp, 0, Infinity}, Assumptions -> b > 0]*(x - xp) /. {b -> g (x - xp)^2} //Simplify Out[]: (x - xp)/Sqrt[g (x - xp)^2] So we have intyp=1/Sqrt[g] as results and it means that Q[g,x] not depends on x. Next step: Integrate[(Exp[-ypp^2/(8 T)])*(-1 + ypp^2/(8 T)) (ypp/(g*(xp - ...


1

To reduce the number of integration in NIntegrate seems reasonable. The effects are somehow dependent on the choices of options for NIntegrate. Choices are Values for lower integral bounds larger than zero. Values replacing the infinite upper bound of the integral to a meaningful numerical integration value. The default method is GlobalAdaptive. This can be ...


0

Clear["Global`*"] Subscript[f, pp1][q_] = (6.254736279890945*(1 - 1.4511668475476842))/ E^(0.2115*q^2); Subscript[F, pp1][q_?NumericQ] := ((4*Pi)/q)* NIntegrate[Subscript[ρ, p][r]*Sin[q*r]*r, {r, 0, Infinity}]; Subscript[ρ, p][r_] = 0.013132593248303927/(1 + \ (E^(1.7543859649122808*(-2.380427976610103 - r)) + E^(1....


3

The integrand can be simplified to a significantly more compact form. First, the floats in f[x,y,z] can be replaced by rationals, which are easier for Mathematica to deal with: f[x_,y_,z_] = {{3/2 - (3*(Cos[x] + Cos[y]))/4, (-I)*Sin[x] - Sin[y], 7/2 - 3/(2*E^(I*z)) - (3*(Cos[x] + Cos[y]))/2, 0}, {I*Sin[x] - Sin[y], 3/2 - (3*(Cos[x] + Cos[y]))/4, 0, ...


1

Hint. n = 3; tmax = 5; E1 = x*y^(1/2)*Exp[-2*x]*Exp[-2*y]*(1 + fg (Exp[-x]*Exp[-y]*x^2))^(-1/2) F = Integrate[Normal[Series[E1, {fg, 0, n}]], {x, 0, Infinity}, {y, 0, Infinity}] /. {fg -> 1/f[z]^2/g[z]} F1 = F - (1/f[z])*(D[f[z], z])^2 G1 = F - (1/g[z])*(D[g[z], z])^2 sol = NDSolve[{D[f[z], {z, 2}] == F1, D[g[z], {z, 2}] == G1, f[0] == 1, f'[0] == 0, g[0] ...


9

The line Table[N[PHI30EQ[u], 30], {u, 0, 0.1, 0.001}] doesn't do what you think it does. You're asking for 30 digits of precision, but you supply u as a machine number. If you mix arbitrary precision and machine precision like that, you'll get machine precision answers. I suspect you instead want is: Table[N[PHI30EQ[u], 30], {u, 0, Rationalize[0.1], ...


1

While Erf can be difficult to work with, a simple, straightforward approach seems to work well: ClearAll[Nfunc]; Nfunc[x_] := E^(-x^2/2)/ Sqrt[2*Pi]*(1/ 2 Erf[(x - 0.256048)/Sqrt[2*1.6^2 + 0.231313^2]/Sqrt[2]] - 1/2 Erf[-Infinity/Sqrt[2]]); (* why not -1 instead of Erf[]? *) nume = NIntegrate[x*Nfunc[x], {x, -Infinity, Infinity}] (* 0....


1

Clear["Global`*"] m2 = 1/2; ℏ = 1; w = 1/2; \[ScriptCapitalO]2 = -ℏ^2/(2 m2) Laplacian[u[x, y], {x, y}] + 1/2 m2 w^2 (x^2 + y^2) u[x, y]; {vals, funs} = NDEigensystem[{\[ScriptCapitalO]2, DirichletCondition[u[x, y] == 0, True]}, u[x, y], {x, -10, 10}, {y, -10, 10}, 28, Method -> {"PDEDiscretization" -> {"...


4

Have a look at the documentation. This is from the ref page of NDEigensystem {vals, funs} = NDEigensystem[-Laplacian[u[x], {x}], u[x], {x, 0, \[Pi]}, 4] NIntegrate[#^2, {x, 0, \[Pi]}] & /@ funs (* {1., 1., 0.999995, 1.} *) Note that the argument u[x] to NDEigensystem tells NDEigensystem that the resulting interpolating functions will also have the ...


0

Alright, I seem to have a solution that works for now, but I do planning on using a loop and this might get a little hard to incorporate into the same, so if anyone has any better ideas please do let me know. I explicitly defined an Integrand as a function of x and y and was able to integrate the functions after that Integrand[x_,y_] = \[Psi]1[x,y] x \[Psi]2[...


3

Note that Last[intFunc[t]] equals t, so the two confusing integrals are just integrating t from 0 to 1, which is why one gets 1/2. What is needed is Indexed[intFunc[t], -1], which extracts a part only when its first argument is a vector. NIntegrate[Indexed[intFunc[t], -1], {t, 0, 1}] (* 0.00378552 *)


9

When I run your code I get a FindRoot warning message: Which makes me suspicious of the result quality. If we assume the result is correct we can speed up the integration by using the FEM for that too. We create a boundary element mesh of the foil: bmeshFoil = ToBoundaryMesh["Coordinates" -> coords[[5 ;; nn]], "BoundaryElements"...


5

Here is a partial non-NIntegrate answer that still needs work but might give you some ideas on how to proceed. I extended the domain so that it would be easier for me to pick line segments related to the airfoil. x1 = -2; x2 = 3; y1 = -1.5; y2 = 1.5;(*domain dimensions*) Then I followed this example from the documentation to grab normals at line segment ...


4

It's a question of numerical precision. Set ky and ky to zero and plot the integrand ii, enlarge working precision and rationalize variable in Plot. ii[p_, u_] = 96*Cos[p/2]*(1 + I*kx*Cos[p] + I*ky*Sin[p])^(-7/2)*1/(6!* Sqrt[Pi])*(((Sqrt[ 2*u]/(Cos[ p/2]*(1 + I*kx*Cos[p] + I*ky*Sin[p])^(1/2)))^6 + 15/2*(Sqrt[ 2*u]/(Cos[...


3

The following problem is hidden by NIntegrate, but which is one of the common sources of the NIntegrate::izero message you get from the OP's code. Exp[-I*k*s]*s^(1/2 + 2)*Exp[-s]*Cos[f/2]* Exp[-Sqrt[8*s*u]/Cos[f/2]] /. {k -> 0, u -> 0.1} /. {f -> 1., s -> 1000.} General::munfl: Exp[-1032.23] is too small to represent as a normalized machine ...


3

Transform Piecewise to UnitStep also helps: S2 = NDSolve[{h''[t] == d[h[t]] // Simplify`PWToUnitStep, h'[-120] == 1, h[-120] == -120}, h, {t, -120, 300}] Plot[h[t] /. S2, {t, -120, 120}, AxesLabel -> {"t", "h[t]"}, PlotLabel -> "tf=300"]


2

Ib[λ_] := (2 0.59552 10^8)/(λ^5 (Exp[14387.75/(λ 385)] - 1)) Plot[Ib[λ], {λ, 0.28, 25}] Integrate[Ib[λ], {λ, 0, ∞}] 396.557 data and dataInt information is missing.


2

This is the typical case for AccuracyGoal (absolute error) and PrecisionGoal (relative error). Define the integrand Iint[q_, p_] = Exp[-I n q] Exp[I n1 q] (rx^2 (Cos[q] - Cos[p])^2 + ry^2 (Sin[q] - Sin[p])^2)^m Sqrt[rx^2 Sin[q]^2 + ry^2 Cos[q]^2] Sqrt[rx^2 Sin[p]^2 +ry^2 Cos[p]^2] (rx^2 Sin[q] Sin[p] + ry^2 Cos[q] Cos[p]); Now do the integrations ...


Top 50 recent answers are included