# Monte Carlos Method: Particular Equation Style Required

• Jul 29th 2011, 01:48 AM
johnact121
Monte Carlos Method: Particular Equation Style Required
I require a numerical program to generate a random equation for me with the requirements listed below. This curve is for a velocity profile of an object, must have positive velocity at all times and has limits on the maximum acceleration and deceleration. I desire to create a slightly different curve each time I generate one, still meeting all the requirements.

The requirements for the function, v(t) for 0 < t < t_max, is:

• global max < b
• global min = 0
• v(0) = v(t_max) = 0
• v ' (t)_max <= a
• v ' (t)_min >= c
• integral ( v(t) ) = d

The following two requirements, as seen on the image, add no new information to the problem:

• v ' (0) <= a
• v ' (t_max) >= c

An image with this information is shown below:

Any help on this problem would be appreciated. If it is very lengthy to explain how to do kindly point me in the right direction as I have some maths textbooks.

For the numerical program I plan to use either MATLAB or Mathematica.

Thanks again!
• Jul 29th 2011, 05:57 AM
Ackbeet
Re: Monte Carlos Method: Particular Equation Style Required
An interesting problem. I have a few comments:

1. $\displaystyle \int_{0}^{t_{\text{max}}}v(t)\,dt = d$ is superfluous if $\displaystyle \min_{0\le t\le t_{\text{max}}}v(t)=0.$ Unless, as I suspect is the case, $\displaystyle d$ is a number you know ahead-of-time. Certainly, if $\displaystyle v(t)\ge 0$ on the interval in question, then $\displaystyle \int_{0}^{t_{\text{max}}}v(t)\,dt$ is the distance traveled, whatever that is.

2. Depending on the maximum slopes of the velocity, a and c, and the required distance traveled, d, you may not be able to produce a function satisfying all those requirements. In particular, if the required distance traveled, d, is large compared with a and c, you will not be able to exhibit the required function.

3. My solution: From the origin, start rising at the maximum slope a. Backwards from the endpoint, descend at the maximum slope c. Then fiddle with what's in-between those two lines to make the area equal to d. You could try to make the area under the extreme lines equal to half of d, and account for the other half of d by various shenanigans that allow you to change things and still keep the same area.
• Jul 29th 2011, 08:44 AM
SpringFan25
Re: Monte Carlos Method: Particular Equation Style Required
Notice that if you have n functions which satisfy the constraints, then any normalised linear combinations of them will also match the constraints

ie, if you have a collection of n functions $\displaystyle f_i(t)$ which satisfy the contraints, then the following does as well:
$\displaystyle g(t) = \sum_{i=1}^{n} w_i f_i(t)$

with the resitrictions $\displaystyle 0 \leq w_i \leq 1$ and $\displaystyle \sum w_i = 1$

So, you can generate your fluctuations by randomising the weights $\displaystyle w_i$. My suggested algorithm for you would be:
(1) Manually find a few (say, 10) functions which satisfy your contraints using any method (eg, trial and error). This should be feasible enough, and youll only have to do it once.
(2) use a random number generator to assign a weight to each function you found in (1)
(3) then plot the function g(t) corresponding to those randomly chosen weights

to create more realisations, you can generate a new set of random numbers for step (2).

.

Notes

• To Generate random weights you can generate $\displaystyle n$ random numbers $\displaystyle X_i \geq 0$ by any method (all maths programs shold be able to do this), then use $\displaystyle w_i=\frac{X_i}{\sum_{i=1}^{n} X_i}$
• a particularly simple implementation is with n=2 since you could then use a single random number generator to get $\displaystyle w_1$ betgween 0 and 1, then use $\displaystyle w_2 = 1- w_1$
• Obviously the n functions you find will need to be dsignificantly different if you want to see significant fluctions
• the types of paths you observe will be heavily dependent on the shapes of the functions you choose, which may be unrealistic for your purpose. i would try to choose functions that look different to each other.

PS: Although i started this post by saying "notice that", you should check the result in case i made a mistake.
• Jul 29th 2011, 06:28 PM
johnact121
Re: Monte Carlos Method: Particular Equation Style Required
Thank you Ackbeet and SpringFan25 for your prompt replies.

