Chief Delphi

Chief Delphi (http://www.chiefdelphi.com/forums/index.php)
-   Programming (http://www.chiefdelphi.com/forums/forumdisplay.php?f=51)
-   -   numerical computation contest (http://www.chiefdelphi.com/forums/showthread.php?t=110034)

Ether 14-12-2012 14:55

numerical computation contest
 
1 Attachment(s)

A diversion, for anyone so inclined.



Jon Stratis 14-12-2012 15:24

Re: numerical computation contest
 
Does the answer have to be completely numerical, or can it contain known variables/equations (for example, could I say pi/8 versus 0.392699..., or 2*sin (1) were either of those the actual answer. Disclaimer: I have done no math at this point, and I highly doubt that I randomly picked the answer :P)

Ether 14-12-2012 15:25

Re: numerical computation contest
 
Quote:

Originally Posted by Jon Stratis (Post 1201879)
Does the answer have to be completely numerical, or can it contain known variables/equations (for example, could I say pi/8 versus 0.392699..., or 2*sin (1) were either of those the actual answer. Disclaimer: I have done no math at this point, and I highly doubt that I randomly picked the answer :P)

Looking for a numerical answer, in decimal form. Winner is person who gets answer correct to the most decimal digits.

edit:
If you can find an exact explicit closed-form analytical solution using only add, subtract, multiply, divide, powers, roots, exponentiation, logarithms, trig functions, and inverse trig functions, you will be declared the winner.


k4mc 14-12-2012 15:45

Re: numerical computation contest
 
I believe the answer is 44.4984550191007992545541 feet

(although I used a "wolfram alpha method":) , not an analytic one to solve the equation)

Christopher149 14-12-2012 15:53

Re: numerical computation contest
 
Is this close?

h = (in ft)
44.49845501910079925455416001673359894996896435222 95612116880595546581686999560519753347236106529250 14422691969000110577445163659978002373921158891281 96107912233956035279094688973264141285038352833979 92808923081071616667853526325563707812597289771290 34979416624799668990733212350734392869697893501318 89111943559177977405431271912121416410363234025409 52753032492871047810149624788135594125012352228840 06058879524106347465534347539833791041297948729285 69084156937833895515653679687331260207084213606326 0

AGPapa 14-12-2012 16:24

Re: numerical computation contest
 
Quote:

Originally Posted by Christopher149 (Post 1201887)
Is this close?

h = (in ft)
44.49845501910079925455416001673359894996896435222 95612116880595546581686999560519753347236106529250 14422691969000110577445163659978002373921158891281 96107912233956035279094688973264141285038352833979 92808923081071616667853526325563707812597289771290 34979416624799668990733212350734392869697893501318 89111943559177977405431271912121416410363234025409 52753032492871047810149624788135594125012352228840 06058879524106347465534347539833791041297948729285 69084156937833895515653679687331260207084213606326 0

The arc is only one foot longer than the chord. How could the answer possibly be larger than one?

EDIT: Anyway, here is my work. SPOILERS!
/http://i.imgur.com/2IYR8.png

So my guess is zero. It's not exact, but close.

Ether 14-12-2012 16:32

Re: numerical computation contest
 
Quote:

Originally Posted by Christopher149 (Post 1201887)
Is this close?

h = (in ft)
44.49845501910079925455416001673359894996896435222 95612116880595546581686999560519753347236106529250 14422691969000110577445163659978002373921158891281 96107912233956035279094688973264141285038352833979 92808923081071616667853526325563707812597289771290 34979416624799668990733212350734392869697893501318 89111943559177977405431271912121416410363234025409 52753032492871047810149624788135594125012352228840 06058879524106347465534347539833791041297948729285 69084156937833895515653679687331260207084213606326 0

Close enough. Tell the folks how you did it.



Jon Stratis 14-12-2012 16:37

Re: numerical computation contest
 
Quote:

Originally Posted by AGPapa (Post 1201895)
The arc is only one foot longer than the chord. How could the answer possibly be larger than one?

Create a simple approximation. Imagine a point C that lies at the top of the line H. You can then create a triangle ABC, where length AB is known (5280 ft) and length BC = AC = 5281/2 = 2640.5.

We know that h bisects line AB - lets call the intersection between AB and h to be H. we know that the length AH = BH = 5280/2 = 2640.

So, now we know two sides to the right triangle AHC - AH and AC. Taking the square of the hypotenuse minus the square of one side gives us the square of the other side. In other words, 2640.5^2 - 2640^2 = h^2 (the pythagorean theorem).

So, in this extremely rough approximation, we get h = 51.383. Given that, it's not hard to imagine that Christopher's answer could be correct, to some number of decimal places.

AGPapa 14-12-2012 16:48

Re: numerical computation contest
 
Quote:

Originally Posted by Jon Stratis (Post 1201900)
Create a simple approximation. Imagine a point C that lies at the top of the line H. You can then create a triangle ABC, where length AB is known (5280 ft) and length BC = AC = 5281/2 = 2640.5.

We know that h bisects line AB - lets call the intersection between AB and h to be H. we know that the length AH = BH = 5280/2 = 2640.

So, now we know two sides to the right triangle AHC - AH and AC. Taking the square of the hypotenuse minus the square of one side gives us the square of the other side. In other words, 2640.5^2 - 2640^2 = h^2 (the pythagorean theorem).

So, in this extremely rough approximation, we get h = 51.383. Given that, it's not hard to imagine that Christopher's answer could be correct, to some number of decimal places.

Thanks!

Ether 14-12-2012 16:51

Re: numerical computation contest
 
Quote:

Originally Posted by Jon Stratis (Post 1201900)
it's not hard to imagine that Christopher's answer could be correct, to some number of decimal places.

498 decimal places :)

I checked them. They appear to be correct.



Christopher149 14-12-2012 16:56

Re: numerical computation contest
 
Quote:

Originally Posted by Ether (Post 1201897)
Close enough. Tell the folks how you did it.



Since this is a circular arc, I start with a circle. Based off of a diagram like this, I draw a triangle and take theta to be as is shown in that picture. From the right triangle, we know radius * sin(theta) = 5280 / 2. From the definition of arc length, we also know that radius * 2 * theta = 5281. We have two equations and two unknowns, so I substituted one into the other and let Mathematica determine theta. From theta, we easily get the radius. Next, we use this equation relating chord length and the height to the radius. Again, I let Mathematica handle the equation to a ridiculous number of decimal places.

Mathematica Code
Code:

Clear[t, h, r, c, a]
c = 5280
a = 5281
t = t /.
  FindRoot[a/(2 t) * Sin[t] == c/2, {t, 1.8},
  AccuracyGoal -> 100000, PrecisionGoal -> 100000,
  WorkingPrecision -> 500]
r = a/(2 t)
FindRoot[r == c^2/(8 h) + h/2, {h, 100}, AccuracyGoal -> 100000,
 PrecisionGoal -> 100000, WorkingPrecision -> 500]

EDIT: I can get more decimal places ... :)

