I am a beginner in mathematica and I am facing a problem that is way bigger than my knowledge..
I am tryng to simulate a semplificated version of Brunel Effect, wich is a Laser-Plasma interaction effect, with an electrostatic Dawson's model.
The theory is simple: the plasma, with a sharp-step density shape, is hit by a sinusoidal laser impulse, some electrons will be pulled outside of the plasma and reach the vacuum region in which they are accelerated and then they re-enter into the plasma with more energy than before.
Mathematically the problem is also quite simple: the differential equations wich governs the electrons motion are different inside and outside the plasma, but for the consistency of the fluid model the electrons trajectories should never overlap and for avoiding this i have to "switch" the label of the trajectories when overlap occur. In a few words, i have to insert some checks while integrating the equations that recognize when an overlap occur and then restart the integration with the velocities swapped from that point, i.e. i have to change the initial conditions of the 2 electrons and then continue the integration till the next overlap.
The problem is that check should act inside the NDSolve procedure and while it is doing the integration.
Here the idea of the code i had been able to develop till now: Initialitation of the Field:
ClearAll;
W = 10^6;
tau = 10*Pi/W;
c = 299792458;
Am = 0.8*W*c;
Ed[t_]:=If[t <((2Pi)/W),Sin[Pi*t/(2*tau)]*Sin[Pi*t/(2*tau)]*10*Am*Sin[W*t], Am*Sin[W*t]];
Plot[Ed[t], {t, 0, (2 Pi/W)*2}]
Define some constants and functions:
t =.
it = 0;
imax = 100;
check = 0;
check2 = 0;
aux = 0;
x0max = 4;
tmax = (2*Pi/W)*7;
V[x0_,t_]:=D[Z[x0,t],t];
x[x0_,t_]:=x0+Z[x0,t];
Initialitation of the eq.:
ndsdata = First@NDSolve`ProcessEquations[{D[V[x0, t], t] ==
If[x0 + Z[x0, t] > 0, -10*W^2*Z[x0, t] - Ed[t],
10*W^2*x0 - Ed[t]], V[x0, 0] == 0, Z[x0, 0] == 0}, Z[x0, t], {x0, 0, x0max},
{t, 0, tmax}, Method -> "ExplicitRungeKutta", MaxSteps -> Infinity]
Where Ed[t] is a Sinusoidal field, W is the pulse frequency, x0 the initial position of the electron, Z[x0,t] the electron displacement (the position of the electron at time t is x=x0+Z[x0,t] and V[x0,t] is the time derivate of Z[x0,t], which is the velocity of the electron
Iteration of the Solver:
For[i = 0, i < imax, i++, {it = it + 1/100000;
NDSolve`Iterate[ndsdata, it];
auxsol = NDSolve`ProcessSolutions[ndsdata];
For[j = 0, j < x0max, j++,
If[Evaluate[x[j, it] /. auxsol] >
Evaluate[x[j + 1, it] /. auxsol],
ndsdata = NDSolve`Reinitialize[ndsdata, {V[j, it] == V[j + 1, it], Z[j, it] ==
Z[j + 1, it]}]];If[it < ((2 Pi/W))*7, imax = imax + 1, i = imax]]
And the process the solution and Plot:
sol = NDSolve`ProcessSolutions[ndsdata];
this code does not work, but if i try to integrate without the swap i am able to plot the (wrong) trajectories...
the code gives me no error, but it seems that the following If cycle neither became true nor false:
If[Evaluate[x[j, it] /. auxsol]>Evaluate[x[j + 1, it] /. auxsol],
{check = 1, V[j, it] == aux,V[j, it] = V[j + 1, it], V[j + 1, it] = aux,ndsdata =NDSolve`Reinitialize[
ndsdata, {V[j, it] == V[j + 1, it], Z[j, it] == Z[j + 1, it]},
check2 = 1]}]
So basically the swap (or what i have tried to do as a swap) never happens.
What i obtain are simply the solution of the integration, that are physically uncorrect:
Plot[Table[Evaluate[x[x0, t] /. sol], {x0, 0, x0max}],{t,0,(2*Pi/W)*7},
PlotRange -> {{0, ((2*Pi/W))*7}, {-700, 300}}]
What i am asking is: How can work inside Ndsolve procedure and change initial conditions of the equations at a certain time, when an event (in this case the overlap) occurs?
Thank you in advance and sorry for the expositive mess..
WhenEvent
documentation. $\endgroup$Z
, it is unlikely that you will receive useful help.x[x0,t]=x0+Z[x0,t]
does not defineZ
, becausex
also is undefined. Please add all this to your question. By the way, use @bbgodfrey in any comments you wish me to see promptly. I should add that you may be making this problem too hard for yourself. $\endgroup$