Approximate Density

Processing math: 0%

Approximate Density

$\def \half {\frac{1}{2}} \def \hieruit {\quad \Longrightarrow \quad} \def \slechts {\quad \Longleftrightarrow \quad} \def \OF {\quad \mbox{or} \quad} \def \MET {\quad \mbox{with} \quad}$ The reverse question. Suppose that the sense

$P(x)$ is given, how can we find then the associated density of points

$\{x_k\}$ ? A crude approximation is to take the value

$P(x_k)$ at the position

$x_k$ and construct the next point by demanding that the rectangle with height

$P(x_k)$ and width

$(x_{k+1}-x_k)$ has an area equal to one:

$P(x_k) (x_{k+1}-x_k) = 1 \hieruit x_{k+1} = x_k + \frac{1}{P(x_k)}$ Herewith it is assumed that

$P(x)$ is positive everywhere. But, as we have seen with the guitar example in the previous section, it is also possible that

$P(x)$ is negative everywhere, giving rise to a monotone decreasing sequence. In that case, the x-axis shall be traversed in the opposite direction, giving:

$P(x_k) (x_k-x_{k-1}) = 1 \hieruit x_{k-1} = x_k - \frac{1}{P(x_k)}$ Anyway, a much better job can be done as follows. Consider a linear approximation of the sense in the neighbourhood of

$x_k$ :

$P_{k+1} = P(x_k) + P'(x_k) (x_{k+1} - x_k)$ Where it is noticed that the estimate

$P_{k+1}$ may be almost, but not exactly equal to

$P(x_{k+1})$ . The density may be approximated with help of the trapezium rule as well:

$\half \left[ P_{k+1} + P(x_k) \right] \left[x_{k+1} - x_k\right] = 1$ Substitution of the first formula into the second gives:

$\half \left[ P(x_k) + P'(x_k) (x_{k+1} - x_k) + P(x_k) \right] \left[x_{k+1} - x_k\right] = 1$

$\hieruit \left[ \half P'(x_k) \right] (x_{k+1} - x_k)^2 + P(x_k)(x_{k+1} - x_k) - 1 = 0$ A quadratic equation. Solve for

$(x_{k+1} - x_k)$ :

$x_{k+1} - x_k = \frac{ - P(x_k) \pm \sqrt{P^2(x_k) + 2 P'(x_k)}}{P'(x_k)}$ Take as an example the linear density:

$P(x) = \frac{1}{\Delta^2}\left(x+\frac{\Delta}{2}\right) \hieruit P'(x)= \frac{1}{\Delta^2} \hieruit \\ x_{k+1} - x_k = \frac{ - (x_k+\Delta/2)/\Delta^2 \pm \sqrt{\left[(x_k+\Delta/2)/\Delta^2\right]^2+ 2/\Delta^2}}{1/\Delta^2} \hieruit \\ x_{k+1} = \pm \sqrt{(x_k+\Delta/2)^2 + 2\Delta^2}-\Delta/2$ Compare with a previously found result to establish that the sign of

$\pm$ must be plus:

$x_k = \left( \sqrt{1/4 + 2 k}-1/2 \right) \Delta \hieruit (x_k+\Delta/2)^2 = \Delta^2(1/4 + 2 k)^2 \hieruit \\ x_{k+1} = \pm \Delta\sqrt{(1/4 + 2 k) + 2}-\Delta/2 = \left( \sqrt{1/4 + 2 (k+1)}-1/2 \right) \Delta$ It should be no surprise that the outcome for a linear density is not approximate but exact.
Anyway, quite in general, we still have:

$x_{k+1} - x_k = \frac{\sqrt{P^2(x_k) + 2 P'(x_k)} - P(x_k)}{P'(x_k)}$ The numerator is of the form

$(a-b)$ . Now multiply with

$(a+b)/(a+b)$ , leading to values of

$(a^2-b^2)/(a+b)$ which are likely more stable numerically:

$x_{k+1} - x_k = \frac{2 P'(x_k)/P'(x_k)}{\sqrt{P^2(x_k) + 2 P'(x_k)} + P(x_k)} \\ \hieruit x_{k+1} = x_k + \frac{1}{\sqrt{\left[P(x_k)/2\right]^2 + P'(x_k)/2} + P(x_k)/2}$ The plus sign (out of

$\pm\sqrt{}$ ) also comes from the fact that, for

$P'(x) = 0$ , the next point has to be found by the (no longer) crude approximation:

$x_{k+1} = x_k + 1/P(x_k)$ . It is noticed, in addition, that the expression under the square root should not become negative. Meaning that, if

$P'(x_k) < 0$ :

$| P'(x_k) | < \half \left[ P(x_k) \right]^2$ It can be argued, though, that the above condition must also hold for positive values of

$P'(x_k)$ . The reason is that the sequence

$x_k$ could be traversed in the opposite direction as well:

$\left. \begin{array}{l} f_{k-1} = P(x_k) + P'(x_k) (x_{k-1} - x_k) \\ \half \left[ P(x_k) + f_{k-1} \right] \left[x_{k-1} - x_k\right] = -1 \end{array} \right\} \hieruit$

$\left[ \half P'(x_k) \right] (x_{k-1} - x_k)^2 + P(x_k)(x_{k-1} - x_k) + 1 = 0$

$\left[ \half P'(x_k) \right] (x_{k+1} - x_k)^2 + P(x_k)(x_{k+1} - x_k) - 1 = 0$ Thus

$(k+1)$ changes into

$(k-1)$ means that

$+1$ changes into

$-1$ . And that's the only difference between traversing in a positive or in a negative direction. A bit of further thinking reveals then that the negative end-result must be:

$x_{k-1} = x_k + \frac{1}{\sqrt{\left[P(x_k)/2\right]^2 - P'(x_k)/2} + P(x_k)/2}$ And therefore, again, the abovementioned condition must hold:

$| P'(x_k) | < \half \left[ P(x_k) \right]^2$ Take for example the linear density, where

$P(x) = (x/\Delta + \half)/\Delta$ ,

$P'(x) = 1/\Delta^2$ . Then our condition translates into:

$| P'(x_k) | < \half \left[ P(x_k) \right]^2 \slechts 1/\Delta^2 < \half \left[ (x_k/\Delta + \half)/\Delta \right]^2 \slechts$

$\left[ (x_k/\Delta + \half) \right]^2 > 2 \slechts x_k/\Delta > - \half + \sqrt{2} \OF x_k/\Delta < - \half - \sqrt{2}$ Only the non-negative solution must be kept. We conclude that the algoritm for constructing the set

$\{x_k\}$ from the linear sense

$P(x) = (x/\Delta + \half)/\Delta$ will certainly work for

$x_k > (\sqrt{2} - 1/2) \Delta \approx 0.9142135 \Delta$ , let's say

$x_k > \Delta$ , which is certainly the case for

$k > 0$ .
Take as an example the sense that is an approximation for the longitudinal wave with small amplitude, eventually without the latter restriction:

$P(x_k) = \frac{1}{\Delta}\left[1 - A\frac{2\pi}{\lambda}\cos\left(\frac{2\pi}{\lambda} x_k\right)\right] \hieruit \\ P'(x_k) = \frac{A}{\Delta}\left(\frac{2\pi}{\lambda}\right)^2\sin\left(\frac{2\pi}{\lambda} x_k\right)$ From which the iteration scheme follows. The result is:

Accompanying software for making the pictures is disclosed as well.