Ackbeet:
Quote:

Unless, as I suspect is the case, http://latex.codecogs.com/png.latex?d is a number you know ahead-of-time.
Yes, d is known and the values a, b, c are also known. Yes, since this is a velocity profile the derivative is acceleration and the integration is the displacement of the object.

For example, let us choose the values of:

• distance (integral) = d = 74m
• max accn = a = 2m/s^2
• max dccn = c = -4m/s^2
• t_max = 20 seconds.

One function which achieves this is a max acceleration for 2 seconds, holding 4m/s for 17 seconds, and max deceleration for 1 second. This generates a function in a trapezium shape with no negative values, meeting all the requirements.

I wish, however, to find a function (continuous) which meets these conditions (notice I have removed the constraint for the global maximum (b) as this is taken care of by the maximum rates of change).

SpringFan25:
Quote:

Notice that if you have n functions which satisfy the constraints, then any normalised linear combinations of them will also match the constraints
Yes, your idea of normalised linear combinations is a good idea. My first question would be: how to I generate the continuous functions in the first place?

I have the following ideas so far:

• A numerical solver to solve y = a.x^(n) + b.x^(n-1) + ... &c. with the various constrains put into the numerical program.
• Analytical method to solve for the equation.
• Simple constant acceleration equations for trapezium shaped functions (not continuous but meet the conditions)
• Others?
• Jul 30th 2011, 03:18 AM
SpringFan25
Re: Monte Carlos Method: Particular Equation Style Required
Quote:

Simple constant acceleration equations for trapezium shaped functions (not continuous but meet the conditions)
Minor point: I think that trapezium functions are continuous but not differentiable. I agree that they meet the conditions as written in post #1.

Major point: if you are going to allow non-differentiable functions then the problem is dramatically simplified. You can just start/finish your graph with two short lines with slope a and c respectively. For simplicity, choose the length of those lines so that the area underneath them is 1.

Now all you need to do is choose some random points in the middle and connect them with straight lines. Scale up/down these middle points by a constant factor so that the area under them is d-1. This will bring the total area under the graph to $\displaystyle d$ as required.

If re-scaling causes the global maximum constraint to be breached, discard your middle points points and randomly choose some more. Alternatively reduce the highest point and then repeat the scaling factor calculation.

Edit: I think ackbeet already suggested the above in post #2.

Many combinations will be discarded but a computer could make thousands/millions of attempts in a very short time.

Poosible Smooth functions
Notice that there is no constraint on v' except at t_max and t_min. Trigonometric functions (specifically, combinations of sine waves) look promising if you run out of simpler options (polynomials etc):

• they can be combined to make many different shapes.
• Because they are periodic they can be easily parameterized to meet the slope restrictions at t_min and t_max.
• An upper bound on the global maximum is also easy to compute. (the global maximum cant be more than the sum of the amplitudes of your waves)
• The integrals should also be tractable as long as you confine yourself to simple functions, eg combinations of sine waves with different periods/phases/amplitudes. The integral calculations will be slightly easier if you chose half-periods that fit into the period (t_max - t_min) a whole number of times.
• One minor challenge will be ensuring that the global minimum is not breached but this can be overcome by trial and error.

PS: I expect the above is tractable but i haven't actually tried it myself! I wont look at it further unless you actually need smooth functions.
• Jul 30th 2011, 08:56 AM
SpringFan25
Re: Monte Carlos Method: Particular Equation Style Required
it was too much fun to ignore so I've attached an excel implementation of the smooth trig functions approach i suggested above.

The main changes from my above post are
• t_min is assumed 0
• i was wrong to say that there was no restriction on v'() except at t_min and t_max. i had misread your first post.
• All sine waves complete an odd number of half periods in the time interval. This makes all graphs symettric so v(0) = v(tmax) =0. We also have v'(0) = -v(t_max).
• To simplify calculations all sine waves have a phase shift of 0 or $\displaystyle \pi$
• Ive used a stronger gradient constraint than you need. I used |v'(t)| < min(|a|,|c|). This is useful because an upper bound on |v'(t)| can be easily calculated from the amplitude and period of the sine waves (derivation below)
• There is a big "calculate me" button to generate a new random graph. it requires macros enabled to run. I dont accept responsibility if this breaks your computer/life etc.

Upper bound on |v'(t)|
I'll derive this. Suppose v(t) is the sum of $\displaystyle n$ sine waves with ammplitude $\displaystyle A_i$, period $\displaystyle B_i$, and phase-shift $\displaystyle C_i$ as follows:

$\displaystyle f_i(t) = A_i \sin(\frac{2 \pi}{B_i} t + C_i)$
$\displaystyle v(t) = \sum f_i(t)$

Differentiate the first equation. This gives:
$\displaystyle f'_i(t) = \frac{2 \pi A_i}{B_i} \cos (2 \pi B_i t + C_i)$

Since $\displaystyle |cos(x)| \leq 1$ it follows that
$\displaystyle |f'_i(t)| \leq \frac{2 \pi A_i}{B_i}$

And we can sum as many sine waves as we want, adding up the upper bounds for each. This is handy since v(t) is the sum of the sine waves:

$\displaystyle |v'(t)| \leq \sum |f'_i(t)| \leq 2 \pi \sum \frac{A_i}{B_i}$
(the left hand inequality is valid. Think "triangle inequality" with lots of variables.)

How does this help?

Note that $\displaystyle |v'(t)| \leq |a| \implies v'(t) \leq a ~~ \text{ since } a > 0$
and, $\displaystyle |v'(t)| \leq |c| \implies v'(t) \geq c ~~ \text{ since } c < 0$

Combine the above pair of results:

$\displaystyle |v'(t)| \leq |a| \text{ and } |v'(t)| \leq |c| \implies c \leq v'(t) \leq a$

it follows that
$\displaystyle |v'(t)| \leq min(|a|,|b|) \implies c \leq v'(t) \leq a$
So testing that $\displaystyle |v'(t)| \leq min(|a|,|c|)$ will ensure that the acceleration constraints are met.

Instead of evaluating |v'(t)| everywhere, It's sufficient to test the upper bound of |v'(t)| that was derived in the previous part of this post. ie the test is:

$\displaystyle 2 \pi \sum \frac{A_i}{B_i} \leq \min{(|a|,|c|)}$

Meeting all the constraints
for the rest of this post assume that v(t) is constructed as follows:

$\displaystyle v(t) = \sum f_i(t)$
$\displaystyle f_i(t) = A_i \sin(\frac{2 \pi}{B_i} t + C_i) ~~~ \text{ with }A_i,B_i \geq 0$

v(0) = v(tmax) = 0
This is the easiest one. Consider time 0:
$\displaystyle f_i(0) = A_i \sin(\frac{2 \pi}{B_i} \cdot 0 + C_i)$
$\displaystyle f_i(0) = A_i \sin(C_i)$

so a sufficient condition for $\displaystyle f_i(0)$ is $\displaystyle \sin(C_i)=0$

so the phase shift should be $\displaystyle 0, \pi, 2\pi, 3\pi$. Since the period is $\displaystyle 2 \pi$ anyway, so only unique phase shifts are 0 and $\displaystyle \pi$.

What about $\displaystyle f_i (t_{max})$?
The sin function returns to its starting point every half-period. So its suffiient to ensure that the interval $\displaystyle [0,t_{max}]$ has a whole number of half periods in it. Then we can use $\displaystyle f_i(t_{max}) = f_i(0) = 0$.

(in fact, we will want an odd number of half-periods later on.).

c < v'(t) < a
I already showed that a sufficient condition is
$\displaystyle 2 \pi \sum \frac{A_i}{B_i} \leq \min{(|a|,|c|)}$

Total displacement = d
This is the tricky one. It's very convenient to restrict the period of each sine wave so that an odd number of half periods fit in the interval $\displaystyle [0,t_{max}]$. This means that the number of completed oscilations will always be 1.5, 2.5, 3.5, 4.5 etc. This is important because the integral of a full oscillation is always zero.

So we only need to compute the integral of the remaining half oscillation of each wave. We already established the the phase shift $\displaystyle C_i$ is 0 or $\displaystyle \pi$. For the case where the phase shift is 0, This area is equal to:

$\displaystyle A_i \int_0^{\pi} \sin(B_i t + 0) dt = A_i \left[ \frac{-\cos(B_i t)}{B_i}\right]_0^{\pi}$

For the case where the phase shift is $\displaystyle \pi$ the half-period area is equal to
$\displaystyle A_i \int_{0}^{\pi} \sin(B_i t + \pi) dt$

But remember that the sin function has this property: $\displaystyle sin(x + \pi) = -sin(x)$

So for the phase shifted case the area is just -1 times the area we just calculated:
$\displaystyle A_i \int_{0}^{\pi} \sin(B_i t + \pi) = -A_i \int_{0}^{pi} \sin(B_i t) dt = A_i \left[ \frac{\cos(B_i t)}{B_i}\right]_0^{\pi}$

Finally, a single formula covering both cases:

$\displaystyle \int_{0}^{\pi} f_i dt = A_i \left[ \frac{-\cos(B_i t)}{B_i}\right]_0^{\pi} \cdot (-1)^{\frac{C_i}{\pi}}$
(if you dont udnerstand what the -1 factor is doing at the end, it evaluates to 0 if there is no phase shift, and -1 if there is a phase shift)

Hence the area under the graph $\displaystyle \Omega$ is:
$\displaystyle \Omega = \int_0^{t_{max}} v(t) dt = \sum_i \left( \int_0^{t_{max}} f_i(t) dt \right)= \sum_i \left( A_i \left[ \frac{-\cos(B_i t)}{B_i}\right]_0^{\pi} \cdot (-1)^{\frac{C_i}{\pi}} \right)$

We require the area under the graph to be exactly d. This is easy (compared to all that integration!). Simply set new coefficients:

$\displaystyle \tilde{A}_i = A_i \frac{d}{\Omega}$

This scales up the entire function so that the total area will be d.

Obviously, this scale up/down may brake the acceleration constraints, since the slopes will also be scaled. The attached workbook checks this, and tries a new set of $\displaystyle A_i, B_i$ if there is a breach.

Algorithm
1) Randomly choose $\displaystyle A_i, B_i.C_i$ subject to the rules derived above, ie:
• $\displaystyle C_i$ = 0 or $\displaystyle \pi$
• $\displaystyle B_i = \frac{t_{max}}{k+0.5} ~~~ : ~ k \in \mathbb{N}$

2) Calcualte adjusted $\displaystyle \tilde{A}_i = A_i \frac{d}{\Omega}$

3) Test the adjusted coefficients meet the accelleration constraints. If they dont, discard and return to step 1. The test is $\displaystyle 2 \pi \sum \frac{A_i}{B_i} \leq \min{(|a|,|c|)}$

The attached workbook does the above. I've chosen ranges on $\displaystyle A_i,B_i$ that produce reasonably varied output.

Sample output
here is some sample output
http://i52.tinypic.com/34o7sit.jpg
• Aug 2nd 2011, 11:41 PM
johnact121
Re: Monte Carlos Method: Particular Equation Style Required
SpringFan25,

Thank you very much for your detailed reply. I will need some time to go through it and understand it.

One thing to point out, it seems that all the functions generated are even and thus symmetrical. That constraint is not required - is there a way to get around that?
• Aug 3rd 2011, 03:45 AM
SpringFan25
Re: Monte Carlos Method: Particular Equation Style Required
The symmetry is convenient because it makes sure that the velocity graph returns to its starting value, ie v(0) = v(tmax) =0.

You can get rid of the symmetry by adding a further asymmetric function to v(t) which also satisfies v(0) = t(tmax) = 0.

One way is to add another pair of sine waves which are 180 degrees out of phase, and with one having double the frequency of the other:

examples of 2 sine waves and their sum are below to illustrate:
plot y&#61;sin&#40;x&#41;,y&#61;0.5sin&#40;2x-pi&#41; on 0,pi - Wolfram|Alpha
plot y&#61;sin&#40;x&#41; &#43; 0.5sin&#40;2x-pi&#41; on 0,pi - Wolfram|Alpha

However these functions are not ideal as the low gradient at the start means that small fluctuations from any other waves you put on top can easily make the overall value of v(t) negative at low values of t.