# Tag Info

2

As an illustration of my initial comment, let's look at various starting points, following @Andrew's method: sols = Table[ NDSolveValue[(2 + Sqrt[2] + s^2 - (-2 + Sqrt[2]) s^4 + s^6)/(1 + s^2)^2 yy[s] - (1 - s^2)^2 yy''[s] == 0 && yy[start] == yy'[start] (start + 1) && yy'[start] == 1, yy[s], {s, start, 1}], {start, -1 + ...

1

As a workaround it is possible to step up a bit from the left end: leftend = -1 + 10^-8; NDSolve[(2 + Sqrt[2] + s^2 - (-2 + Sqrt[2]) s^4 + s^6)/(1 + s^2)^2 yy[s] - (1 - s^2)^2 yy''[s] == 0 && yy[leftend] == 0 && yy'[leftend] == 1, yy[s], {s, leftend, 1}][[1]]

5

The problem is that Mathematica prematurely threads r[t] - p not knowing r[t] is actually in $\mathbb{R}^2$ In[]:= r[t]-{0,0} Out[]= {r[t],r[t]} Which is not what you want. A quick fix for these types of issues is to create a function that only evaluates for numerical values (Changed to NDSolve since I only have v8): dummy[r_?(VectorQ[#, NumericQ] ...

7

I think your initial condition is singular. In order to solve the ODEs, consider $(x(t), y(t))$ as a planar curve, we may try changing the parameter $t$ to the arc length parameter $s$: \left\{\begin{split} \frac{\mathrm{d}x}{\mathrm{d}s}=\frac{x'(t)}{\sqrt{x'(t)^2+y'(t)^2}}\\ \frac{\mathrm{d}y}{\mathrm{d}s}=\frac{y'(t)}{\sqrt{x'(t)^2+y'(t)^2}}\\ ...

5

- update - @whuber gives insightful comment and I agree I should mention here the following. Below we consider a perturbed form (solvable exactly in Bessel functions) of original equation up to linear term. This allows us to understand behavior of the system around t ~ 0. - original - The following shows that under some general assumptions there is no ...

2

Your functions for y and x are single values and not lists, and thus doing a ParametricPlot for the sum of these is not going to work. However, to solve the equations and get a plot you can run the following: alpha = 0.5 beta = 1 solf[b_, c_, tp_] := Module[{sol}, sol = NDSolve[{x'[t] == -x[t]*a[t] - 3*x[t] + alpha*y[t]^2 + 1/2*beta*(r1[t]^2 - ...

0

Could be interpreted as a nonlinear least squares problem. pts = {{70.26, 45.78}, {71.04, 46.32}, {37.23, 24.67}, {37.91, 28.78}}; ndsoln[{k1_, k2_, k3_}, x0_] := x[720] /. NDSolve[{x'[t] == -k1 x[t]^2 + k2 x[t] y[t], y'[t] == k1 x[t]^2 - k2 x[t] y[t] - k3 y[t], x[0] == x0, y[0] == 0}, {x, y}, {t, 0, 720}][[1]] ssfun[{k1_?NumericQ, ...

13

Some frames from my version of the animation: Here's the code I used: orbit[posStart_?VectorQ, derStart_?VectorQ] := Block[{c = -Rationalize[6.672*^-11*7*^17], x, y, z, t}, {x, y, z} /. First @ NDSolve[ Join[Thread[{x''[t], y''[t], z''[t]} == c {x[t], y[t], z[t]}/Norm[{x[t], y[t], z[t]}]^3], ...

2

This is too long for a comment, but it isn't an answer. Perhaps you would like to consider a more compact way to write down your equations: m = -(6.672*10^-11) (7*10^17) ; st = {{1000, 1000, 1000, 0, -100, 0}, {500, -1000, -1000, -110, 100, 0}, {0, 100, 500, 350, -100, 0}}; r = {x @ t, y @ t, z @ t}; o[n_] := NDSolve[Join[{Equal @@ ...

6

As has been noted by ruebenko in the comments, there does seem to be a bug in the handling of infinite-range Bessel function integrals when MinRecursion and MaxRecursion are both set to non-default values. For instance, even the simple NIntegrate[BesselJ[0, x], {x, 0, ∞}, MinRecursion -> 10, MaxRecursion -> 15] chokes with a NIntegrate::minmax ...

3

@ruebenko confirmed that this is a bug. Here's a workaround: currentTimeStep = 0; monitor[t_?NumericQ] := (currentTimeStep = t; 1) pfun = ParametricNDSolveValue[{y'[t] == a monitor[t] y[t], y[0] == 1}, y, {t, 0, 10000}, {a}, MaxSteps -> 100000]; Dynamic[currentTimeStep] pfun[1]

0

One option is to reduce the PrecisionGoal and AccuraryGoal. Those are usually set quite high and perhaps this is a viable option here. Setting , PrecisionGoal -> 4, AccuracyGoal -> 4 for the Ns=10 case runs in 12 seconds, where the original code runs in 20 seconds. It were good if you removed the calls to N and replaced the real numbers with exact ...

4

This integral equation is solvable using the LaplaceTransform technique: Clear[s, t]; eqn = y'[t] == -Integrate[y[t1] Exp[t1 - t], {t1, 0, t}] LaplaceTransform[eqn, t, s] (* ==> s LaplaceTransform[y[t], t, s] - y[0] == -( LaplaceTransform[y[t], t, s]/(1 + s)) *) Solve[%, LaplaceTransform[y[t], t, s]] (* ==> {{LaplaceTransform[y[t], t, s] -> ...

6

This is a classical shock-tube problem in which a initially diaphragm separates a hi-pressure, high-density region from one of lower pressure and density. The classical exact solution has multiple discontinuities, a shock wave and a contact-surface (density discontinuity) that propagate to the right, and a continuous rarefaction wave traveling into the ...

4

The angle signs mean that you're looking at a shorthanded expression. Take a look at FullForm[s] to see the real thing. The curly braces are there because NDSolve is going to give you consistently formatted results whether you look for only one function with only one solution or your solution involves more than one function, and more than one possible ...

3

By default, NDSolve uses a "shooting" method to satisfy the boundary conditions: it picks an initial condition and then evolves the equation to see what boundary values are produced by the initial conditions. It then rejiggers these initial conditions to produce the boundary values specified by the boundary conditions. Here is the Mathematica documentation ...

3

"LevinRule" should work splendidly here, I think: NIntegrate[-m Exp[-m] BesselJ[1, m]^2, {m, 0, Infinity}, Method -> "LevinRule", WorkingPrecision -> 20] -0.18196415067209554877 ruebenko's answer has given a closed form for this particular definite integral. Personally, I prefer it when the parameters of the elliptic integrals are ...

4

This is a bug. As a workaround for this specific integral you could use a symbolic solution: Integrate[-m*Exp[-m]*BesselJ[1, m]^2, {m, 0, Infinity}] (* (-3*EllipticE[-4] + 5*EllipticK[-4])/(5*Pi) *)

1

Here is a way: sol = ParametricNDSolveValue[{-y*a'[x] == (a[x]^2 - 1)*x + a[x], a[-10] == 0.6}, a, {x, -10, 10}, {y}] sol[1] NIntegrate[sol[y][t] /. y -> 1, {t, 0, 1}] or alternatively via a Table: ListLinePlot[ Table[NIntegrate[sol[y][t], {t, 0, 1}], {y, 0.1, 1, 0.1}]]

6

This will work: \[Beta] = \[Omega]0/4; \[Omega] = 2 \[Pi]; \[Omega]0 = 3/2 \[Omega]; tend = 500; s = ParametricNDSolveValue[{x''[t] + 2 \[Beta] x'[t] + \[Omega]0^2 Sin[ x[t]] == \[Gamma] \[Omega]0^2 Cos[\[Omega] t], x[0] == 0, x'[0] == 0}, x, {t, 0, tend}, \[Gamma], MaxSteps -> Infinity]; Manipulate[ Plot[Evaluate[D[s[\[Gamma]][t], ...

1

The issue is not with NDSolve but plotting: \[Beta]=\[Omega]0/4;\[Gamma]=1/5;\[Omega]=2 \[Pi];\[Omega]0=3/2 \[Omega]; tend = 500; s = NDSolveValue[{x''[t] + 2 \[Beta] x'[t] + \[Omega]0^2 Sin[ x[t]] == \[Gamma] \[Omega]0^2 Cos[\[Omega] t], x[0] == 0, x'[0] == 0}, x, {t, 0, tend}, MaxSteps -> Infinity] Plotting a smaller subrange helps: ...

3

Using exact values (e.g. setting γ = 1/5 instead of γ = 0.2) and increasing the MaxSteps option should yield a reasonable result. Below, I add the filigree of also using a method I tend to prefer for problems with oscillatory solutions (although I've found that just setting Method -> "StiffnessSwitching" without increasing MaxSteps still works well): xp ...

1

If you just want to do a simple cartoon (as opposed to evaluating the solution components of the Lorenz equations at particular values), you can just directly extract the points generated by NDSolve[]. Here's one way to go about it: sol[σ_] := Transpose[Through[{x, y, z}["ValuesOnGrid"] /. First @ NDSolve[{x'[t] == σ (y[t] - x[t]), ...

3

I got it to work - had to redefine the boundary conditions with Derivative - I'll leave this up here if it's useful to others - working code block is ic = {Ef[r, 0] == kd, Derivative[1, 0][Ef][rn, t] == 0, Derivative[1, 0][Ef][ro, t] == 0};

3

You can define a solution function depending on your parameter and then use it for the animation : sol[sigma_] := {x, y, z} /. NDSolve[{x'[t] == sigma (y[t] - x[t]), y'[t] == 28 x[t] - y[t] - x[t] z[t], z'[t] == x[t]*y[t] - 8/3 z[t], x[0] == z[0] == 0, y[0] == 2}, {x, y, z}, {t, 0, 25}][[1]] Animate[ With[{f = sol[sigma]}, ...

0

The problem may be due to confusion over numerical values and dummy variables. The documentation states that NInverseFourierTransform[expr,w,t] gives a numerical approximation to the inverse Fourier transform of expr at the numerical value t, where expr is a function of w. In your question, you have reversed t and w, which is not a real problem since these ...

12

Here is a code that is about 2 orders of magnitude faster. We will use a finite element method to solve the issue at hand. Before we start, note however, that the transition between the Dirichlet values should be smooth. We use the finite element method because that works for general domains and some meshing utilities exist here and in the links there in. ...

3

An alternative approach is to form an approximate interpolating function from your actual function, which is (often) cheaper to integrate. To wit, nsin = FunctionInterpolation[Sin[x], {x, 0, 10}]; ni = Derivative[-1][nsin]; Plot[ni[t], {t, 0, 10}] In this case, we know that the integral of sine can be expressed simply, so we can compare the exact ...

0

Another possibility is to sieze control over the evaluation of f[t]. For instance: ListLinePlot[Table[f[t], {t, 0, 10, 0.1}]] evaluates f[t] at exactly the points you choose. In this case, it gives substantially the same plot as the original, and is about 10 times faster.

9

You can express your integral in terms of a differential equation and use NDSolve. Since NDSolve builds up the solution as it goes, this is typically much faster. Clear[y]; y[x_] = y[x] /. First[ NDSolve[{y'[x] == Sin[x], y[0] == 0}, y[x], {x, 0, 10}] ]; t = AbsoluteTime[]; Plot[y[t], {t, 0, 10}] AbsoluteTime[] - t

1

If you increase the working precision this will work: Needs["DifferentialEquationsNDSolveProblems"]; Needs["DifferentialEquationsNDSolveUtilities"]; k1 = 1*10^13; k2 = 1*10^6; k3 = 2000; s = NDSolve[{x'[t] == -k1*x[t]*y[t], y'[t] == -k3*y[t] + k1*x[t]*y[t], x[0] == 1*10^(-8), y[0] == 1*10^(-14)}, {x, y}, {t, 0, 1/1000}, WorkingPrecision -> ...

Top 50 recent answers are included