Mathematica Asked on June 27, 2021
I’ve looked up a few examples of similar errors for ODE’s with singular points, but have yet to figure out a way to solve this problem for my case. A simple patch up where the singular point is ‘skipped’ would be sufficient for me, but I would really appreciate some help.
Edit: my point being that ParametricNDSolve
is giving the ParametricNDSolve::ndsz: At t == 0.13801826916665527`, step size is effectively zero; singularity or stiff system suspected.
error
So I have a velocity equation V(x,t) for a particular as a function of (x,t) and some other parameters:
Edit: in my original unedited post, I think the code was copying across incorrectly. The following should show the issues with the ParametricNDSolve
hitting the singularities.
ClearAll["Global`*"]
v = 0.3;
k0r = k0;
k0l = k0;
[Sigma]r = [Sigma];
[Sigma]l = [Sigma];
A[x_, t_, k0_, [Sigma]_, [Alpha]_] =
Sqrt[2/[Pi]] [Sigma]r Sqrt[(1 - v)/(
1 + v)] [Alpha] E^(-2 ((1 - v)/(1 + v)) (t - x)^2 [Sigma]r^2) -
Sqrt[2/[Pi]] [Sigma]l Sqrt[(1 + v)/(
1 - v)] (1 - [Alpha]) E^(-2 ((1 + v)/(
1 - v)) (t + x)^2 [Sigma]l^2) +
Sqrt[(1 + v)/(
1 - v)] (-Sqrt[(2/[Pi])] Sqrt[[Sigma]r [Sigma]l] Sqrt[k0l/k0r]
Sqrt[[Alpha] (1 - [Alpha])]
E^(-((1 + v)/(1 - v)) (t + x)^2 [Sigma]l^2 - ((1 - v)/(
1 + v)) (t - x)^2 [Sigma]r^2) (Cos[
k0r (-Sqrt[((1 - v)/(1 + v))] (t - x)) +
k0l Sqrt[(1 + v)/(1 - v)] (t + x)] -
2*[Sigma]l^2/k0l Sqrt[(1 + v)/(
1 - v)] (t + x) Sin[
k0r (-Sqrt[((1 - v)/(1 + v))] (t - x)) +
k0l Sqrt[(1 + v)/(1 - v)] (t + x)])) +
Sqrt[(1 - v)/(
1 + v)] (Sqrt[2/[Pi]] Sqrt[[Sigma]r [Sigma]l] Sqrt[k0r/k0l]
Sqrt[[Alpha] (1 - [Alpha])]
E^(-((1 + v)/(1 - v)) (t + x)^2 [Sigma]l^2 - ((1 - v)/(
1 + v)) (t - x)^2 [Sigma]r^2) (Cos[
k0r (-Sqrt[((1 - v)/(1 + v))] (t - x)) +
k0l (Sqrt[(1 + v)/(1 - v)] (t + x))] +
2*[Sigma]r^2/k0r Sqrt[(1 - v)/(
1 + v)] (t - x) Sin[
k0r (-Sqrt[((1 - v)/(1 + v))] (t - x)) +
k0l Sqrt[(1 + v)/(1 - v)] (t + x)]));
B[x_, t_, k0_, [Sigma]_, [Alpha]_] =
Sqrt[2/[Pi]] [Sigma]r Sqrt[(1 - v)/(
1 + v)] [Alpha] E^(-2 ((1 - v)/(1 + v)) (t - x)^2 [Sigma]r^2) +
Sqrt[2/[Pi]] [Sigma]l Sqrt[(1 + v)/(
1 - v)] (1 - [Alpha]) E^(-2 ((1 + v)/(
1 - v)) (t + x)^2 [Sigma]l^2) +
Sqrt[(1 + v)/(
1 - v)] (Sqrt[2/[Pi]] Sqrt[[Sigma]r [Sigma]l] Sqrt[k0l/k0r]
Sqrt[[Alpha] (1 - [Alpha])]
E^(-((1 + v)/(1 - v)) (t + x)^2 [Sigma]l^2 - ((1 - v)/(
1 + v)) (t - x)^2 [Sigma]r^2) (Cos[
k0r (-Sqrt[((1 - v)/(1 + v))] (t - x)) +
k0l Sqrt[(1 + v)/(1 - v)] (t + x)] -
2*[Sigma]l^2/k0l Sqrt[(1 + v)/(
1 - v)] (t + x) Sin[
k0r (-Sqrt[((1 - v)/(1 + v))] (t - x)) +
k0l Sqrt[(1 + v)/(1 - v)] (t + x)])) +
Sqrt[(1 - v)/(
1 + v)] (Sqrt[2/[Pi]] Sqrt[[Sigma]r [Sigma]l] Sqrt[k0r/k0l]
Sqrt[[Alpha] (1 - [Alpha])]
E^(-((1 + v)/(1 - v)) (t + x)^2 [Sigma]l^2 - ((1 - v)/(
1 + v)) (t - x)^2 [Sigma]r^2) (Cos[
k0r (-Sqrt[((1 - v)/(1 + v))] (t - x)) +
k0l (Sqrt[(1 + v)/(1 - v)] (t + x))] +
2*[Sigma]r^2/k0r Sqrt[(1 - v)/(
1 + v)] (t - x) Sin[
k0r (-Sqrt[((1 - v)/(1 + v))] (t - x)) +
k0l Sqrt[(1 + v)/(1 - v)] (t + x)]));
which I am using ParametricNDSolve
to solve:
params = {x, t, k0, [Sigma] , [Alpha]};
V[x_, t_,
k0_, [Sigma]_, [Alpha]_] = (A[x, t, k0, [Sigma], [Alpha]])/(B[
x, t, k0, [Sigma], [Alpha]]);
Startpoints[k0i_, [Sigma]i_, vi_, [Alpha]i_, fr_, to_, st_,
nu_] := (samps =
Table[B @@ (params /. {k0 -> k0i, [Sigma] -> [Sigma]i,
v -> vi, [Alpha] -> [Alpha]i} /. t -> -2), {x, fr, to, st}];
tab = Table[Plus @@ #[[1 ;; n]]*st, {n, 1, Length[#]}] &[samps];
Flatten[
Table[FirstPosition[tab, x_ /; x > n, {Length[tab]}] - 1, {n, 0,
1, 1/(nu - 1)}]]*st + fr)
subs = {k0 -> 5, [Sigma] -> 1, v -> 0, [Alpha] -> Sqrt[1]/2};
de = {xtraj'[t] == V @@ (params /. subs /. {x -> xtraj[t]}),
xtraj[-2] == x0};
parasols = ParametricNDSolve[de, xtraj[t], {t, -2, 2}, x0];
startpoints =
Startpoints @@ ({k0, [Sigma], v, [Alpha], -3, 3, 0.01, 50} /.
subs);
trajectories =
Table[{xtraj[t][x0], t} /. parasols, {x0,
startpoints}]; ParametricPlot[trajectories, {t, -2, 2},
AspectRatio -> 1/2, PlotStyle -> Red]
There are singular points in the V(x,t) equation, most obviously when the denominator equals zero.
Is there a way for me to tell Mathematica to ‘skip’ the singular points? At the singularity itself I expect the trajectory to go horizontal (the velocity being infinite for that point) before going backwards for a while, turning around and then hitting another singular point and continuing on. If one runs the code above, you can see that it Mathematica doesn’t want to plot the initial conditions for which there are singular points.
Get help from others!
Recent Answers
Recent Questions
© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP