A message from our CEO about the future of Stack Overflow and Stack Exchange. Read now.

New answers tagged

1

You need to control the element mesh used to generate the solution. Add the option Method -> {"FiniteElement", "MeshOptions" -> {"MaxCellMeasure" -> 0.001}} The error is reduced: ContourPlot[psixx == 0, {x, y} \[Element] regPl, PlotRange -> All, PlotLegends -> Automatic, AspectRatio -> 0.8] Here is the mesh created with the ...


2

This is more an extended comment than a formal answer. With the code and parameters given in the question, the first plot in the question should be something like Plot[Evaluate@Table[(i[k] /. sol)[t], {k, 1.5, 2.5, 0.2}], {t, 0, 10}, PlotRange -> All, ImageSize -> Large, AxesLabel -> {t, i}, LabelStyle -> {15, Bold, Black}] The maximum of ...


2

To plot the solution for a range of delc, make these changes to your code: (1) Remove delc = 1; Leave delc undefined. It will be your parameter. (2) Keep the ParametricDSolveValue command and its first argument, but change the other arguments to get s = ParametricNDSolveValue[ ... , {1/2*(V11[t] + V22[t] - 2*V12[t])^(-1)}, {t, 0, 100}, delc]; ...


2

Here is a working serial version of your problem: AbsoluteTiming[ myDataSerial = {}; initialConditions = Flatten[#, 1] &@ Table[{i, j}, {i, 0., 3(*10*), 0.1}, {j, 0., 3(*10*), 0.1}]; myEvent = WhenEvent[x[t] == 0.8, AppendTo[myDataSerial, t]]; Map[NDSolve[{x''[t] == -x[t], {x[0], x'[0]} == #, myEvent}, x, {t, 0, 10}] &, ...


1

You can always supply the statement as an argument for NDSolve, Method -> {"TimeIntegration" -> "ExplicitEuler"}, or Method -> {"TimeIntegration" -> {"ExplicitRungeKutta", "DifferenceOrder" -> d}} where d is the chosen order you want.


4

This is the famous Lotka-Volterra predator-prey model. I don't think there is a closed-form solution for the time dynamics but there is one for the orbits in the phase plane, which is described on that wikipedia page. You can derive it in Mathematica with DSolve[{y'[x] == (-c*y[x] + d*y[x]*x)/(a*x - b*x*y[x])}, y, x]


5

