Tri-48

January 30, 2018

There are two possibilities. If 48 is the hypotenuse of the triangle, the problem is to find two integers with squares that sum to a target value, which we examined in a previous exercise; in this case, there are no solutions in triangles:

(define (squares n)
  (let loop ((x (isqrt n)) (y 0) (zs '()))
    (cond ((< x y) zs)
          ((< (+ (* x x) (* y y)) n) (loop x (+ y 1) zs))
          ((< n (+ (* x x) (* y y))) (loop (- x 1) y zs))
          (else (loop (- x 1) (+ y 1) (cons (list x y) zs))))))
> (squares (* 48 48))
((48 0))

The other possibility is that 48 is one of the other legs of the triangle, in which case we count up from 1, testing at each step if the hypotenuse is a square:

(define (tri n)
  (let loop ((a 1))
    (let ((c (+ (* a a) (* n n))))
      (when (square? c)
        (display (sort < (list a (isqrt c) n)))
        (newline)))
    (loop (+ a 1))))
> (tri 48)
(14 48 50)
(20 48 52)
(36 48 60)
(48 55 73)
(48 64 80)
(48 90 102)
(48 140 148)
(48 189 195)
(48 286 290)
(48 575 577)
CTRL-C
break>

This is intentionally an infinite loop. Beiler shows that any number expressible as 16k (such as 48) always has exactly ten solutions, so the list is exhaustive.

You can run the program at https://ideone.com/yX11Fm:

Advertisements

Pages: 1 2

7 Responses to “Tri-48”

  1. Paul said

    In Python. Look here here for solving the second part.

    from ma.primee import divisor_gen # generate all divisors
    
    z = 48
    z2 = z ** 2
    
    # try x^2 + y^2 = z^2
    x, y = z, 0
    while x >= y > 0:
        test = x ** 2 + y ** 2 - z2
        if test == 0:
            print(x, y)
            x -= 1
            y += 1
        elif test < 0:
            y += 1
        else: x -= 1
    
    # try x^2 - y^2 = z^2
    for div in divisor_gen(z2):
        if div & 1:  # as z = 4 * k
            continue
        if div >= z:
            continue
        u, v = div, z2 // div
        x = (u + v) // 2
        y = (v - u) // 2
        print(x, y)
    
  2. A Perl 1-liner…. we know the largest value of a side is (48^2/4 + 1) as {(n+1)^2-(n-1)^2 = 4n = 48^2} so we loop through all values – work out what the other side is… remove entries for which the other side is 0 or non-integer… finally display the results…

    print map  {"@{[sort{$a<=>$b}@{$_}]}\n"}
          grep {$_->[2]&&$_->[2]==int $_->[2]}
          map  {[48,$_,sqrt abs$_*$_-2304]}
          1..577;
    

    This gives:

    14 48 50
    20 48 52
    36 48 60
    48 55 73
    48 64 80
    48 90 102
    48 140 148
    48 189 195
    48 286 290
    48 575 577

  3. Mark said
    [(a, b) for a in range(1,1000) for b in range(1,1000) if a**2 + 48**2 == b**2]
    
  4. Just realised there is a minor optimisation in mine to avoid duplicates… – can change the bottom of the loop to 34….

    if you want to do for arbitrary N…

    $N=shift;
    print
      map { "@{[sort{$a<=>$b}@{$_}]}\n"}
      grep{$_->[2]&&$_->[2]==int$_->[2]}
      map {[$N,$_,sqrt abs$_*$_-$N*$N]}
      ($N/sqrt 2)..($N*$N/2+1);
    
  5. kernelbob said

    This was a fun exercise, thanks! To keep it interesting, i avoided a brute force search. Instead, I relied on Euclid’s formula. https://en.wikipedia.org/wiki/Pythagorean_triple#Generating_a_triple

    There is still optimization to be done, but this iterates the innermost loop 71 times.

    def primes():
        """Generate all prime numbers."""
        yield 2
        yield 3
        yield 5
        raise RuntimeError('How many primes are there, anyway?')
    
    def factors(n):
        if n == 1:
            yield 1
        else:
            for p in primes():
                count = 0
                while n % p == 0:
                    n //= p
                    count += 1
                if count:
                    for f in factors(n):
                        for i in range(count + 1):
                            yield p**i * f
                    break
    
    def pythagorean_triples(num):
        for f in factors(num):
            k = num // f
            if f % 2 == 0:
                for n in factors(f):
                    m = f // 2 // n
                    if m <= n:
                        continue
                    a = m**2 - n**2
                    b = 2 * m * n
                    assert(b == f)
                    c = m**2 + n**2
                    yield tuple((k * a, k * b, k * c))
    
            for d in range(1, f*f):
                md2 = f + d**2
                if md2 % (2 * d) == 0:
                    m = md2 // (2 * d)
                    n = m - d
                    if n <= 0:
                        break
                    a = m**2 - n**2
                    assert(a == f)
                    b = 2 * m * n
                    c = m**2 + n**2
                    yield (k * a, k * b, k * c)
    
    
    for t in sorted(set(tuple(sorted(t)) for t in pythagorean_triples(48))):
        print(t)
        a, b, c = t
        assert(a**2 + b**2 == c**2)
    
  6. matthew said

    Let’s see: any k,a,b, with a and b coprime and of different parity (so a+b is odd) define a Pythagorean triple (k(a²+b²), 2kab, k(a²-b²)). If n = kj with j = a²+b², j must be an odd divisor of n, expressible as the sum of two squares; if n = 2kab, a and b are any coprime divisors of n/2 of different parity; if n = kj with j = a²-b² = (a+b)(a-b), j is product of two odd factors c,d of n, c > d, with a = (c+d)/2, b = (c-d)/2.

    For n = 48, odd factors are just 1 and 3, even factors are 2,4,6,8,12,24,48.

    Neither 1 or 3 are the sum of (non-zero) squares.
    Even factors of 24 coprime to 1 are: 2,4,8,6,12,24 gives a,b,k = (2,1,12),(4,1,6),(8,1,3),(6,1,4),(12,1,2),(24,1,1)
    Even factors of 24 coprime to 3 are: 2,4,8 gives a,b,k = (2,3,4),(4,3,3),(8,3,1)
    Single pair of odd factors 1 and 3 give a,b,k = ((3+1)/2,(3-1)/2,12) = (2,1,16)

    10 triples total.

    Can argue the same way for any number of form 16p with p an odd prime to get 10 different legs. Might also get a hypotenuse if p is the sum of two squares (eg. if n = 80, p = 5, n = 16(1²+2²) – see https://en.wikipedia.org/wiki/Fermat%27s_theorem_on_sums_of_two_squares).

  7. Daniel said

    Here’s a brute force Python 2.7 solution that also works for other sized sides.

    import itertools
    
    # https://en.wikipedia.org/wiki/Integer_square_root#Using_only_integer_division
    def isqrt(n):
        # Works for integers n >= 1
        x = n
        while True:
            y = (x + (n / x)) / 2
            if y - x in (0, 1): return x
            x = y
    
    def triples(x):
        """Returns the other sides of Pythagorean triples that include x."""
        x2 = x ** 2
        for a in itertools.count(1):
            a2 = a ** 2
            b2 = x2 - a2         # treats x as the hypotenuse
            if b2 > a2:          # avoids duplicates:
                b = isqrt(b2)    #   e.g., (3,4) and (4,3) for x == 5
                if b * b == b2:
                    yield (a, b)
            c2 = a2 + x2         # treats x as a leg
            c = isqrt(c2)
            if c * c == c2:
                yield (a, c)
            if a == c: return
    
    for x in triples(48):
        print x
    

    Output:

    (14, 50)
    (20, 52)
    (36, 60)
    (55, 73)
    (64, 80)
    (90, 102)
    (140, 148)
    (189, 195)
    (286, 290)
    (575, 577)
    

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: