Services & Resources / Wolfram Forums / MathGroup Archive
-----

MathGroup Archive 2009

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Length distribution of random secants on a unit

  • To: mathgroup at smc.vnet.net
  • Subject: [mg95737] Re: [mg95712] Length distribution of random secants on a unit
  • From: Bob Hanlon <hanlonr at cox.net>
  • Date: Sun, 25 Jan 2009 21:50:29 -0500 (EST)
  • Reply-to: hanlonr at cox.net

When the two points lie on the same side (1/4 of the time)

With[{d = RandomReal[]},
 Show[
  RegionPlot[Tooltip[Abs[x1 - x2] <= d],
   {x1, 0, 1}, {x2, 0, 1}],
  Plot[{x1 - d, d + x1}, {x1, 0, 1},
   PlotStyle -> {Darker[Blue], Red},
   PlotRange -> {0, 1}]]]

Integrating over the valid region to obtain the CDF

Assuming[{0 <= d <= 1}, 
 FullSimplify[Integrate[1, {x, 0, 1}, {y, Max[0, x - d], Min[1, d + x]}]]]

Piecewise[{{3/4, 2*d == 1}, {1, d == 1}, {(-(-2 + d))*d, 
   0 < d < 1/2 || 1/2 < d < 1}}, 0]

The CDF for this case is then

dist1CDF[d_] = Piecewise[{{1, d >= 1}, {((2 - d))*d, 
        0 < d < 1}}, 0]

Piecewise[{{1, d >= 1}, {(2 - d)*d, 0 < d < 1}}, 0]

Similarly, when the two points lie on adjacent sides (1/2 of the time)

With[{d = Sqrt[2]*RandomReal[]},
 Show[
  RegionPlot[Tooltip[Sqrt[x^2 + y^2] < d],
   {x, 0, 1}, {y, 0, 1}],
  Plot[Sqrt[d^2 - x^2], {x, 0, 1},
   PlotStyle -> Red,
   PlotRange -> {0, 1}]
  ]]

dist2CDF[d_] = FullSimplify[Piecewise[{
    {1, d >= Sqrt[2]},
    {Integrate[1, {x, 0, d}, {y, 0, Sqrt[d^2 - x^2]},
      Assumptions -> {0 < d <= 1}], 0 < d <= 1},
    {Integrate[1, {x, 0, 1}, {y, 0, Min[1, Sqrt[d^2 - x^2]]},
      Assumptions -> {1 < d < Sqrt[2]}], 1 < d < Sqrt[2]}}]]

Piecewise[{{1, d >= Sqrt[2]}, {(d^2*Pi)/4, 
   Inequality[0, Less, d, LessEqual, 1]}, 
     {Sqrt[-1 + d^2] + (1/2)*d^2*(ArcCsc[d] - ArcCsc[d/Sqrt[-1 + d^2]]), 
   1 < d < Sqrt[2]}}, 0]

When the two points lie on opposite sides (1/4 of the time)

With[{d = (Sqrt[2] - 1)*RandomReal[] + 1},
 Show[
  RegionPlot[Tooltip[Sqrt[(x2 - x1)^2 + 1] < d],
   {x1, 0, 1}, {x2, 0, 1}],
  Plot[{x1 - Sqrt[d^2 - 1], Sqrt[d^2 - 1] + x1}, {x1, 0, 1},
   PlotStyle -> {Darker[Blue], Red},
   PlotRange -> {0, 1}]]]

Assuming[{1 <= d <= Sqrt[2]}, FullSimplify[Integrate[1, {x1, 0, 1},
   {x2, Max[0, x1 - Sqrt[d^2 - 1]], Min[1, Sqrt[d^2 - 1] + x1]}]]]

Piecewise[{{3/4, 2*d == Sqrt[5]}, {1, 
   d == Sqrt[2]}, {1 - d^2 + 2*Sqrt[-1 + d^2], 
       1 < d < Sqrt[5]/2 || Sqrt[5]/2 < d < Sqrt[2]}}, 0]

dist3CDF[d_] = Piecewise[{
   {1, d >= Sqrt[2]},
   {1 - d^2 + 2*Sqrt[-1 + d^2], 1 < d < Sqrt[2]}}]

Piecewise[{{1, d >= Sqrt[2]}, {1 - d^2 + 2*Sqrt[-1 + d^2], 
   1 < d < Sqrt[2]}}, 0]

The CDF for the distribution is then

distCDF[d_] = FullSimplify[PiecewiseExpand[
   dist1CDF[d]/4  + dist2CDF[d]/2  + dist3CDF[d]/4]]

Piecewise[{{1, d >= Sqrt[2]}, {(2 + Pi)/8, 
   d == 1}, {(1/8)*d*(4 + d*(-2 + Pi)), 0 < d < 1}, 
     {1/2 + 
    Sqrt[-1 + d^2] + (1/4)*d^2*(-1 + ArcCsc[d] - ArcCsc[d/Sqrt[-1 + d^2]]), 
   1 < d < Sqrt[2]}}, 0]

Comparing with the experimental data

len = Norm[(First[#] - Last[#])] &;
corners = {{0, 0}, {1, 0}, {1, 1}, {0, 1}};
dir = {{1, 0}, {0, 1}, {-1, 0}, {0, -1}};

p[t_] := Block[{n, r}, n = Mod[IntegerPart[t], 4];
  r = FractionalPart[t];
  corners[[n + 1]] + r dir[[n + 1]]]

Show[Histogram[Table[len[{p[RandomReal[{0, 4}]],
     p[RandomReal[{0, 4}]]}], {100000}],
  Automatic, "ProbabilityDensity"],
 Plot[distPDF[d], {d, 0, Sqrt[2]},
  PlotStyle -> {AbsoluteThickness[2], Red},
  PlotRange -> {0, 2.4}]]


Bob Hanlon

---- andreas.kohlmajer at gmx.de wrote: 

=============
I need to work with the length distribution of random secants (of two
random points on the perimeter) on a unit square. It's easy to
generate some random data and a histogram. I used the following code
(Mathematica 7.0):

len = Norm[(First[#] - Last[#])] &;
corners = {{0, 0}, {1, 0}, {1, 1}, {0, 1}};
dir = {{1, 0}, {0, 1}, {-1, 0}, {0, -1}};

p[t_] := Block[{n, r},
  n = Mod[IntegerPart[t], 4];
  r = FractionalPart[t];
  corners[[n + 1]] + r dir[[n + 1]]
  ]

Histogram[
 Table[len[{p[RandomReal[{0, 4}]], p[RandomReal[{0, 4}]]}], {100000}],
  PlotRange -> All]

The histogram shows a small increase close to 1, a big peak at 1 and
some kind of exponential decay to Sqrt[2] (= maximum).

Does anybody know how to calculate this distribution exactly? What
about moving from a unit square to a random rectangle or a random
polygon? Thanks!



  • Prev by Date: Re: Re: Re: Mathematica Animation Drives Me Crazy!
  • Next by Date: Re: iterative method?
  • Previous by thread: VectorPlot3D and VectorStyle
  • Next by thread: Re: Length distribution of random secants on a unit