BigJ 14-12-2012 16:57

Re: numerical computation contest
 
I was bored and tried to find a closed form by hand, it didn't work out by the time I refreshed it and it was over :p


Ether 14-12-2012 17:16

Re: numerical computation contest
 
Quote:

Originally Posted by Christopher149 (Post 1201908)
Mathematica Code
Code:

Clear[t, h, r, c, a]
c = 5280
a = 5281
t = t /.
  FindRoot[a/(2 t) * Sin[t] == c/2, {t, 1.8},
  AccuracyGoal -> 100000, PrecisionGoal -> 100000,
  WorkingPrecision -> 500]
r = a/(2 t)
FindRoot[r == c^2/(8 h) + h/2, {h, 100}, AccuracyGoal -> 100000,
 PrecisionGoal -> 100000, WorkingPrecision -> 500]


Here's Maxima code:
Code:

fpprec: 600$
L:5280$ eps:1$
y: x/sin(x)-(L+eps)/L$
a: bf_find_root(y,x,1/1000,4/100);
h: (L/2)*(1-cos(a))/sin(a);

Quote:

EDIT: I can get more decimal places ... :)
So can I :)
Code:

Maxima 5.27.0 http://maxima.sourceforge.net
using Lisp GNU Common Lisp (GCL) GCL 2.6.8 (a.k.a. GCL)
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.

(%i1) fpprec: 600$
L:5280$ eps:1$
y: x/sin(x)-(L+eps)/L$
a: bf_find_root(y,x,1/1000,4/100);
h: (L/2)*(1-cos(a))/sin(a);

4.4498455019100799254554160016733598949968964352229561211688059554658168\
699956051975334723610652925014422691969000110577445163659978002373921158891281\
961079122339560352790946889732641412850383528339799280892308107161666785352632\
556370781259728977129034979416624799668990733212350734392869697893501318891119\
435591779774054312719121214164103632340254095275303249287104781014962478813559\
412501235222884006058879524106347465534347539833791041297948729285690841569378\
338955156536796873312602070842136063260258182278201500825791387457906343098887\
9258590281352783951539842610751922916884732140489864979154825b1
(%i7)


Christopher149 14-12-2012 17:28

Re: numerical computation contest
 
Quote:

Originally Posted by Ether (Post 1201911)
Code:

4.4498455019100799254554160016733598949968964352229561211688059554658168\
699956051975334723610652925014422691969000110577445163659978002373921158891281\
961079122339560352790946889732641412850383528339799280892308107161666785352632\
556370781259728977129034979416624799668990733212350734392869697893501318891119\
435591779774054312719121214164103632340254095275303249287104781014962478813559\
412501235222884006058879524106347465534347539833791041297948729285690841569378\
338955156536796873312602070842136063260258182278201500825791387457906343098887\
9258590281352783951539842610751922916884732140489864979154825b1
(%i7)


4.4498 ... I thought it was 44.498 ... ;)

Ether 14-12-2012 17:53

Re: numerical computation contest
 
Quote:

Originally Posted by Christopher149 (Post 1201913)
4.4498 ... I thought it was 44.498 ... ;)

Um... The exponent is at the end of those 600 digits :)



Richard Wallace 14-12-2012 20:34

Re: numerical computation contest
 
Russ,

I understand you intended this one as a programming exercise -- hence you posted it in the programming forum.

However, it can also serve as an example of engineering approximation. Jon showed earlier how the Pythagorean theorem gives 51.48 ft as a rough approximation.

For a slightly better one, I used the Taylor series expansion for the sinc function, defined as sinc(x) = sin(x)/x. The first two terms of that expansion are:

sinc(x) = 1 - (x^2)/6 + ....

The arc length is 2rx and the chord length is 2rsin(x), where 2x is the angle subtended by the arc. Since the chord is one mile long and the arc is one foot longer, we have:

sinc(x) = 5280/5281

and the truncated series approximation above gives

x = sqrt(6/5281) = 0.03371 radians, or about 1.931 degrees.

Then the arc is part of a circle with perimeter

(5281/5280) x (180/1.931) = 93.23 miles,

so the radius r is 14.84 miles [this is a very flat arc!]

and the height h is 14.84 x (1 - cos(x)) miles = 44.49 ft.

Back in the day, engineers without benefit of cheap computers could have gotten this close using slide rules; i.e., four figures.

Ether 14-12-2012 21:21

Re: numerical computation contest
 
3 Attachment(s)

Excellent post Richard.

You anticipated where I was going to take this.



DampRobot 14-12-2012 23:23

Re: numerical computation contest
 
Ether,

I'm not a programmer, but I love a good geometrical challenge.

The equations I came up with were: (r-h)sin(l/2r)=x/2 and 2rsin(l/2r)=(x^2/4+h^2)^1/2
where r was the radius of the curve, h was the height specified in your diagram, l was the length of the curve, and x was the length of the chord. To find the numerical solution, input values for l and x, and solve for h or r.

Was the "challenge" to write a program that found the solutions of these functions through guess and check? Or did I come up with a different solution? Plugging in the values on this thread with my trusty graphing calculator, I got a reasonable approximation of the answers posted on this thread.

Tom Line 14-12-2012 23:28

