Fifty Shades of J/Chapter 35

From J Wiki
Jump to: navigation, search

Table of Contents ... Glossary ... Previous Chapter ... Next Chapter

How many Obtuse angled triangles are there?

Principal Topics

 ?. (roll) o. (circle functions) : (grade down) ~ (reflex), simulation, random numbers, random angles, accept-reject technique

In the world, that is. Lots, I know, but let me ask rather what proportion of all triangles in the world are obtuse-angled? Since any two of the angles, say a and b, of a triangle can be determined independently, the every possible triangle corresponds to a point within the large right-angled triangle in the diagram below

Fsoj35.png

Acute-angled triangles must have a+b>90o and so their points all lie within the inner triangle, from which it follows that 75% of all triangles are obtuse-angled. Well let’s ask J and see what he/she/it says using statistical Monte-Carlo methods. Assume that a point is represented by a 2-list of numbers such as 4 7. The verb `py` calculates the square of the distance between two such points

   py=.(+/@:*:@:-/)"2    NB. square of dist between 2 points
   py 4 7,:1 3           NB. (3,4,5) triangle
25

A triangle is a 3-list of 2-lists (points). To calculate the squares of the three sides, take the points in pairs and apply `py` :

   sides =: 2 [\ 4&$
   ]t=: 3 2$_3 _2 _1 6 0 2
_3 _2
_1  6
 0  2
   <"2 sides t
┌─────┬────┬─────┐
│_3 _2│_1 6│ 0  2│
│_1  6│ 0 2│_3 _2│
└─────┴────┴─────┘

A triangle is obtuse if the square of the largest side exceeds the sum of the squares of the other two, so define ifbigone which uses reflexive grade down to test whether the largest of a set of numbers is greater than the sum of all the others

   ifbigone =: ({. > +/@}.)@\:~     NB. 1 if biggest > +/rest
   ifbigone 49 43 100
1

This leads naturally to

   obtest =: ifbigone@:(py"2@sides)"2 NB. 1 if obtuse
   obtest t
1

randtri1 generates y random triangles with vertices in the unit square:

   randtri1 =: 0 ?@$~ ,&3 2

   <"2 randtri1 3
┌─────────────────┬─────────────────┬─────────────────┐
│0.577101 0.430063│0.356508 0.879054│0.190733 0.629572│
│0.689446 0.610395│0.898293  0.11287│ 0.45149 0.964217│
│0.401262 0.267585│0.184275 0.444405│0.482228 0.170644│
└─────────────────┴─────────────────┴─────────────────┘

Everything is now in place to test a single random triangle with given resolution :

   trial1=.obtest@:randtri1     NB. 1 if obtuse
   trial1 8
0 0 1 1 1 0 0 1

Of the triangles in this particular trial four of eight were obtuse-angled. Here are a thousand trials:

   +/ trial1 1000
740

... and a few more :

   +/ trial1 1000
739
   +/ trial1 1000
708
   +/ trial1 1000
736

By now you should be getting a feel for where the true value lies, somewhere between 72 and 75 per cent. Another way to conduct the simulation is to generate random angles rather than random sides. Start by writing a verb which generates a number of random values (angles) between 0 and 2π to use as angles from the origin to the vertex:

   rnd =: ?@$&0    NB. y rand uniform nos. in (0,1)
   TAU =: 2p1             NB. tauday.com
   randang =: TAU * rnd    NB. random angles in radians
   randang 3
5.66985 3.78605 2.69206

The sines and cosines of each angle generate a point within the unit circle, and so random angles can be transformed to random points within the unit circle by

   cs =: 2 1 o.~/~ randang

   cs 8
 0.515802 _0.856708
  0.19751 _0.980301
 0.590711  0.806883
_0.659526 _0.751682
_0.451618 _0.892211
_0.956942  _0.29028
_0.980289 _0.197569
  _0.4976  0.867407
   

and then to y random triangles in a unit circle by

   randtri2 =: [: cs ,&3

   <"2 randtri2 4
┌──────────────────┬──────────────────┬─────────────────────┬────────────────────┐
│_0.718971  0.69504│ 0.503937 _0.86374│0.0327604   _0.999463│_0.998892 _0.0470605│
│  0.96262 0.270857│_0.795225 0.606314│ _0.91815   _0.396232│_0.554286   0.832326│
│ 0.442591 0.896724│_0.401744 0.915752│ 0.999965 _0.00831741│ 0.727516  _0.686091│
└──────────────────┴──────────────────┴─────────────────────┴────────────────────┘

All is in place for determining whether these triangles are obtuse or not, and for doing some more experiments

   ([: +/ trial2)@:(1000"_) ^:(>: i. 4) ''
763 710 754 736

A little bit higher than our previous estimate, but close to the theoretical 75%. Maybe sampling from a square frame includes disproportionately ‘skinny’ acute-angled triangles which have one vertex tucked into one of the far corners. To test this conjecture, sample the coordinates within a square, but filter only those triangles all of whose points lie inside the unit circle. The polynomial transforms the triangle vertices to have absolute value <: 1 . All of these must have length <: 1.

   ssq =: +/@:*:    NB. sum of squares
   randtri3 =: (#~ ([: (*./"1) 1 >: ssq"1))@:(_1 2 p. randtri1)

   NB. graphics shows the remaining points are within the unit circle
   require 'plot'
   'point'plot _2 j./\ , randtri3 999

As before everything is in place for the next experiment. Since randtri3 generates fewer triangles than requested, include the average :

      (+/ % #) trial3 10000
0.729905

Not much different from `trial1` – ah well, that’s one conjecture gone for a bang!

Another accept-reject technique for random triangles consists of generating the values of three random sides in order, and rejecting the result if the sides fail to form a triangle because one of them is bigger than the sum of the other two :

   randtri4 =: (#~ (>./ < +/ - >./)"1)@:([: rnd ,&3)

   randtri4 9
0.482861 0.631554 0.575322
0.790549  0.74182 0.983299
 0.21859 0.814304 0.746637
0.758181 0.448893 0.836001

Since randtri4 returns lengths of sides the obtuseness test is applied to the squares of the sides :

   trial4 =: ifbigone"1@:*:@randtri4

   (+/ % #) trial4 99999
0.57199

Seems I’ve mislaid nearly 20% of the world’s obtuse-angled triangles (careless of me!) Another possibility based on sides is to choose two sides at random, say these were 6 and 10. The range of acceptable values for a third side is then from 10-6 to 10+6, or more generally from m-n to m+n where m and n are the lengths of the sides in decreasing order. The third side can then be taken as (m-n) plus a number drawn at random from 0 to 2n:

   avg =: +/ % #
   rnd =: ?@$&0
   third =: , (>./ + <./ * _1 2 p. ?@:0:)
   rantri5 =: [: third"1 [: rnd ,&2
   trial5 =: ifbigone"1@:*:@rantri5   
   avg trial5 100000
0.73331

Let’s get back to the original question . How many obtuse-angled triangles did I say there were? ...

Code Summary

trial1=: obtest@randtri1                NB. first random trial
obtest =: ifbigone@:(py"2@sides)"2      NB. 1 if obtuse
ifbigone =: ({. > +/@}.)@\:~            NB. 1 if biggest > +/rest
py =: (+/@:*:@:-/)"2                    NB. sq of dist between 2 pts
sides =: 2 [\ 4&$                       NB. sides from coords
randtri1 =: 0 ?@$~ ,&3 2                NB. generate y triangles in unit square

trial2 =: obtest@randtri2               NB. second random trial
randtri2 =: [: cs ,&3
cs =: 2 1 o.~/~ randang                 NB. cos/sin of rand angle
TAU =: 2p1
randang =: TAU * rnd                    NB. random angles in radians
rnd =: ?@$&0                            NB. random values on (0,1) with shape y

trial3 =: obtest@randtri3               NB. third random trial

randtri3 =: (#~ ([: (*./"1) 1 >: ssq"1))@:(_1 2 p. randtri1)
                                        NB. generate at most y triangles in unit square,
                                        NB. discarding those outside the inscribed circle.

ssq =: +/@:*:                           NB. sum of squares

trial4 =: ifbigone@:*:@randtri4         NB. fourth random trial
randtri4 =: (#~ (>./ < +/ - >./)"1)@:([: rnd ,&3)
                                        NB. generate at most y triangles by random length sides
                                        NB. discarding triples which can not form triangles.

trial5 =: ifbigone@:*:@rantri5          NB. fifth random trial
rantri5 =: [: third"1 [: rnd ,&2        NB. choose 2 sides, 3rd side randomly chosen as plausible
third =: , (>./ + <./ * _1 2 p. ?@:0:)  NB. join feasible third side

Script

File:Fsojc35.ijs