x[0] = 0;(*Supplementary conditions*) y[0] = 4.; h = 0.1; f[x_, y_] := x^2 - y^2 x[n_] := n h + x[0] y[n_] := y[n] = y[n - 1] + h f[x[n - 1], y[n - 1]] data = Table[{x[i], y[i]}, {i, 0, 10}]; Show[Plot[ Evaluate[yy[x] /. First[NDSolve[{yy'[x] == x^2 - yy[x]^2, yy[0] == 4}, yy[x], {x, 0, 1}]]], {x, 0, 1}, PlotStyle -> {Red, Dashed}], ...


2

Let DSolve do the job. eq = F'[x] + Sqrt[ a^3 x F[x]] == 0; dsol = DSolve[{eq, F[1] == 0}, F, x] (* {{F -> Function[{x}, 1/9 (a^3 - 2 a^3 x^(3/2) + a^3 x^3)]}} *) sol = Solve[1/a == 2 Integrate[F[x] /. First@dsol, {x, 0, 1}], a, Reals] (* {{a -> -10^(1/4)}, {a -> 10^(1/4)}} *) Fsol[x_] = F[x] /. First@dsol /. sol[[2]] // FullSimplify ...


2

Try DSolve (without timerange): DSolve[{(v1[t] - t) n1'[t] + n1[t] n1'[t] == 0, 3 sigma1 n1[t] n1'[t] + (v1[t] - t) v1'[t] + phi'[t] == 0, (v2[t] - t) n2'[t] + n2[t] v2'[t] == 0, 3 Q sigma2 n2[t] n2'[t] + (v2[t] - t) v2'[t] + Q phi'[t] == 0, phi'[t] - (1/ne[t]) ne'[t] == 0, phi'[t] - (sigmas1/ns1[t]) ns1'[t] == 0, (vs20 + vs2[t] - t) ns2'[t] + ns2[t] ...


0

It turns out that the problem isn't about precision, but about the setup of StopIntegration, where I set to be WhenEvent[z[ρ] == 10^-3, "StopIntegration"], WhenEvent[z[ρ] == 1.1` Sqrt[zt^2 - ρc^2], "StopIntegration"], however, if we have a closer look on the solution $z[ρ]$, we would find it's increasing at first, then decreasing until hitting the ρ-axis. ...


1

You need to specify initial conditions: FirEqn = -D[IPF[z, t], z] - D[IPF[z, t], t] - IPF[z, t] SecEqn = -D[IPB[z, t], z] - D[IPB[z, t], t] - IPB[z, t] (* a wild guess *) ROCP = 1; RICP = 1; solIntEqn = NDSolve[{ FirEqn == 0, SecEqn == 0, IPF[0, t] == 4 Exp[-(t - 1)^2], IPF[5, t] == IPB[5, t]/ROCP, IPB[0, t] == 4 Exp[-(t - 1)^2]/RICP, IPB[...


7

This problem can be solved as follows. First, solve diff with h[z] treated as a parameter and the constants as listed in the queston to obtain an explicit expression for f'[z]. diff /. {h[z] -> h, n -> 1/2, i -> 1/2, L -> 2/5} (* (f'[z] (1/2 + 1/(2 (1 + (4 f'[z]^2)/(25 h^2))^(1/4))))/h == -h z *) Simplify[(h # - i f'[z]) & /@ %]; solf = ...


4

With little modifications it seems to wotk: s = ParametricNDSolveValue[{Derivative[2][y][x] + (a + b*(2 + (2/Pi)*ArcTan[x]))*y[x] == 0, y[-10] == Exp[I*10*Sqrt[a + b]], Derivative[1][y][-10] == (-I)*Sqrt[a + b]*Exp[I*10*Sqrt[a + b]]}, y, {x, -10, 10}, {a, b}] u = ParametricNDSolveValue[{Derivative[2][z][x] + (a + b*(2 + (2/Pi)*ArcTan[x]))*z[x] == 0, z[...


2

Let me give a little example. Since i am working with MMa version 8.0, where ParametricNDSolve is not jet implemented, here a little workaround. Please adapt it to ParametricNDSolve. I add intial conditions and look for minimum, not maximum, here. eqs = {alpha* Derivative[2, 0, 0][x][t, alpha, beta] + y[t, alpha, beta]*Derivative[1, 0, 0][x][t, alpha, ...


6

That's weird considering that the system can be solved one-by-one with Mathematica itself. Clear[x, y, z]; eqns = {x''[t] + x[t] == a, y''[t] == x[t]^2 - y[t], z[t] == 2 x[t]*y[t] - z''[t]}; x[t_] = x[t] /. DSolve[eqns[[1]], x[t], t][[1]] y[t_] = y[t] /. DSolve[eqns[[2]], y[t], t][[1]] z[t_] = z[t] /. DSolve[eqns[[3]], z[t], t][[1]]


1

Here is something to get you started: BCs and ICs need to be consistent. You need ICs and derivatives of ICs (you have second order time derivatives). Mixing u[x,t] and v[t] is most likely not going to work. eqn1 = D[u[x, t], x] + 5 D[u[x, t], x, x] + D[v[x, t], t, t] - 4 == 0; eqn2 = D[v[x, t], t] + D[u[x, t], t, t] + v[x, t] == 0; BCE1 = u[0, t] == 0;...


6

Use a decimal point for your coordinate values so that regions are converted to real numbers. ri = Rectangle[{40., 40.}, {60., 60.}]; rbi = Element[{x, y}, RegionBoundary@ri]; ro = Rectangle[{0., 0.}, {100., 100.}]; rdiff = RegionDifference[ro, ri]; rbo = Element[{x, y}, RegionBoundary@ro]; eqn = Laplacian[V[x, y], {x, y}] == 0; dci = DirichletCondition[V[x,...


2

Here's a proof of concept. Problems with code usually require the code for the problem to be analyzed, so while it accomplishes what is described, I don't know if it can be adapted to the OP's case. {sol} = NDSolve[{x'[t] == -y[t] - x[t]^3, y'[t] == x[t] - y[t]^3, x[0] == 1, y[0] == 0}, {x, y}, {t, 0, 20}]; xp[x_, y_, e_] := -y - e x^3; (* RHS for x'[t]...


6

Since this can be readily solved analytically using seperation of variables, another option is to use the analytical solution, and replace the paramaters in the solution without the need to solve it each time. Mathematica can't solve this analytically, but the analytical solution is $$ u \left( x,t \right) =\sum _{n=1}^{\infty }-4\,{\frac {\sin \left( nx ...


8

Look carefully at the documentation for NDSolve and Manipulate and see if you can understand how every character of this is working. You can even click on the orange Details and Options for each of those to get additional information. And then Manipulate[ sol=u/.NDSolve[{D[u[t,x],t]==α D[u[t,x],{x,2}]+a u[t,x], u[t,0]==0,u[t,Pi]==0,u[0,x]==x(Pi-x)},u,{...


1

Figured it out. Easiest way is to use IF to specify values for particular regions. Needs["NDSolve`FEM`"] outerR = 0.002; (*outer position*) innerR = 0.001; (*inner position*) WalinPos = 0.0011;(*inner position of wall*) WalOutPos = 0.0013; (*outer position of wall*) vel1 = 0.001; vel2 = 0; vel3 = -0.001; DT1 = 0.000000143; DT2 = 0.0143; DT3 = 0.000000143; ...


3

We can normalize the equations by 10^5. So that the oscillation period is of the order of 1, we make the change t1->t*10^5. Then it is possible to integrate up to 10^5, since the dynamics are also visible on this scale. \[Theta] = 34.5*\[Pi]/180; m2 = 1/2*7.3*10^-5*{{Cos[ 2 \[Theta]], -Sin[2 \[Theta]]}, {-Sin[2 \[Theta]], -Cos[ 2 \[Theta]]}};...


4

You can find examples of this for the numerical solver NDSolve in the FEMDocumentation/tutorial/NonlinearFiniteElementVerificationTests Specifically the following three tests are taken from the source you give: Diffusion—FEM-NL-Transient-1D-Single-Diffusion-0001 Diffusion—FEM-NL-Transient-1D-Single-Diffusion-0002 Diffusion—FEM-NL-Transient-1D-Single-...


1

The Helmholtz equation is a reaction-diffusion equation. There are many examples in the documentation. Maybe something like this: region = RegionDifference[Rectangle[{0, 0}, {5, 10}], Disk[{5, 5}, 3]]; op = -Laplacian[u[x, y], {x, y}] - u[x, y]; bcs = {DirichletCondition[u[x, y] == 0, x == 0 && 8 <= y <= 10], DirichletCondition[u[x, y] == ...


3

This is not an answer to your question but the solution you get is not what you think it is: When you specify Method -> {"FiniteElement", "MeshOptions" -> {MaxCellMeasure -> 0.005}} You are solving this a 2D spatial problem but you want to solve this as a 1D spatial and time dependent problem. The easiest way to do so is by using: heatdist = ...


4

Decreasing the MaxCellMeasure to $5\times10^{-5}$ improves the quality of the solution considerably while still taking only a few seconds to calculate, so I would suggest that as a viable compromise: heatdist = NDSolve[{HEATIMP == 0, BCTot, ic}, u, {r, 1, 2}, {t, 0.00000, 1}, Method -> {"FiniteElement", "MeshOptions" -> {MaxCellMeasure -> 5*...


4

Here is a partial solution. Begin by noting that pde1 does not involve T and so can be solved for F and the solution then inserted into pde2 to solve for T. The boundary condition F[0, y] + y == -y*Derivative[0, 1][F][0, y] should be integrated to yield -y/2. The boundary condition Derivative[1, 0][F][lb, y] == 1 is inconsistent with the initial condition ...


3

To get the first-order state-space equations you can use an internal function that NonlinearStateSpaceModel uses under the hood. NonlinearStateSpaceModel cannot handle this at the moment because it has algebraic equations, but the internal function has no problem with such equations. {{f, h, e}, x} = Control`DEqns`nonaffinestatespaceForm[ eqns, {PX[t], ...


2

NDSolveValue would allow you to get the requested expressions directly: sValue = NDSolveValue[eqns, {Abs[a1'[t]]^2, Abs[a1[t]]^2}, {t, 0, 0.00001}]; ParametricPlot[ sValue, {t, 0.000008, 0.000009}, PlotRange -> Automatic, AspectRatio -> 0.5 ]


6

Here is a possible direction to get started. I have based this example off of the Navier-Stokes example in the documentation from the Fluid Flow section. I have re-written the equations in inactive form: op = {Inactive[ Div][{{-0.2, 0}, {0, 0}}.Inactive[Grad][u[t, x, y], {x, y}], {x, y}] + Inactive[ Div][{{0, 0}, {0, -0.2}}.Inactive[Grad]...


0

I introduced two boundary conditions on $x$ which are $$ \frac{\partial f}{\partial x}f(0,t) = \frac{\partial f}{\partial x}f(L,t) = 0 $$ to facilitate the Laplace Transform technique. Using the Laplace transformation Clear[f, c, t, x] eqn = {D[f[x, t], t] - b D[f[x, t], {x, 2}] + k f[x, t] - g c[x, t] == 0, D[c[x, t], t] - k f[x, t] + g c[x, t] == 0}; ic ...


1

Since the first equation is totally independent of the second variable P, you can solve it analytically first. Saying "the second starts at t = 20", I think you mean P[20]==0 and for lower t. im = Kb = Yab = Kdm = Kt = Kdp = 1; mRNAsol = mRNA /. First@ DSolve[{Derivative[1][mRNA][t] == 1 - 1/E - mRNA[t], mRNA[0] == 0}, mRNA, t] (* Function[{t}, (-1 +...


1

I'm not sure I understood, but if you only need to add some lines to existing plot from cited question: T = 10; Y = ParametricNDSolveValue[{X'[t] == Boole[X[t] > 0], X[0] == x}, X, {t, 0, T}, {x}]; Show[Plot[{x, x + Range[10]}, {x, 0, 10}, PlotRange -> {{-10, 10}, {0, 10}}, PlotStyle -> {Blue, Green}], Table[ParametricPlot[{Y[x][t], t}, {t,...


1

Here's one way: im = Kb = Yab = Kdm = Kt = Kdp = 1; NDSolve[{mRNA'[t] == im (1 - Exp[-Kb Yab] - Kdm mRNA[t]), P'[t] == If[t > 20, Kt mRNA[t] - Kdp P[t], 0], mRNA[0] == 0, P[0] == 0}, {mRNA, P}, {t, 0, 400}]


2

This is a nonlinear equation. To be able to use FEM you have to get the equation into the coefficient form: Since the form you present is not in this coefficient form strange things happen. The easiest way around this is probably to use Activate and have the FEM solver figure it out - which seems to work in this case but is not a general cure. ClearAll[...


11

The following example makes use of my package, IGraph/M. Needs["IGraphM`"] First we generate some random points to use as Voronoi centres. We apply a few steps of Lloyd relaxation to make them more equi-distant. pts = RandomReal[{-1, 1}, {100, 2}]; Do[ pts = PropertyValue[{VoronoiMesh[pts, {{-1, 1}, {-1, 1}}], {2, All}}, MeshCellCentroid], {2} ] This ...


1

Regarding the initial conditions ic = {u[0, x] == t u0 x/lr, s[0, x] == 1}; there is an inconsistency because it is taken at $t = 0$. Making ic = {u[0, x] == 0, s[0, x] == 1}; it works fine. sols = NDSolve[{PDE1, PDE2, bcs1, bcs2, ic}, {u, s}, {x, ll, lr}, {t, 0, 10}] utx = Evaluate[u[t, x] /. sols] stx = Evaluate[s[t, x] /. sols] Show[Plot3D[utx, {x, ...


3

DSolve often cannot solve differential equations for which symbolic solutions actually do exist. So, out of curiosity, I attempted to solve these equations and made some progress, obtaining v'[l] in terms of u[l], and u[l] as an implicit function of l. Here is the calculation. Obtain v'[l]: Equal @@ (DSolve[eq2, u[l], l][[1, 1]]) /. Exp[2 C[1]] -> c ...


Top 50 recent answers are included