Tag Info

New answers tagged differential-equations

0

The following may be slightly more direct than the approach given in the answer. I don't split up the wave function into real and imaginary parts, so I can work with the original equations, making it easier to spot errors. The equation I want to check is eqn, and the equation of motion is nLS. I use two different function names for $\psi$ and its complex ...

1

No, your Fao[ ] doesn't work as posted. Ro = .007; Caorta = 1/.48; k = 110; ω = 2 π; x = 1000000; Rsystemic = 3.1; Vo = 108; Cheart = 6.67; m[ω_, t_] := 1/2*k*(1 + Cos[ω t]) + 10 PaoMwt[ω_, t_] := m[ω, t] - Paorta[t] pw[fun_] := Piecewise[{{Ro, fun > 0}}, x*Ro] VhMwt[ω_, t_] := m[ω, t] - ((Vheart[t] - Vo)*m[ω, t] + Vheart[t]/Cheart) sol = ...

1

Solving this problem is probably going to get you in to the deep weeds of Method specifications in NDSolve. When Mathematica encounters the code as originally stated, it defaults to a finite element method. The problem is that Mathematica's finite element methods can't handle equations with higher than second-order derivatives in the spatial variables, ...

6

Alternatively, this can be treated as a boundary-value problem. a = 1.1; b = 1.2; tf = 20; eqn = {x''[t] + (a + I*b)*x[t] == 0}; inits = {x[0] == 1, x[tf] == 0}; sys = Join[eqn, inits]; sol = NDSolve[sys, x[t], {t, 0, tf}]; Plot[Evaluate[ReIm[x[t]] /. sol], {t, 0, tf}] MichaelE2's caveat applies here too, "It might be significantly more difficult in a ...

1