Re: numerical computation contest
 
Quote:

Originally Posted by DampRobot (Post 1202018)
Ether,

I'm not a programmer, but I love a good geometrical challenge.

The equations I came up with were: (r-h)sin(l/2r)=x/2 and 2rsin(l/2r)=(x^2/4+h^2)^1/2
where r was the radius of the curve, h was the height specified in your diagram, l was the length of the curve, and x was the length of the chord. To find the numerical solution, input values for l and x, and solve for h or r.

Was the "challenge" to write a program that found the solutions of these functions through guess and check? Or did I come up with a different solution? Plugging in the values on this thread with my trusty graphing calculator, I got a reasonable approximation of the answers posted on this thread.

Even getting a reasonable approimation is good enough, because there is no exact answer. I simplified until I ended up with a formula that looked something like this:

sin (theta/2) = (theta * Arc) / (2 * chord length)

With theta being the angle inside the two radii that would have gone from the center of the circle to the end of the cords.

That can only be solved iteratively, graphically, or by a lookup table (for those of us who remember not having calculators....)

Once you have theta, the height is trivial.

Ether 15-12-2012 01:15

Re: numerical computation contest
 
3 Attachment(s)
Quote:

Originally Posted by DampRobot (Post 1202018)
Was the "challenge" to write a program that found the solutions of these functions through guess and check?

There is no explicit solution so you have to use iteration.

For folks like Christopher149 who know how to use modern math tools like Mathematica this is very simple to do, and you can get answers of arbitrary precision.

Mathematica is very expensive, but there are free CAS (Computer Algebra System) tools like Maxima which anybody can download and install for free. Maxima can also solve this problem with arbitrary precision (see my earlier post where I used Maxima to get 600 digits).

Attached are:

- The derivation of the equations

- A simple Delphi console app code which uses a binary search to find theta


If you replace sin(theta)/theta - 5280/5281 with a 12th degree Taylor expansion (see Richard's earlier post about Taylor expansion), you can squeeze out 16 decimal digits of accuracy with this simple code using Intel's X86 Extended Precision floats



IKE 15-12-2012 12:50

Re: numerical computation contest
 
This reminds me of a chapter in one of my favorite books:
The History of PI

Specifically, you guys are reminding me of chapter 10, "The Digit Hunters".

I am very impressed. Happy hunting.

Ether 15-12-2012 13:51

Re: numerical computation contest
 
Quote:

Originally Posted by DampRobot (Post 1202018)
The equations I came up with were: (r-h)sin(l/2r)=x/2 and 2rsin(l/2r)=(x^2/4+h^2)^1/2
where r was the radius of the curve, h was the height specified in your diagram, l was the length of the curve, and x was the length of the chord.

I think you made a mistake somewhere along the way. Those equations don't look correct.

x/2 = rsin(l/2r), not (r-h)sin(l/2r)

and

2rsin(l/r) = x, not (x^2/4+h^2)^1/2



Ether 15-12-2012 17:30

Re: numerical computation contest
 
Quote:

Originally Posted by k4mc (Post 1201885)
I believe the answer is 44.4984550191007992545541 feet

(although I used a "wolfram alpha method":) , not an analytic one to solve the equation)

Sorry, I completely overlooked your post !

You were the first to get a correct answer.

Would you please post you work?



Alpha Beta 15-12-2012 18:20

Re: numerical computation contest
 
Quote:

Originally Posted by Ether (Post 1202099)
Sorry, I completely overlooked your post !

You were the first to get a correct answer.

I don't think you overlooked it. There seems to be a delay in people getting their first couple of posts published when they are new to the forum. The awkward thing is those posts then sneak back into the list of posts in the spot they were written at. Easy to miss them completely if a lot of posts have been written in the mean time.

Ether 15-12-2012 18:51

Re: numerical computation contest
 
Quote:

Originally Posted by Alpha Beta (Post 1202103)
...those posts then sneak back into the list of posts in the spot they were written at.

Aha. So I'm not losing my mind (yet).



Ether 15-12-2012 21:41

Re: numerical computation contest
 
Quote:

Originally Posted by Tom Line (Post 1202019)
sin (theta/2) = (theta * Arc) / (2 * chord length)

That should be sin (theta/2) = (theta * chord) / (2 * arc)



Ether 16-12-2012 18:21

Re: numerical computation contest
 
2 Attachment(s)

Just for fun, here's a graph showing the behavior of the 2nd, 4th, 6th, and 8th degree Taylor approximations for the sinc function from zero to 3pi/2 radians.

Here's Maxima code to generate the graph.

Maxima is free software with an active and helpful mailing list.



Richard Wallace 16-12-2012 18:38

Re: numerical computation contest
 
Quote:

Originally Posted by Ether (Post 1202342)

Just for fun, here's a graph ...

It would be interesting to see a graph comparing the 2nd, 4th, 6th, and 8th order expansions of 1 - sinc(x) plotted vs. x, for 0 < x < 0.05. The range of most interest will be near 1/5280 .... ;)

I'd do it myself, but have been unable to locate Maxima as an iPad app (yet).

Ether 16-12-2012 20:19

Re: numerical computation contest
 
3 Attachment(s)
Quote:

Originally Posted by Richard Wallace (Post 1202350)
It would be interesting to see a graph comparing the 2nd, 4th, 6th, and 8th order expansions of 1 - sinc(x) plotted vs. x, for 0 < x < 0.05. The range of most interest will be near 1/5280 .... ;)

... or 1/5281 :)

sin(θ)/θ = 5280/5281 = 1 - 1/5281 => 1 - sin(θ)/θ = 1/5281

The Taylor expansions (even the 2nd degree one) are so accurate over the range 0<x<0.05 that you can't see the difference on the graph at that scale. Which maybe was your point :)

Also attached is the 64 decimal digit precision computation of θ using the zero of the function y = 5280/5281 - sin(θ)/θ and its Taylor expansions.

Finally, notice that the Taylor8 expansion actually gives a more accurate answer than the original function itself, when using only 16 digit precision. Presumably, this is because Taylor8 is numerically more stable.



Bill_B 16-12-2012 21:11

Re: numerical computation contest
 
This whole thread reminds me of an easy (for you all) geometry question. Imagine a sphere the diameter of the Earth. Add to that a belt around the equator that has no space between it and the surface of the sphere. By how much must the length of that belt be increased to allow a playing card to be slipped underneath the belt at all points along it?

Ether 16-12-2012 21:19

Re: numerical computation contest
 
Quote:

Originally Posted by Bill_B (Post 1202425)
This whole thread reminds me of an easy (for you all) geometry question. Imagine a sphere the diameter of the Earth. Add to that a belt around the equator that has no space between it and the surface of the sphere. By how much must the length of that belt be increased to allow a playing card to be slipped underneath the belt at all points along it?

That problem has a simple explicit closed-form solution.

C = pi*d

C' = pi*d' = pi*(d+2t)

C'-C = pi*(d+2t) - pi*d = pi*2t

... where "t" is the thickness of the card.



Richard Wallace 16-12-2012 21:23

Re: numerical computation contest
 
Actually this problem is considerably more practical:

Imagine a steel rail one mile long, pinned to earth at points A and B one mile apart. The rail is warmed by about thirty Fahrenheit degrees, so that it elongates by one foot, causing it to deform from its original straightness into a circular arc. How far does the arc deviate from the original line, at a point halfway between A and B?

Bill_B 16-12-2012 22:19

Re: numerical computation contest
 
Quote:

Originally Posted by Ether (Post 1202427)
That problem has a simple explicit closed-form solution.

C = pi*d

C' = pi*d' = pi*(d+2t)

C'-C = pi*(d+2t) - pi*d = pi*2t

... where "t" is the thickness of the card.



I told you it was easy! What usually hangs fledgling geometers is the "2t" to allow the card under on both sides of the sphere. Everyone else just guesses "big" numbers because the sphere is so big.


All times are GMT -5. The time now is 21:39.

Powered by vBulletin® Version 3.6.4
Copyright ©2000 - 2017, Jelsoft Enterprises Ltd.
Copyright © Chief Delphi