Problem39

30分プログラム、その312。Problem39 - Project Euler

辺の長さが{a,b,c}と整数の3つ組である直角三角形を考え, その周囲の長さをpとする. p = 120のときには3つの解が存在する:
{20,48,52}, {24,45,51}, {30,40,50}
p < 1000 で解の数が最大になる p を求めよ.

数式をがんがん変形したら、かなり高速に解けた。うん、キレイに解けると気持ちいいね。

解法

まず、制約条件はx^2+y^2-z^2=0x+y+z=pの2つ。
これを連立して解くと2(x-p)(y-p)-p^2=0になる。これをちょっと変形すると(x-p)(y-p)=p^2/2になる。

あとは、x,y,pが全て整数なので、(x-p)と(y-p)も整数になる。なので解の個数はp^2/2の約数の個数から求められる。xとyは交換可能なことに注意すると、約数の個数の半分が解の数になる。

で、約数の個数の半分が解の数になるので、約数の前半部分だけを求めたい。前半ということは求めたい数の平方根までで十分。(x < n/xからx^2 < nが求まる。)

使い方

$ time perl problem39.pl
2 1
4 2
6 3
12 6
24 9
36 10
48 12
60 18
120 27
180 30
240 36
360 45
420 54
720 60
840 81
perl problem39.pl  0.24s user 0.02s system 76% cpu 0.349 total

ソースコード

#! /usr/bin/perl
# -*- mode:perl; coding:utf-8 -*-
#
# problem39.pl -
#
# Copyright(C) 2008 by mzp
# Author: MIZUNO Hiroki / mzpppp at gmail dot com
# http://howdyworld.org
#
# Timestamp: 2008/05/27 21:06:25
#
# This program is free software; you can redistribute it and/or
# modify it under MIT Lincence.
#

use strict;
use warnings;

sub divisor($){
    my ($n) = @_;
    grep { $n % $_ == 0 } (1..sqrt($n));
}

sub answer_count($){
    my ($p) = @_;
    return () unless ($p*$p % 2 == 0);

    divisor($p*$p/2);
}

my $max = 0;
for my $p (1..1000){
    my @xs = answer_count($p);
    my $n = @xs;
    if($n > $max){
	print "$p $n\n";
	$max = $n;
    }
}