Although this Answer does not use Floquet analysis, it does cast the equations into a more useful form and provide a sample numerical solution. X[t] = {{x1[t], x2[t]}}\[Transpose]; b[t] = b1 E^(-I δ t) + b2 E^(I δ t); a[t] = a1 E^(-I δ t) + a2 E^(I δ t); M[t] = {{0, FullSimplify[Conjugate[b[t]], t ∈ Reals && δ ∈ Reals && b1 ∈ ...

5

If you can define an objective function that measures the size of the solution, you could optimize it. This is simple to do on the simple test case. It's a linear system, so the convergence/divergence will depend only on the ratio x'[0]/x[0]. One can optimize varying x'[t] for x[0] == 1 and test x[0] == 0 separately (best to do x[0] first, but I omit the ...

2

I rewrote your code and it works for me so I suspect there is a bug in your code. Some tips: when you write out the equation as you have done it is very easy to make mistakes. What I recomment is that you define $P=\{P_1,P_2,P_3\}$ and $Q=\{Q_1,Q_2,Q_3\}$ only once (instead of retyping it every time you need it) and then write the equation directly in ...

1

I can reproduce this behavior like this: y' = y^2; s = NDSolve[{y'[x] == y[x], y[0] == 1}, y, {x, 0, 1}] NDSolve::dvnoarg: The function y appears with no arguments. >> (* NDSolve[{(y^2)[x] == y[x], y[0] == 1}, y, {x, 0, 1}] *) Executing Clear[y, x] or even ClearAll[y, x] won't work because the problem is stored in the SubValues for Derivative: ...

1

This seems like a cell/formatting issue. Try using a new notebook. Btw, om my pc (math 10) your example gives the correct result. Cheers

2

And we Go: ClearAll[b, c, k, x, t T0]; SOL = First @ Assuming[{b, c, k, T0} ∈ Reals, DSolve[{x'[t] == b x[t] - (c - x[t])^2/(4 k), x[T0] == 0}, x[t], t]]; // Quiet FullSimplify[SOL] $\left\{x(t)\to 2 \sqrt{b} \sqrt{k} \sqrt{b k+c} \tanh \left(\frac{\sqrt{b} (t-\text{T0}) \sqrt{b k+c}}{2 \sqrt{k}}-\tanh ^{-1}\left(\frac{2 b k+c}{2 ... 4 As noted by @bbgodfrey, the "shooting" algorithms that Mathematica tends to use are not well-adapted to this particular equation. Better would be some kind of relaxation method, which is what Mathematica uses (I think) for solving PDEs on a mesh. And an ODE is just a PDE in one dimension, so let's try solving this equation on a one-dimensional mesh: ... 4 To see why NDSolve has difficulty with this problem for very small e, consider that NDSolve solves this two-point boundary value problem by some form of shooting. In other words, it varies y'[min] until one is found that yields the desired y[max]. However, as e becomes very small, the sensitivity of y[max] to y'[min] becomes great, because the differential ... 6 DSolve can handle this. Clear[y]; y[x_, e_] = y[x] /. DSolve[{ e y''[x] + y'[x] x == -e Pi^2 Cos[Pi x] - Pi x Sin[Pi x], y[-1] == -2, y[1] == 0}, y[x], x][[1]] // Simplify Manipulate[ Plot[y[x, e], {x, -1, 1}], {{e, 0.01}, 0.0001, 0.1, Appearance -> "Labeled"}] 3 Replace ,DirichletCondition[cA[t, r, z] == 0, True] ,DirichletCondition[cB[t, r, z] == 0, True] ,DirichletCondition[cC[t, r, z] == 0, True] by , DirichletCondition[cA[t, r, z] == 0, z == z1 || z == z2] , DirichletCondition[cB[t, r, z] == 0, z == z1 || z == z2] , DirichletCondition[cC[t, r, z] == 0, z == z1 || z == z2] , DirichletCondition[cA[t, r, z] == ... 2 Since you are specifically asking about versions below 10, it may be useful to point out that this problem is equivalent to the electrostatics problem of finding the potential in a region bounded by conductors held at fixed voltages. This can be solved, e.g., with the simple relaxation method I implemented in this answer, where I actually allow for lots of ... 2 This is the outline of an answer. To begin, assume that R u v is meant instead of Ruv. As noted by MarcoB, attempting to solve the equations by brute force, as suggested in the Question, takes forever. Instead, first simplify the three linear ODEs by adding constants to the dependent variables. eq = {x'[t] == R x[t] v + R y[t] u + R u v - x[t] - u, ... 1 If you have v.10 you can explicitly use the finite element method: Needs["NDSolveFEM"] mesh = ToElementMesh[Rectangle[{0, 0}, {10, 10}]] sol = First@NDSolveValue[{Laplacian[w[x, y], {x, y}] == 0, DirichletCondition[w[x, y] == 100, y == 0], DirichletCondition[w[x, y] == 400, y == 10], DirichletCondition[w[x, y] == 0, x == 0], ... 0 First: Getting the right solution Using InternalWithLocalSettings to set the Method -> Reduce option on Solve, and adding the assumption x ∈ Reals, we get the correct answer: sys = {y'[x] == (1 + 2 x) Sqrt[y[x]], y[0] == 1}; With[{opts = Options[Solve]}, InternalWithLocalSettings[ SetOptions[Solve, Method -> Reduce], s = Assuming[x ∈ Reals, ... 1 If we run WolframAlpha["y'[x] == (1 + 2 x)Sqrt[y[x]], y[0] == 1"], we can see the step by step solution, which exposes the issue at hand. WolframAlpha["y'[x] == (1 + 2 x)Sqrt[y[x]], y[0] == 1"] (* click the Step-by-step solution *) Now, I cut off the image here on purpose. Notice to solve for y[x], we had to square both sides. This could give ... 2 [Edit: Rewrote and expanded explanation of how to set up DSolve and why it fails.] You were almost there.... First, I want to show how such a problem should be specified, if you want DSolve to handle it all for you. DSolve can do it, provided such a system is solvable by Mathematica to begin with. Understanding how DSolve works with such equations helps ... 3 EDIT: added handling special case of x == -1/2 Expanding on comment by@murray eqns = {y'[x] == (1 + 2 x) Sqrt[y[x]], y[0] == 1}; sol = DSolve[eqns, y, x] {{y -> Function[{x}, (1/4)*(4 - 4*x - 3*x^2 + 2*x^3 + x^4)]}, {y -> Function[{x}, (1/4)*(4 + 4*x + 5*x^2 + 2*x^3 + x^4)]}} eqns /. sol ... 2 First, your use of v[t] doesn't add any new information that isn't already in s[t]. It's not really a variable we solve for, so just use s[t]. Next, you are overspecifying the problem like @b.gatessucks says. The differential equation has one constant of integration, because it is first order. We can only apply one boundary condition in this case. Lets ... 3 Building up from DSolveValue... sol = DSolveValue[{s'[t] == A Exp[k t], s[0] == 0}, s, t]; s0 = Rationalize[2.37 10^16, 0]; v0 = 299792458; NSolve[{sol[6000] == s0, sol'[6000] == v0}, {A, k}, Reals] (* {{A -> 4.72755*10^13, k -> -0.00199474}} *) 3 sol = First@DSolve[{s'[t] == a*Exp[k t], s[0] == 0}, s[t], t]; Displacement: S[t] = s[t] /. sol $$\frac{a \left(e^{k t}-1\right)}{k}$$ Speed: V[t] = D[s[t] /. sol, t] $$a e^{k t}$$ Solution. NSolve[{Rationalize[2.37*10^16, 0] == (a (-1 + Exp[k t]))/k, 299792458 == a Exp[k t]} /. t -> 6000, {k, a}, Reals]$$\{\{k\to ... 2 If you look at the output of your original data, you'll notice that you're getting small imaginary parts in your results. This implies to me that Mathematica is trying to take square roots, which is why it's failing when$E_z$is going to zero; and the main culprit, I believe, is your constraint that$x'(s)^2 + z'(s)^2 = 1$. Let's differentiate that ... 3 As noted in the comment above, k1 etc must be called with its argument list. With this change the r1 etc equations become r1 = k1[r, ϕ] dens[r, ϕ]^2 nH[r] nOH[r]; r2 = k2[r, ϕ] dens[r, ϕ]^2 nH[r] nH2O[r]; r17 = k17[r, ϕ] dens[r, ϕ]^2 nH2[r] nO[r]; r18 = k18[r, ϕ] dens[r, ϕ]^2 nH2[r] nOH[r]; r62 = k62[r, ϕ] dens[r, ϕ]^2 nO[r] nOH[r]; r63 = k63[r, ϕ] ... 7 First thing: I have very limited experience with PDEs but what you present doesn't look like the Laplace equation. AFAIK the r and r^2 need to be in the denominator. You try to solve this in polar coordinates. Maybe it is possible to use NDSolve if you can do your real problem in Cartesian coordinates, because it seems the new region functionality within ... 2 Unfortunately the "shooting" is done by NDSolveProcessEquations, which converts the BVP into an IVP. Thus NDSolveReinitialize[ndssdata] is basically operating on an IVP, and you won't be able to approach your problem in this way. If we examine the state data, we see that the initial conditions are already there. {sdb, sdf} = ndssdata@"SolutionData" (* ... 1 I have not been able to obtain a solution in a general way. Nonetheless, some progress can be made. If we happen to know the value of V at phi == 0 and phi == 2 Pi, then using it of course gives the desired solution, s = NDSolveValue[{r^2*D[V[r, phi], r, r] + r*D[V[r, phi], r] + D[V[r, phi], phi, phi] == 0, V[1, phi] == 0, V[2, phi] == ... 5 WhenEvent[cond, act] works in time, i.e., an event happens only when a time step causes the condition cond to change from False to True, save the special cases such as f == c described in the documentation. Those subtleties aside, the main thing to understand in using WhenEvent in the method of lines is what is substituted for the dependent variables such ... 1 These are the transitions I usually use: (both are$C^{\infty}\$) f[x_] = Piecewise[{{(Erf[Sqrt[2 π] ArcTanh[x]] + 1)/2, -1 < x < 1}}, UnitStep[x]] g[x_] = Piecewise[{{(Tanh[Sqrt[2] Tan[π/2 x]] + 1)/2, -1 < x < 1}}, UnitStep[x]] Plot[{UnitStep[x], f[x], g[x]}, {x, -2, 2}]

5

To get a fixed step size with the BDF method you can lower the AccuracyGoal and PrecisionGoal to increase the adaptive step sizes and then use MaxStepSize to limit the step size to any value you want. Get an example stiff system from the documentation: Needs["DifferentialEquationsNDSolveProblems"]; Needs["DifferentialEquationsNDSolveUtilities"]; system ...

8

You can make a change of variable to solve the problem. Here I'll use dChange for this task: r0 = 0.5; eqn = {k[z] f[z] + f'[z] == 0, k'[z] == 0, f[-1] == Exp[2], f[1] == 1}; c = Piecewise[{{r0, z > 0}}, 1]; neweqn = dChange[eqn, f[z] == c g[z]]; {solg, solk} = NDSolveValue[neweqn, {g, k}, {z, -1, 1}]; Plot[c solg[z], {z, -1, 1}]

11

Working solution One can manually implement the shooting method with ParametricNDSolveValue and FindRoot: psol = ParametricNDSolveValue[{k f[z] + f'[z] == 0, f[-1] == Exp[2], WhenEvent[z == 0, f[z] -> r0 f[z]]}, f, {z, -1, 1}, {k, r0}]; k0 = k /. FindRoot[psol[k, 0.5][1] == 1, {k, 1}] (* 0.653426 *) Plot[psol[k0, 0.5][z], {z, -1, 1}] ...

9

Analytical Solution A corrected version of the now deleted answer by Nasser, that uses a trick explained in this answer by Jens: Clear[f, if, y, x, g]; f /: Integrate[f[x_], x_] := if[x]; SetAttributes[if, {NumericFunction}]; sol = Integrate[f[x], {x, 0, 2}]; bc2 = sol == 5 /. if -> g -g[0] + g[2] == 5 With g beeing the antiderivative of y. ode = ...

10

s = DSolveValue[{y''[x] - y[x] == (x^2), y[0] == 1}, y[x], x]; First@Solve[Integrate[s, {x, 0, 2}] == 5, C[1]]; s /. % {* E^-x (3 - 2 E^x - (9 + 26 E^2)/(3 (-1 + E^2)^2) + (E^(2 x) (9 + 26 E^2))/ (3 (-1 + E^2)^2) - E^x x^2) *} Alternative solution The answer by Karsten7 suggests the following: With y[x] replaced by g'[x], D[DSolveValue[{g'''[x] - ...

3

Because the equations are nonlinear, there is no assurance that equilibrium solutions even exist for a given set of k and initial conditions. However, if they do exist, an alternative, and perhaps more informative, approach is to compute them directly: xs = Solve[Thread[Table[0, {i, Length[des[[1]]]}] == des[[2]] /. x[n_][t] -> x[n]], Table[x[i], ...

0

I think for a Matlab user this would probably be something relatively simple to understand but still making decent use of Mathematica: k = 1.0; solution = Table[ NDSolveValue[{ (1/Pi)^2*(1 + k*Sin[th[zr]]^2)*th''[zr] + (1/Pi)^2*k*Sin[th[zr]]* Cos[th[zr]]*(th'[zr])^2 + vr^2*Sin[th[zr]]*Cos[th[zr]] == 0, th[0] == 0.0, th[1] == Pi/2 }, th, ...

2

You can check that a differential equation is solved correctly by plugging it in to the original equation and seeing that it holds. For your case: sol = DSolve[z''[x] == λ*x^(3/4)*z[x]/Sqrt[x], z[x], x] f[x_] := sol[[1, 1, 2]] so f[x] is the candidate function. Now verify that: FullSimplify[D[f[x], {x, 2}]] == FullSimplify[λ*x^(3/4)*f[x]/Sqrt[x]] Since ...

4

Based on the OP's original "set-up": des = {x[1]'[t], x[2]'[t], x[3]'[t], x[4]'[t], x[5]'[t], x[6]'[t], x[7]'[t], x[8]'[t], x[9]'[t]} == {-k[1] x[1][t] x[3][t] + k[2] x[5][t] + k[3] x[5][t] - k[7] x[1][t] x[7][t] + k[8](*[t]*) x[8][t], -k[4] x[2][t] x[4][t] + k[5] x[6][t] + k[6] x[6][t] - k[9] x[2][t] x[7][t] + k[10] x[9][t], ...

2

It would help a lot if you write each of the equations separately and name them appropriately; e.g. eq1, eq2, eq3... etc. Then you can use NDSolve routine to try and solve the system (along with the boundary/initial conditions). Also, try formatting your code properly so that the stack-exchange community can help you more easily. From a quick view of the ...

3

An example perturbing a simple RL circuit when it's about to reach the steady state. R = 1; L = 1; V = 1; SeedRandom[43]; sol = NDSolve[{R i[t] + L i'[t] == V, i[0] == 0, WhenEvent[{i'[t] == .01}, i[t] -> i[t] + RandomReal[{-.5, .5}]], WhenEvent[{i'[t] == -.01}, i[t] -> i[t] + RandomReal[{-.5, .5}]]}, i, {t, ...

2

p = ParametricNDSolve[{y''[x] == (l x^(3/4) y[x])/Sqrt[1 - x], y[0] == 1, y'[0] == 1/2}, y, {x, 0, 1}, l] Show@Quiet@Table[Plot[(y /. p)[l][x], {x, 0.1, 1}], {l, 1, 3}]

3

Take a look at X'[t]/Z'[t]: This is not the angle of rotation that you need. What you need is ArcTan[Z'[t],X'[t]], which looks like this If you change that in your code, you get the desired result.

0

Maybe you can disregard that region of numerical solving by specifying another interval for NDSolve that does not include the region where you observe the oscillations. Another solution would be to impose such boundary/initial conditions that make this irrelevant oscillations disappear. But this will also affect your solution, so unless these new conditions ...

5

The following isn't quite an answer to the question—rather than having DSolve remind you that you might not have all the solutions, the procedure below tries to find other ones (and succeeds in the case you've outlined.) As I mentioned in my comment above, the "missing" solution is the envelope of the general family of curves that you found: In ...

2

Please behave like a good citizen of the site and read the answers you receive more carefully. The following is a direct application of my answer to your previous question here. w[k_, ω_, t_] := 1/2*k*(1 + Cos[ω t]) + 10; pnd = ParametricNDSolve[{ Paorta'[t] == 1/Caorta ((w[k,ω,t] - Paorta@t)/ Piecewise[{{ρ, w[k,ω,t] - Paorta@t ...

4

At x = 0 we have a singular point. Calculation of analytical solutions will take eternity. sol = NDSolve[{(4*y[x] - 1 - 4 x - 6 x^2) (x^2 y'[x] + (1 - x) y[x]) == 3*y[x] (1 + 3 x (1 + 2 x + 2 (x^2))), y[0.00001] == 1}, y[x], {x, 0, 10}] // Quiet; Plot[y[x] /. sol, {x, 0, 10}] If we want to express solution using Taylor Series: initconds = {y[0] == ...

3

This is mostly because your code is completely unrelated to your equation. What is w? Why do you use = when you want ==? Where does the 5 come from when you want 1.5? Why are all of your initial values 0? Here is something to start with but first go to Evaluation -> Quit Kernel sol = NDSolveValue[{ u'[t] == 1.5 u[t] v[t] - u[t], v'[t] == -1.5 ...

1

In the absence of a functioning code in the Question, I used the straightforward L = 1.5; s = NDSolveValue[{D[u[x, t], t] + u[x, t] D[u[x, t], x] == 0, u[L, t] == u[-L, t], u[x, 0] == Max[1 - x^2, 0]}, u, {x, -L, L}, {t, 0, 1}, MaxStepSize -> 0.001]; Plot3D[s[x, t], {x, -L, L}, {t, 0, 1}] Despite a warning message about error estimates, the ...

Top 50 recent answers are included