Math Quiz 9

[GeeTwo]

Originally Posted by GeeTwo

What is the average length of all the line segments which can be drawn within the unit circle (1 unit in radius, 2 units in diameter)?

This can actually be solved in closed form with much less intricate calculus than the case of the square, with proper selection of coordinate system.

Some warmup questions to this if you can't figure out where to start:

• What is the average of the squares of the lengths of all the line segments which can be drawn within the unit circle? [taking out the square root simplifies the problem greatly, and helps provide an upper bound to the original question]
• What is the average of the lengths of all the chords of the unit circle? (That is, line segments with both endpoints on the border of the circle)
• What is the average of the lengths of all the line segments between a point on the [border of] a circle of radius R, and a point within the same circle?

OBTW, with proper scaling, Greg's answer for the 100000-gon of area 1 provides an answer to the original problem good to within 1 part per 10,000:

Originally Posted by Greg Woelki

Here is a generalized Monte Carlo simulation ....

And, OBTW2: if you can find the closed form, each of the four questions can be calculated in fewer than 10 keystrokes on a calculator with the following buttons:

0 1 2 3 4 5 6 7 8 9 + - / * π √ =

[GeeTwo]

I have posted this figure which I created to describe variables in my calculation of the length of the average line segment from a point on the edge of a circle to a point on the interior. I used this, or a very similar coordinate system, to solve all the "circle questions" except for the one about the mean square of the segments. I am also confident that the same questions can be answered for the sphere using a similar point of view.

[GeeTwo]

Originally Posted by GeeTwo

I am also confident that the same questions can be answered for the sphere using a similar point of view.

I did the integral for the average length of a line segment within a sphere, and found that the answer is rational! This problem was actually simpler to integrate (though just a bit trickier to set up) than any of the two-dimensional problems. If no one posts anything (yes, even a "hold on, I'm working it!") on the circle or sphere questions by Wednesday at 1800 US-CDT, I'll post my solutions so we can close out these quiz questions.

[AMendenhall]

0.9525383819936485
Sample size: 60 million

[Ether]

Originally Posted by AMendenhall

0.9525383819936485
Sample size: 60 million

Post your code and we'll show you where you went wrong.

[Ether]

Originally Posted by GeeTwo

I did the integral for the average length of a line segment within a sphere, and found that the answer is rational!

The average length of sphere chords (line segments with both endpoints on sphere surface) is also rational.

[GeeTwo]

Originally Posted by Ether

The average length of sphere chords (line segments with both endpoints on sphere surface) is also rational.

Yes, as well as the intermediate case where one endpoint is on the sphere and the other is internal. I don't recall whether the mean square length was rational or not; I'll check this evening.

[Ether]

Originally Posted by GeeTwo

I don't recall whether the mean square length was rational or not.

It is. (mean square length with one endpoint on surface).

[GeeTwo]

Originally Posted by Ether

It is. (mean square length with one endpoint on surface).

I was referring to both endpoints in the volume, but I see that for all three cases (both exterior, one exterior, both interior), the mean length as well as the mean square length are all rational.

Given AM's interest, I will hold off until at least Friday at 6pm to post solutions.

When posting a candidate solution, please show your work, and specify whether you are solving for:

• Circle or sphere?
• Mean length or mean square length?
• Both, one, or none of the endpoints constrained to the boundary?

[Greg Woelki]

A Monte Carlo simulation written in python and running on my laptop isn't a great solution if 1 part in a million accuracy is required, but since I had essentially already written the program, I figured that I might as well give it a shot. For the average length of line segments in the interior of an r=1 circle, I got 0.9054124568 after 4e10 trials. I'm not certain whether or not this meets the accuracy requirement. I also made minor modifications so that I could simulate the chord problem, and I got 1.2732310914 after 7e09 trials. This problem requires far fewer trials to get an accurate result and I believe this solution meets the accuracy spec. I can post my code if anyone is interested.

[Ether]

Originally Posted by Greg Woelki

0.9054124568 after 4e10 trials

Wrong by ~3 ppm.

1.2732310914 after 7e09 trials

Wrong by ~7 ppm.

[Ether]

Originally Posted by Greg Woelki

0.9054124568 after 4e10 trials

1) How long did 4e10 trials take to run?

2) Yes, please post your code.

edit:

3) What's the period of the PRNG you used?

[Greg Woelki]

Originally Posted by Ether

1) How long did 4e10 trials take to run?

About 16 hours. My laptop has an Intel Core i7-6700HQ Quad Core and I ran 4 copies of my program at once instead of implementing multithreading. This seemed to work as desired based on the cpu usage and left me with enough overhead to not impact my normal computer use.

Originally Posted by Ether

2) Yes, please post your code.

Internal points code:

import numpy as np
import time

def genPt(r_cir):
    
    r = r_cir*(np.random.random()**0.5) #Square root correction ensures even probability distribution over area
    theta = 2*np.pi*np.random.random()
    
    x = r*np.cos(theta)
    y = r*np.sin(theta)
    
    return x,y


t = time.time()

trials = 10**9

r = 1.0

total = 0.0

for i in range(trials):
    x1,y1 = genPt(r)
    x2,y2 = genPt(r)

    total += ((x2-x1)**2 + (y2-y1)**2)**0.5
    
    if not i%10**7:
        print(i)
        print(total/(i+1))
        print('')

    
print(total/trials)
print(time.time()-t)

Chord problem code:

import numpy as np
import time

t = time.time()

trials = 10**10 #Number of random samples

total = 0.0

r = 1.0

x1,y1 = 0.0,r

for i in range(trials):
    
    theta = 2*np.pi*np.random.random()
    
    x2 = r*np.cos(theta)
    y2 = r*np.sin(theta)

    total += ((x2-x1)**2 + (y2-y1)**2)**0.5
    
    if not i%10**7:
        print(i)
        print(total/(i+1))
        print('')

    
print(total/trials)
print(time.time()-t)

I have since made both programs about 2.5x faster by using the polar distance formula to cut down on trig functions (really add up after a few billion trials), by cutting the number of angles generated from 2 to 1 in each trial of the line segment problem, and by switching to using the math module trig functions instead of the NumPy ones. I was somewhat surprised to find that the math module ones are about twice as fast as NumPy's. I've often been impressed by the speed of doing things in NumPy, but it turns out that the Python developers also have some idea what they're doing .

Originally Posted by Ether

3) What's the period of the PRNG you used?

I used numpy.random.random(), which uses the Mersenne Twister sequence and has a period of 2^19937 − 1. Unfortunately I wasn't able to run quite that many trials.

[GeeTwo]

Originally Posted by Greg Woelki

I have since made both programs about 2.5x faster by using the polar distance formula to cut down on trig functions (really add up after a few billion trials), by cutting the number of angles generated from 2 to 1 in each trial of the line segment problem, and by switching to using the math module trig functions instead of the NumPy ones.

You can shave a few more CPU cycles by storing the squares of the radii rather than the radii, so you only need two square roots per segment, rather than three. I am not familiar with this detail in python, but in most languages, the dedicated sqrt(x) function is faster than x ** 0.5.

[Ether]

Greg, try this for the internal points problem: taking advantage of symmetry, generate 2 random points in Quadrant1, then reflect one of the points into Quadrants 2 thru 4 so you get 4 lengths per iteration. Then you can use 4 times fewer iterations. Also, doing this allows you to significantly optimize the code to eliminate a lot of repetitious floating point operations. I think you'll find that approach to be much faster. I was able to get 1 ppm accuracy in 7 minutes of runtime on a 10-year-old Pentium D desktop machine (using 32 bit compiled code). Here's the pseudocode for it:

pio2=pi/2; // initialize constant

// iterate:
r1=random; r2=random; dq=pio2*(random-random); 
tmp1=r1+r2; tmp2=2*sqrt(r1*r2);
cosdq=cos(dq); sindq=sin(dq);
sum +=
  sqrt(tmp1-tmp2*cosdq)	
+ sqrt(tmp1+tmp2*sindq)	
+ sqrt(tmp1+tmp2*cosdq)	
+ sqrt(tmp1-tmp2*sindq);	

// don't forget to divide by 4*iterations when done iterating

Using symmetry for the chords problem, you can reduce the length computation to:

sum += sin(pi*abs(random-random));

... and when the iterations are complete, multiply the sum by 2.