Hot answers tagged

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


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

The following is not particularly fast and could be more accurate but does make progress toward the goals set in the Question. To begin, consider the analytical properties of f, as defined in the Question. By observation, it has branch points at {I Sqrt[6/5] ξ, I Sqrt[2/5] ξ} and their conjugates. Poles are obtained by p /. ...


7

You can first discretize the region as follows: reg = DiscretizeRegion[ ImplicitRegion[0 < Sinh[u]/Cosh[v] < 1 && 0 < Sinh[v]/Cosh[u] < 1, {{u, 0, 3}, {v, 0, 3}}]]; Now NIntegrate gives us better values: NIntegrate[1, Element[{u, v}, reg]] 1.23371 I don't know why putting Infinity as the limits ...


5

Some of what is going on is that the region is numericized to the extent that the precision of the numbers in R is set to 7.. Then NIintegrate discretizes the region with dR = DiscretizeRegion[ ImplicitRegion[ 0 < Sech[v] Sinh[u] < 1.000000 && 0 < Sech[u] Sinh[v] < 1.000000 && u >= 0 && v >= 0, {u, v}], ...


4

It appears that by far most time is spent on symbolic precomputations during NIntegrate. You can avoid this by calling the NIntegratecommand with the option Method -> {Automatic, "SymbolicProcessing" -> 0}. fs[25., 25., 0.01] // AbsoluteTiming fsold[25., 25., 0.01] // AbsoluteTiming (* {0., 695.762} *) (* {0., 695.762} *)


4

Here is a simple example for the default method (LSODA) that shows a couple of the issues related to the situation in the question: $$y' = \exp(-10 \, x^2), \quad y(-20) = 1, \quad -20 \le x \le 20$$ The exact solution is a scaled error function (Erf), which has a sigmoid graph going from $y(-20)=1$ to $y(20) \approx 1.5605$. For $|x|$ large, $y'$ is ...


4

I am not sure I understand the question 100% but here is what I think you are looking for: {vals3DL0, funs3DL0} = NDEigensystem[-(R''[r]/2) - (R'[r]/r) - (R[r]/r), R[r], {r, 0, 200}, 3, Method -> {"Eigensystem" -> {"Arnoldi", "Criteria" -> "RealPart"}, "SpatialDiscretization" -> {"FiniteElement", {"MeshOptions" -> \ ...


4

Changing parameter values during integration works better with DiscreteVariables. But I think the problem with OP's code, in the question and the OP's answer, has more to do with Mathematica numerics. My solution Clear[bind]; zdot = 1/2 (1 - z[t]); ydot = 1/20*y[t] + z[t] - x[t]; xdotbind = D[Solve[-ydot - zdot == 0, x[t]][[1, 1, 2]], t] /. {y'[t] -> ...


3

The trace shows what happens, but not why, which reasons remain obscure to me. NIntegrate internally evaluates the integrand, so ultimately the outside integral is computed on the integrand 0.5, which is the value of both a1 and a2[x]. The slower way has essentially three evaluations of NIntegrate. The faster way has only two, but it also has several extra ...


3

TestCDF2 has a number of singularities, each of which might cause a convergence warning message: (* Data prepared with SeedRandom[0] for reproducibility *) Plot[TestCDF2, {x, Min[Flatten[TestData1]], Max[Flatten[TestData1]]}, PlotPoints -> 200] There is probably a limit in both time and number of singularities that NIntegrate will detect on its ...


3

Your system is singular and cannot be fixed with Adams's method (see also J.M.'s comment). You can solve it with DSolve. nsol = NDSolveValue[{y'[x] == y[x]*y[x]*Exp[x] - 2*y[x], y[0] == 2}, y, {x, 0, 10}] At x == 0.6931469688260267`, step size is effectively zero; \ singularity or stiff system suspected. {x1, x2} = nsol["Domain"][[1]] (* {0., ...


3

TechSupport acknowledged and proposed a simple workaround by putting inert expression, e.g. empty string, inside the { }: NDSolve[{x'[t] == x[t], x[0] == 1, WhenEvent[Mod[t, 1] == 0, {""}]}, x, {t, 0, 3}]


3

The problem is that you are not using large enough precision goal for an oscillatory integrand with very small absolute values over the integration range. Look at the log plot -- there are more oscillations than the ones visible with Plot: The remedy is to use higher precision goal. Try this: NIntegrate[ Re[(x^4 Sin[x b])/(9 + 4 d^2 x^2)^6 1/ 2 ...


3

Short answer: One workaround is to use Method -> "LevinRule" instead. Long answer: As mentioned by Xavier in a comment above, changing BesselJ[0, x] to BesselJ[0, Re[x]] resolves the issue: NIntegrate[{1, -1} BesselJ[0, Re@x], {x, 0, ∞}, WorkingPrecision -> 32, Method -> "ExtrapolatingOscillatory"] Precision /@ % ...


3

Below is a workaround for the simple case. The OP can say whether it works more general. I haven't quite tracked down yet why the system is set up incorrectly with the default Method -> {"EquationSimplification" -> "Solve"} and with Method -> {"EquationSimplification" -> "Residual"}. But it works in this case with Method -> ...


2

As I noted in a comment above, it is quite possible for the solution of a nonlinear PDE to become singular at finite t, and that appears to be occurring here. And, as noted by Dr. Wolfgang Hintze, the right side of the PDE, when integrated over {x, -6, 6} is zero. So, the integral of f must be a constant, and indeed it is. For the parameters given in the ...


2

I'd suggest a simpler approach that seems to work just fine: Use arbitrary precision input if you can. In your case, instead of d -> 0.001/5, use d -> 5/1000 integrand = Re[(x^4 Sin[x b])/(9 + 4 d^2 x^2)^6 1/2 E^(-(1/2) x^2 - I x c - 1/2 a^2) (Erfc[(I x - c)/Sqrt[2]] + E^(2 I x c) Erfc[(I x + c)/Sqrt[2]]) /. {b -> 5, ...


2

As noted in the question, the first code segment does run to completion. However, NDSolve also issues two warning messages. NDSolve::ibcinc: Warning: boundary and initial conditions are inconsistent. >> NDSolve::eerr: Warning: scaled local spatial error estimate of 37.88105237757757at t = 100. in the direction of independent variable x is much ...


2

I should confess that I am still slightly perplexed by what you actually want to do, but I'll proceed with the interpretation that you're starting from a Cesàro equation implicitly defined by a second-order ODE (which, BTW, can be solved in terms of elliptic functions, but let's pretend we don't know that), and you want to play with the boundary conditions. ...


2

So, after playing around with this problem all day, I finally stumbled upon a version of the code that gives me my desired result, although I'm not at all clear on why it works where other versions fail. First, here is the working version: zdot=.5*(1-z[t]); ydot=.05*y[t]+z[t]-x[t]; ...


2

the culprit here is AccuracyGoal: a = 3; c = 6; d = 0.00033; b = 2; NIntegrate[ x^3 (SphericalBesselJ[0, b x] + SphericalBesselJ[2, b x])/(4 d^2 x^2 + 9)^6 1/ 2 E^(-(1/2) x^2 - I x c - 1/2 a^2) (1 + E^(2 I x c) - I E^(2 I x c) Erfi[(x - I c)/Sqrt[2]] - I Erfi[(x + I c)/Sqrt[2]]), {x, 0, 8}, MaxRecursion -> 22, AccuracyGoal ...


2

I was just inspired by @kglr to illustrate: f[x_] := 1/1000 x^4 - 280/1000 x^2 + 25; sim[a_, b_, n_] := With[{u = RandomVariate[UniformDistribution[{a, b}], n]}, {#, f@#} & /@ u]; Manipulate[ Module[{res = sim[a, b, n], mn = Integrate[f[x], {x, a, b}]/(b - a)}, Show[Plot[f[x], {x, a, b}], ListPlot[res, PlotStyle -> Red], GridLines ...


2

You can also use the built-in functions Expectation or NExpectation taking the function argument x to be uniformly distributed over the interval {a,b}. f[x_] := 1/1000 x^4 - 280/1000 x^2 + 25 Expectation[f[x], Distributed[x, UniformDistribution[{a, b}]]] NExpectation[f[x], Distributed[x, UniformDistribution[{-12, 12}]]] 15.7072 ...


2

See also Average Function and Average Value of a Function. Say you have messured values and express them with the Function f(x): f[x_] := 1/1000 x^4 - 280/1000 x^2 + 25 You can Plot that Function: Plot[f[x], {x, -15, 15}] And you like to find the average value as of -12 ... 12 a = -12; b = 12; solNI = NIntegrate[f[x], {x, a, b}]/(b - a) ...


2

There are two different methods for solving the OP"s problem : The Method of Lines with the option "SpatialDiscretization" -> {"TensorProductGrid"... The Method of Lines with the option "SpatialDiscretization" -> {"FiniteElement"}. This solution is the Mathematica 10 implementation of the Finite Element Method for transcient PDEs. In both cases ...


1

NDEigensystem was added in version 10.2 (or was it 10.3?) but version 10.0 is not going to work unless you use this.


1

α = 3; δ = 2/α; int = Assuming[{A > 0, km > 0, r > 0, λ > 0}, Integrate[(1 - Exp[-x*h*r^(-1/δ)])*km*λ*π, {r, A, Infinity}, GenerateConditions -> False]]; The condition suppressed by use of GenerateConditions -> False is Re[h x] > 0, so include assumption that x > 0 and FullSimplify (this is quite slow) the expression for ...



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