46 people like it.

Primitive Pythagorean triples

Primitive Pythagorean triples generator. It uses an Algorithm found on Wolfram MathWorld and the F# PowerPack matrix library.

Primitive Pythagorean triples generator

 1: 
 2: 
 3: 
 4: 
 5: 
 6: 
 7: 
 8: 
 9: 
10: 
11: 
12: 
13: 
14: 
15: 
16: 
17: 
18: 
19: 
20: 
21: 
22: 
23: 
24: 
25: 
26: 
27: 
28: 
29: 
30: 
31: 
32: 
#r "FSharp.PowerPack.dll"
open Microsoft.FSharp.Math

// Pythagorean triples : http://en.wikipedia.org/wiki/Pythagorean_triple
// Algorithm: http://mathworld.wolfram.com/PythagoreanTriple.html Equations(2 .. 9)

// Primitive Pythagorean triples generator
let primitives take =
    let org = rowvec [ 3.; 4.; 5. ]
    
    let U = matrix [[  1.; 2.; 2. ];
                    [ -2.;-1.;-2. ];
                    [  2.; 2.; 3. ];]

    let A = matrix [[  1.; 2.; 2. ];
                    [  2.; 1.; 2. ];
                    [  2.; 2.; 3. ];]

    let D = matrix [[ -1.;-2.;-2. ];
                    [  2.; 1.; 2. ];
                    [  2.; 2.; 3. ];]

    let triplets (p:RowVector<float>) = (p*U,p*A,p*D)
    
    let rec primitives' next cont acc = 
        if take next then
            let u,a,d = triplets next
            next::acc |> primitives' u (primitives' a (primitives' d cont))
        else
            cont acc
    
    primitives' org id []

Example: Project Euler 75

 1: 
 2: 
 3: 
 4: 
 5: 
 6: 
 7: 
 8: 
 9: 
10: 
11: 
12: 
13: 
14: 
15: 
16: 
// http://projecteuler.net/index.php?section=problems&id=75
let limit = 1500000.

let perimeter (t :RowVector<_>) = t.[0] + t.[1] + t.[2]

let perimeterUnder1500000 t = perimeter t <= limit

let solution = primitives perimeterUnder1500000
               |> Seq.map perimeter
               // include the perimeter of the primitive and all the multiples
               |> Seq.map(fun p -> seq { p .. p .. limit })
               |> Seq.concat
               |> Seq.countBy id
               |> Seq.filter(snd >> (=) 1)
               |> Seq.length
solution |> printfn "solution: %d"
namespace Microsoft
namespace Microsoft.FSharp
val primitives : take:('a -> bool) -> 'a list

Full name: Script.primitives
val take : ('a -> bool)
val org : 'a
val U : obj
val A : obj
val D : obj
val triplets : ('a -> 'a * 'a * 'a)
val p : 'a
Multiple items
val float : value:'T -> float (requires member op_Explicit)

Full name: Microsoft.FSharp.Core.Operators.float

--------------------
type float = System.Double

Full name: Microsoft.FSharp.Core.float

--------------------
type float<'Measure> = float

Full name: Microsoft.FSharp.Core.float<_>
val primitives' : ('a -> ('a list -> 'b) -> 'a list -> 'b)
val next : 'a
val cont : ('a list -> 'b)
val acc : 'a list
val u : 'a
val a : 'a
val d : 'a
val id : x:'T -> 'T

Full name: Microsoft.FSharp.Core.Operators.id
val limit : float

Full name: Script.limit
val perimeter : t:'a -> float

Full name: Script.perimeter
val t : 'a
val perimeterUnder1500000 : t:'a -> bool

Full name: Script.perimeterUnder1500000
val solution : int

Full name: Script.solution
module Seq

from Microsoft.FSharp.Collections
val map : mapping:('T -> 'U) -> source:seq<'T> -> seq<'U>

Full name: Microsoft.FSharp.Collections.Seq.map
val p : float
Multiple items
val seq : sequence:seq<'T> -> seq<'T>

Full name: Microsoft.FSharp.Core.Operators.seq

--------------------
type seq<'T> = System.Collections.Generic.IEnumerable<'T>

Full name: Microsoft.FSharp.Collections.seq<_>
val concat : sources:seq<#seq<'T>> -> seq<'T>

Full name: Microsoft.FSharp.Collections.Seq.concat
val countBy : projection:('T -> 'Key) -> source:seq<'T> -> seq<'Key * int> (requires equality)

Full name: Microsoft.FSharp.Collections.Seq.countBy
val filter : predicate:('T -> bool) -> source:seq<'T> -> seq<'T>

Full name: Microsoft.FSharp.Collections.Seq.filter
val snd : tuple:('T1 * 'T2) -> 'T2

Full name: Microsoft.FSharp.Core.Operators.snd
val length : source:seq<'T> -> int

Full name: Microsoft.FSharp.Collections.Seq.length
val printfn : format:Printf.TextWriterFormat<'T> -> 'T

Full name: Microsoft.FSharp.Core.ExtraTopLevelOperators.printfn

More information

Link:http://fssnip.net/39
Posted:13 years ago
Author:Cesar Mendoza
Tags: mathematics , euler problem , powerpack , algorithms