2
$\begingroup$

I'have a 2d-system of differential equations, analitycally solved, depending on a parameter. I know that, by continuity, there exist a minimum value of the parameter such that trajectory passes through a specific point. My problem is to find such value (also approximatively).

This is the system-equation:

ti = 0; yi = 0; zi = -.75; zf = -.5; eps=.01; sol = FullSimplify[DSolve[{y'[t] == 1/(2 zf) y[t] - u z[t], z'[t] == -1 + 1/zf z[t] + u y[t], y[ti] == yi, z[ti] == zi},{y[t],z[t]},t]];

The parameter is u.

I have to find the first u such that trajectorie intersect the point: z=zf; y=Sqrt[-2*(zf)*eps + eps^2]

Note that this point is done by the intersection between z=zf and the circle of radius |zf|+eps. In the following, you can see an animation of the system:

Manipulate[ Module[{sol, y, z, t}, sol = First@DSolve[{y'[t] == 1/(2 zf) y[t] - u z[t], z'[t] == -1 + 1/zf z[t] + u y[t], y[ti] == yi, z[ti] == zi}, {y[t], z[t]}, t]; Show[p1, p2, Graphics[{Red, Line[{{0, -1}, {0, 1}}]}], Graphics[{DotDashed, Red, Thickness[.006], Line[{{-1, zf}, {1, zf}}]}], ParametricPlot[{y[t] /. sol, z[t] /. sol}, {t, 0, tend}, PlotStyle -> Thickness[.004]]]], {{tend, .1, "t"}, .01, 20, .1, Appearance -> "Labeled"}, {{u, 10, "u"}, -300, 300, .5, Appearance -> "Labeled"}, {{zi, -.75, "zi"}, -1, 1, Appearance -> "Labeled"}, {{yi, 0, "yi"}, -1, 1, Appearance -> "Labeled"}, Initialization :> (deltaA[y_, z_] := 1/(2 zf) y^2 - z + 1/zf z^2; circ[y_, z_] := y^2 + z^2; p2 = ContourPlot[{circ[y, z] == (-zf + eps)^2}, {y, -1, 1}, {z, -1, 1}, ContourStyle -> Yellow, GridLines -> Automatic, Frame -> True, FrameLabel -> {"y", "z"}, RotateLabel -> False, LabelStyle -> {FontSize -> 20}]; p1 = ContourPlot[{deltaA[y, z] == 0}, {y, -1, 1}, {z, -1, 1}, ContourStyle -> Green, GridLines -> Automatic, Frame -> True, FrameLabel -> {"y", "z"}, RotateLabel -> False, LabelStyle -> {FontSize -> 20}];)]
Improve this question
$\endgroup$
6
  • $\begingroup$ y depend on t, the time, AND u, a free parameter ... $\endgroup$ Commented Jul 16, 2014 at 17:44
  • $\begingroup$ I'have just added a Manipulate to show the behaviour of the system. Anyway, maybe I've not understand your question...the definitions of y[t],z[t] come from DSolve. They are analytic.. $\endgroup$ Commented Jul 16, 2014 at 18:18
  • $\begingroup$ Maybe I am the one not understanding :) Your code doesn't compute, brackets mismatch. $\endgroup$ Commented Jul 16, 2014 at 18:21
  • $\begingroup$ ??? :) Try again...If your problem was in DSolve, I'have rewrited it $\endgroup$ Commented Jul 16, 2014 at 18:29
  • $\begingroup$ I'm afraid that it doesn't work. $\endgroup$ Commented Jul 16, 2014 at 18:38

1 Answer 1

Reset to default
2
$\begingroup$

Here is a stab at what I think you are asking:

 ti = 0; yi = 0; zi = -.75; zf = -.5; eps = .01; sol = {y[t], z[t]} /. First@FullSimplify[DSolve[{ y'[t] == 1/(2 zf) y[t] - u z[t], z'[t] == -1 + 1/zf z[t] + u y[t], y[ti] == yi, z[ti] == zi}, {y[t], z[t]}, t]]; target = {Sqrt[-2*(zf)*eps + eps^2], zf};

brute force discretize the solution and find the minimum distance to the target point (for such a highly nonlinear function this is much faster than NMinimise, and guarantees finding a global min, within the discretization approximation of course )

 dis[u0_] := Min@Table[Norm[ (sol /. u -> u0) - target ] , {t, 0, 10, .05}];

note the range and increment on t here are important tuning parameters to play with.

 Plot[dis[u], {u, -1, 1}]

enter image description here

you come very close to your target point around 0.21.

Armed with a good guess now we can use FindMinimum:

 FindMinimum[Norm[sol - target], {u, .21 } , {t, 1}]

{3.23748*10^-9, {u -> 0.230625, t -> 1.62219}}

 ParametricPlot[ Table[Chop[sol /. u -> u0 ], {u0, {.1, .230625, .3, 1, 5}}] , {t, 0, 10}, Epilog -> Point[target], PlotRange -> All, AspectRatio -> 1]

enter image description here

Improve this answer
$\endgroup$
2
  • $\begingroup$ thank you very much for this very useful and didactic answer. I have a question: why do you need chop? Furthermore I try to make a numerical procedure, using a Do and the function FindAllCrossing. But it doesn't work and I don't know why. Can I post an answer with my wrong-code to understand my error? $\endgroup$ Commented Jul 17, 2014 at 0:04
  • $\begingroup$ Chop removes a numerically negligible complex part. Hitting a point in 2d isn't formulated in terms of 'zero crossings' so i cant see how you could use the methods in the link. You need to find minimum distance to the point and if you need to automate it you need to put in a check that the minimum you find is really zero. $\endgroup$ Commented Jul 17, 2014 at 12:49

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.