Differential Equations with Discontinuities

Translate the entire conversation x

I am solving a system of odes with imposed discontinuites on the right hand side. I can do this by calling a solver up to the time of the discontinuity, imposing the jump condition, and then restarting the calculation (as in the attached worksheet). But I want to do this for a series of discontinuities. Since all the jumps are the same magnitude and occur at regular time intervals, it seemed like there should be a way to program it rather than doing it manually, but I could not see how to do it. Can anyone help?

8 REPLIES 8

(To:JohnRudnicki)

Attached is a version using Odesolve (Mathcad 15). For more than two sets of discontinuous jumps you should be able to modify it in an obvious way.

It produces the following sort of results

Alan

(To:AlanStevens)

Thanks for the reply. It seems like that should work but your plot does not show a discontinuity in theta (t) at t = 1000. I did not show the plot of this in my previous post, but it is shown below:

Ultimately, I would like to do this for a staircase of discontinuties (with the same stepsize and time interval). BTW, I used dk=0 even though k is constant so I could change the value in the initial condition for different cases (though I think this is easier to do with Odesolve). I have not been using Odesolve because in looking at similar systems I have had a lot of trouble with convergence, the error "return value does not match...", and instabilities that is sometimes helped by fiddling with different solvers and tolerances.

(To:JohnRudnicki)

A little more like this then, where I've kept k constant and only solved for v and theta (you should be able to extend the method to the others easily enough). I've also added another step:

Alan

(To:AlanStevens)

Thanks. Looks good, but I will have to study a bit to make sure I understand.

(To:JohnRudnicki)

@JohnRudnicki wrote:

Thanks. Looks good, but I will have to study a bit to make sure I understand.

The conditions for the solution of the problem:

$d x d t = F (t, x (t)), x (t 0) = x 0,$

to exist and to be unique is that $F \in C_{t, x}^{0}$ and Lipschitz continuous for $t > t_{0}$ . Then, the solution of the problem is given by:

$x (t) = x 0 + \int t t 0 F (t, x (t)) d t$

Note that for the function $F$ to be integrable in $t > t_{0}$ , it must contain finitely many discontinuities so $x (t)$ remains continuous in $t > t_{0}$ whilst, otherwise, $x^{'} (t)$ may not, in general.

As an example, consider the problem:

$x' (t) = H (t - 1), x (0) = 0, t > 0.$

Its solution is then given by:

$x (t) = (t - 1) H (t - 1)$

which is continuous but not differentiable at $t = 1$ .

Cheers!

(To:AlanStevens)

Hi Alan,

There still seems to be a discrepancy in how the jumps are implemented as shown below. JWRCalc uses my (laborious) method in worksheet Jumps04_jwr. Plot is added to the bottom of your worksheet in Jumps04_as. I could not detect where the problem is.

(To:JohnRudnicki)

I looked at this a little more but could not see how to modify your program to capture the jumps accurately (see bottom of attached worksheet). I do not know that it is possible using odesolve because it interpolates between the solution points. The manual way I did it, the same time value is assigned to the indices N and N+1, but the state (theta) has different values for N and N +1. So it gets the discontinuity exactly. I could not see how to adapt your program to do this by calling the solvers directly. So, unless you have some other suggestions, I guess I am stuck with solving between the discontinuites, implementing the jump condition and repeating. But thanks for your help.

(To:JohnRudnicki)

After having to do some other things for a few weeks, I got back to this. I managed to program it (program is at bottom of the worksheet). It seems to be working correctly, though there is probably a better way to do it. If anyone has suggestions, I would be happy to hear them.

Top Tags