Special pythagorean triplet

From Rosetta Code
Special pythagorean triplet is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.
Task
The following problem is taken from Project Euler


Related task



ALGOL 68

Using Euclid's formula, as in the XPL0 sample. <lang algol68>BEGIN # find the product of the of the Pythagorian triplet a, b, c where: # 

     #                            a + b + c = 1000, a2 + b2 = c2, a < b < c #
   INT  sqrt 1000 = ENTIER sqrt( 1000 );
   FOR n TO sqrt 1000 DO # m and must have different parity, i.e. one even, one odd       #
       FOR m FROM n + 1 BY 2 TO sqrt 1000 DO
           # a = m^2 - n^2, b = 2mn, c = m^2 + n^2 ( Euclid's formula ), so               #
           # a + b + c = m^2 - n^2 + 2mn + m^2 + n^2 = 2( m^2 + mn ) = 2m( m + n )        # 
           IF m * ( m + n ) = 500 THEN
               INT m2 = m * m, n2 = n * n;
               INT a = m2 - n2;
               INT b = 2 * m * n;
               INT c = m2 + n2;
               print( ( "a = ", whole( a, 0 ), ", b = ", whole( b, 0 ), ", c = ", whole( c, 0 ), newline ) );
               print( ( "a * b * c = ", whole( a * b * c, 0 ), newline ) )
           FI
       OD
   OD

END</lang>

Output:
a = 375, b = 200, c = 425
a * b * c = 31875000

ALGOL W

Translation of: Wren

...but doesn't stop on the first solution (thus verifying there is only one). <lang algolw>% find the Pythagorian triplet a, b, c where a + b + c = 1000 % for a := 1 until 1000 div 3 do begin 

   integer a2, b;
   a2 := a * a;
   for b := a + 1 until 1000 do begin
       integer c;
       c := 1000 - ( a + b );
       if c <= b then goto endB;
       if a2 + b*b = c*c then begin
           write( i_w := 1, s_w := 0, "a = ", a, ", b = ", b, ", c = ", c );
           write( i_w := 1, s_w := 0, "a + b + c = ", a + b + c );
           write( i_w := 1, s_w := 0, "a * b * c = ", a * b * c )
       end if_a2_plus_b2_e_c2 ;
       b := b + 1
   end for_b ;

endB: end for_a .</lang>

Output:
a = 200, b = 375, c = 425
a + b + c = 1000
a * b * c = 31875000

F#

Here I present a solution based on the ideas on the discussion page. It finds all Pythagorean triplets whose elements sum to a given value. It runs in O[n] time. Normally I would exclude triplets with a common factor but for this demonstration it is better to leave them. <lang fsharp> // Special pythagorean triplet. Nigel Galloway: August 31st., 2021 let fG n g=let i=(n-g)/2L match (n+g)%2L with 0L->if (g*g)%(4L*i)=0L then Some(g,i-(g*g)/(4L*i),i+(g*g)/(4L*i)) else None 

                                            |_->if (g*g-2L*i-1L)%(4L*i+2L)=0L then Some(g,i-(g*g)/(4L*i+2L),i+1L+(g*g)/(4L*i+2L)) else None

let E9 n=let fN=fG n in seq{1L..(n-2L)/3L}|>Seq.choose fN|>Seq.iter(fun(n,g,l)->printfn $"%d{n*n}(%d{n})+%d{g*g}(%d{g})=%d{l*l}(%d{l})") [1L..260L]|>List.iter(fun n->printfn "Sum = %d" n; E9 n) </lang>

Output:
Sum = 1
Sum = 2
Sum = 3
Sum = 4
Sum = 5
Sum = 6
Sum = 7
Sum = 8
Sum = 9
Sum = 10
Sum = 11
Sum = 12
9(3)+16(4)=25(5)
Sum = 13
Sum = 14
Sum = 15
Sum = 16
Sum = 17
Sum = 18
Sum = 19
Sum = 20
Sum = 21
Sum = 22
Sum = 23
Sum = 24
36(6)+64(8)=100(10)
Sum = 25
Sum = 26
Sum = 27
Sum = 28
Sum = 29
Sum = 30
25(5)+144(12)=169(13)
Sum = 31
Sum = 32
Sum = 33
Sum = 34
Sum = 35
Sum = 36
81(9)+144(12)=225(15)
Sum = 37
Sum = 38
Sum = 39
Sum = 40
64(8)+225(15)=289(17)
Sum = 41
Sum = 42
Sum = 43
Sum = 44
Sum = 45
Sum = 46
Sum = 47
Sum = 48
144(12)+256(16)=400(20)
Sum = 49
Sum = 50
Sum = 51
Sum = 52
Sum = 53
Sum = 54
Sum = 55
Sum = 56
49(7)+576(24)=625(25)
Sum = 57
Sum = 58
Sum = 59
Sum = 60
100(10)+576(24)=676(26)
225(15)+400(20)=625(25)
Sum = 61
Sum = 62
Sum = 63
Sum = 64
Sum = 65
Sum = 66
Sum = 67
Sum = 68
Sum = 69
Sum = 70
400(20)+441(21)=841(29)
441(21)+400(20)=841(29)
Sum = 71
Sum = 72
324(18)+576(24)=900(30)
Sum = 73
Sum = 74
Sum = 75
Sum = 76
Sum = 77
Sum = 78
Sum = 79
Sum = 80
256(16)+900(30)=1156(34)
Sum = 81
Sum = 82
Sum = 83
Sum = 84
144(12)+1225(35)=1369(37)
441(21)+784(28)=1225(35)
Sum = 85
Sum = 86
Sum = 87
Sum = 88
Sum = 89
Sum = 90
81(9)+1600(40)=1681(41)
225(15)+1296(36)=1521(39)
Sum = 91
Sum = 92
Sum = 93
Sum = 94
Sum = 95
Sum = 96
576(24)+1024(32)=1600(40)
Sum = 97
Sum = 98
Sum = 99
Sum = 100
Sum = 101
Sum = 102
Sum = 103
Sum = 104
Sum = 105
Sum = 106
Sum = 107
Sum = 108
729(27)+1296(36)=2025(45)
Sum = 109
Sum = 110
Sum = 111
Sum = 112
196(14)+2304(48)=2500(50)
Sum = 113
Sum = 114
Sum = 115
Sum = 116
Sum = 117
Sum = 118
Sum = 119
Sum = 120
400(20)+2304(48)=2704(52)
576(24)+2025(45)=2601(51)
900(30)+1600(40)=2500(50)
Sum = 121
Sum = 122
Sum = 123
Sum = 124
Sum = 125
Sum = 126
784(28)+2025(45)=2809(53)
Sum = 127
Sum = 128
Sum = 129
Sum = 130
Sum = 131
Sum = 132
121(11)+3600(60)=3721(61)
1089(33)+1936(44)=3025(55)
Sum = 133
Sum = 134
Sum = 135
Sum = 136
Sum = 137
Sum = 138
Sum = 139
Sum = 140
1600(40)+1764(42)=3364(58)
1764(42)+1600(40)=3364(58)
Sum = 141
Sum = 142
Sum = 143
Sum = 144
256(16)+3969(63)=4225(65)
1296(36)+2304(48)=3600(60)
Sum = 145
Sum = 146
Sum = 147
Sum = 148
Sum = 149
Sum = 150
625(25)+3600(60)=4225(65)
Sum = 151
Sum = 152
Sum = 153
Sum = 154
1089(33)+3136(56)=4225(65)
Sum = 155
Sum = 156
1521(39)+2704(52)=4225(65)
Sum = 157
Sum = 158
Sum = 159
Sum = 160
1024(32)+3600(60)=4624(68)
Sum = 161
Sum = 162
Sum = 163
Sum = 164
Sum = 165
Sum = 166
Sum = 167
Sum = 168
441(21)+5184(72)=5625(75)
576(24)+4900(70)=5476(74)
1764(42)+3136(56)=4900(70)
Sum = 169
Sum = 170
Sum = 171
Sum = 172
Sum = 173
Sum = 174
Sum = 175
Sum = 176
2304(48)+3025(55)=5329(73)
3025(55)+2304(48)=5329(73)
Sum = 177
Sum = 178
Sum = 179
Sum = 180
324(18)+6400(80)=6724(82)
900(30)+5184(72)=6084(78)
2025(45)+3600(60)=5625(75)
Sum = 181
Sum = 182
169(13)+7056(84)=7225(85)
Sum = 183
Sum = 184
Sum = 185
Sum = 186
Sum = 187
Sum = 188
Sum = 189
Sum = 190
Sum = 191
Sum = 192
2304(48)+4096(64)=6400(80)
Sum = 193
Sum = 194
Sum = 195
Sum = 196
Sum = 197
Sum = 198
1296(36)+5929(77)=7225(85)
Sum = 199
Sum = 200
1600(40)+5625(75)=7225(85)
Sum = 201
Sum = 202
Sum = 203
Sum = 204
2601(51)+4624(68)=7225(85)
Sum = 205
Sum = 206
Sum = 207
Sum = 208
1521(39)+6400(80)=7921(89)
Sum = 209
Sum = 210
1225(35)+7056(84)=8281(91)
3600(60)+3969(63)=7569(87)
3969(63)+3600(60)=7569(87)
Sum = 211
Sum = 212
Sum = 213
Sum = 214
Sum = 215
Sum = 216
2916(54)+5184(72)=8100(90)
Sum = 217
Sum = 218
Sum = 219
Sum = 220
400(20)+9801(99)=10201(101)
Sum = 221
Sum = 222
Sum = 223
Sum = 224
784(28)+9216(96)=10000(100)
Sum = 225
Sum = 226
Sum = 227
Sum = 228
3249(57)+5776(76)=9025(95)
Sum = 229
Sum = 230
Sum = 231
Sum = 232
Sum = 233
Sum = 234
4225(65)+5184(72)=9409(97)
5184(72)+4225(65)=9409(97)
Sum = 235
Sum = 236
Sum = 237
Sum = 238
Sum = 239
Sum = 240
225(15)+12544(112)=12769(113)
1600(40)+9216(96)=10816(104)
2304(48)+8100(90)=10404(102)
3600(60)+6400(80)=10000(100)
Sum = 241
Sum = 242
Sum = 243
Sum = 244
Sum = 245
Sum = 246
Sum = 247
Sum = 248
Sum = 249
Sum = 250
Sum = 251
Sum = 252
1296(36)+11025(105)=12321(111)
3136(56)+8100(90)=11236(106)
3969(63)+7056(84)=11025(105)
Sum = 253
Sum = 254
Sum = 255
Sum = 256
Sum = 257
Sum = 258
Sum = 259
Sum = 260
3600(60)+8281(91)=11881(109)

I present results with timing for increasing powers of 10 to demonstrate its O[n] running. 

E9 1000L
40000(200)+140625(375)=180625(425)
Real: 00:00:00.001

E9 10000L
4000000(2000)+14062500(3750)=18062500(4250)
Real: 00:00:00.000

E9 100000L
400000000(20000)+1406250000(37500)=1806250000(42500)
478515625(21875)+1296000000(36000)=1774515625(42125)
Real: 00:00:00.001

E9 1000000L
40000000000(200000)+140625000000(375000)=180625000000(425000)
47851562500(218750)+129600000000(360000)=177451562500(421250)
Real: 00:00:00.005

E9 10000000L
54931640625(234375)+23814400000000(4880000)=23869331640625(4885625)
4000000000000(2000000)+14062500000000(3750000)=18062500000000(4250000)
4785156250000(2187500)+12960000000000(3600000)=17745156250000(4212500)
Real: 00:00:00.040

E9 100000000L
5493164062500(2343750)+2381440000000000(48800000)=2386933164062500(48856250)
400000000000000(20000000)+1406250000000000(37500000)=1806250000000000(42500000)
478515625000000(21875000)+1296000000000000(36000000)=1774515625000000(42125000)
Real: 00:00:00.382

E9 1000000000L
549316406250000(23437500)+238144000000000000(488000000)=238693316406250000(488562500)
40000000000000000(200000000)+140625000000000000(375000000)=180625000000000000(425000000)
47851562500000000(218750000)+129600000000000000(360000000)=177451562500000000(421250000)
Real: 00:00:03.704

Go

Translation of: Wren

<lang go>package main

import ( 

   "fmt"
   "time"

)

func main() { 

   start := time.Now()
   for a := 3; ; a++ {
       for b := a + 1; ; b++ {
           c := 1000 - a - b
           if c <= b {
               break
           }
           if a*a+b*b == c*c {
               fmt.Printf("a = %d, b = %d, c = %d\n", a, b, c)
               fmt.Println("a + b + c =", a+b+c)
               fmt.Println("a * b * c =", a*b*c)
               fmt.Println("\nTook", time.Since(start))
               return
           }
       }
   }

}</lang>

Output:
a = 200, b = 375, c = 425
a + b + c = 1000
a * b * c = 31875000

Took 77.664µs

jq

Works with: jq

Works with gojq, the Go implementation of jq <lang jq>range(1;1000) as $a | range($a+1;1000) as $b | (1000 - $a - $b) as $c | select($a*$a + $b*$b == $c*$c) | {$a, $b, $c, product: ($a*$b*$c)}</lang>

Output:
{"a":200,"b":375,"c":425,"product":31875000}

Or, with a tiny bit of thought: <lang jq>range(1;1000/3) as $a | range($a+1;1000/2) as $b | (1000 - $a - $b) as $c | select($a*$a + $b*$b == $c*$c) | {$a, $b, $c, product: ($a*$b*$c)}}</lang>

Julia

julia> [(a, b, c) for a in 1:1000, b in 1:1000, c in 1:1000 if a < b < c && a + b + c == 1000 && a^2 + b^2 == c^2]
1-element Vector{Tuple{Int64, Int64, Int64}}:
 (200, 375, 425)

or, with attention to timing: 

julia> @time for a in 1:1000
           for b in a+1:1000
               c = 1000 - a - b; a^2 + b^2 == c^2 && @show a, b, c
           end
       end
(a, b, c) = (200, 375, 425)
  0.001073 seconds (20 allocations: 752 bytes)

Nim

My solution from Project Euler:

<lang Nim>import strformat from math import floor, sqrt

var 

 p, s, c : int
 r: float

for i in countdown(499, 1): 

 s = 1000 - i
 p = 1000 * (500 - i)
 let delta = float(s * s - 4 * p)
 r = sqrt(delta)
 if floor(r) == r:
   c = i
   break

echo fmt"Product: {p * c}" echo fmt"a: {(s - int(r)) div 2}" echo fmt"b: {(s + int(r)) div 2}" echo fmt"c: {c}"</lang>

Output:
Product: 31875000
a: 200
b: 375
c: 425

Phix

brute force (83000 iterations)

Not that this is in any way slow (0.1s, or 0s with the displays removed), or deliberately avoids using sensible loop limits, you understand. 

with javascript_semantics
constant n = 1000
integer count = 0
for a=1 to floor(n/3) do
    for b=a+1 to floor((n-a)/2) do
        count += 1
        integer c = n-(a+b)
        if a*a+b*b=c*c then
            printf(1,"a=%d, b=%d, c=%d, a*b*c=%d\n",{a,b,c,a*b*c})
        end if
    end for
end for
printf(1,"%d iterations\n",count)
Output:
a=200, b=375, c=425, a*b*c=31875000
83000 iterations

smarter (166 iterations)

It would of course be 100 iterations if we quit once found. 

with javascript_semantics
constant n = 1000
integer count = 0
for a=2 to floor(n/3) by 2 do
    count += 1
    integer nn2a = n*(n/2-a),
            na = n-a
    if remainder(nn2a,na)=0 then
        integer b = nn2a/na,
                c = n-(a+b)
        printf(1,"a=%d, b=%d, c=%d, a*b*c=%d\n",{a,b,c,a*b*c})
    end if
end for
printf(1,"%d iterations\n",count)
Output:
a=200, b=375, c=425, a*b*c=31875000
166 iterations

PL/M

Based on the XPL0 solution.
As the original 8080 PL/M compiler only has unsigned 8 and 16 bit integer arithmetic, the PL/M long multiplication routines and also a square root routine based that in the PL/M sample for the Frobenius Numbers task are used - which makes this somewhat longer than it would otherwose be... <lang pli>100H: /* FIND THE PYTHAGOREAN TRIPLET A, B, C WHERE A + B + C = 1000 */ 

  /* CP/M BDOS SYSTEM CALL                                                 */
  BDOS: PROCEDURE( FN, ARG ); DECLARE FN BYTE, ARG ADDRESS; GOTO 5; END;
  /* I/O ROUTINES                                                          */
  PRINT$CHAR:   PROCEDURE( C ); DECLARE C BYTE;    CALL BDOS( 2, C ); END;
  PRINT$STRING: PROCEDURE( S ); DECLARE S ADDRESS; CALL BDOS( 9, S ); END;
  PRINT$NL:     PROCEDURE;   CALL PRINT$STRING( .( 0DH, 0AH, '$' ) ); END;

  /* LONG MULTIPLICATION                                                   */
  /* LARGE INTEGERS ARE REPRESENTED BY ARRAYS OF BYTES WHOSE VALUES ARE    */
  /* A SINGLE DECIMAL DIGIT OF THE NUMBER                                  */
  /* THE LEAST SIGNIFICANT DIGIT OF THE LARGE INTEGER IS IN ELEMENT 1      */
  /* ELEMENT 0 CONTAINS THE NUMBER OF DIGITS THE NUMBER HAS                */
  DECLARE LONG$INTEGER  LITERALLY '(21)BYTE';
  DECLARE DIGIT$BASE    LITERALLY '10';
  /* PRINTS A LONG INTEGER                                                 */
  PRINT$LONG$INTEGER: PROCEDURE( N$PTR );
     DECLARE N$PTR ADDRESS;
     DECLARE N BASED N$PTR LONG$INTEGER;
     DECLARE ( D, F ) BYTE;
     F = N( 0 );
     DO D = 1 TO N( 0 );
        CALL PRINT$CHAR( N( F ) + '0' );
        F = F - 1;
     END;
  END PRINT$LONG$INTEGER;
  /* SETS A LONG$INTEGER TO A 16-BIT VALUE                                 */
  SET$LONG$INTEGER: PROCEDURE( LN, N );
     DECLARE ( LN, N ) ADDRESS;
     DECLARE V ADDRESS;
     DECLARE LN$PTR ADDRESS, LN$BYTE  BASED LN$PTR BYTE;
     DECLARE LN$0   ADDRESS, LN$0BYTE BASED LN$0   BYTE;
     LN$0, LN$PTR = LN;
     LN$0BYTE     = 1;
     LN$PTR       = LN$PTR + 1;
     LN$BYTE      = ( V := N ) MOD DIGIT$BASE;
     DO WHILE( ( V := V / DIGIT$BASE ) > 0 );
        LN$PTR   = LN$PTR   + 1;
        LN$BYTE  = V MOD DIGIT$BASE;
        LN$0BYTE = LN$0BYTE + 1;
     END;
  END SET$LONG$INTEGER;
  /* IMPLEMENTS LONG MULTIPLICATION, C IS SET TO A * B                     */
  /*     C CAN BE THE SAME LONG$INTEGER AS A OR B                          */
  LONG$MULTIPLY: PROCEDURE( A$PTR, B$PTR, C$PTR );
     DECLARE ( A$PTR, B$PTR, C$PTR ) ADDRESS;
     DECLARE ( A BASED A$PTR, B BASED B$PTR, C BASED C$PTR ) LONG$INTEGER;
     DECLARE MRESULT LONG$INTEGER;
     DECLARE RPOS    BYTE;

     /* MULTIPLIES THE LONG INTEGER IN B BY THE INTEGER A, THE RESULT      */
     /*     IS ADDED TO C, STARTING FROM DIGIT START                       */
     /*     OVERFLOW IS IGNORED                                            */
     MULTIPLY$ELEMENT: PROCEDURE( A, B$PTR, C$PTR, START );
        DECLARE ( B$PTR, C$PTR )                 ADDRESS;
        DECLARE ( A, START )                     BYTE;
        DECLARE ( B BASED B$PTR, C BASED C$PTR ) LONG$INTEGER;
        DECLARE ( CDIGIT, D$CARRY, BPOS, CPOS )  BYTE;
        D$CARRY = 0;
        CPOS    = START;
        DO BPOS = 1 TO B( 0 );
           CDIGIT = C( CPOS ) + ( A * B( BPOS ) ) + D$CARRY;
           IF CDIGIT < DIGIT$BASE THEN D$CARRY = 0;
           ELSE DO;
              /* HAVE DIGITS TO CARRY                                      */
              D$CARRY = CDIGIT  /  DIGIT$BASE;
              CDIGIT  = CDIGIT MOD DIGIT$BASE;
           END;
           C( CPOS ) = CDIGIT;
           CPOS = CPOS + 1;
        END;
        C( CPOS ) = D$CARRY;
        /* REMOVE LEADING ZEROS BUT IF THE NUMBER IS 0, KEEP THE FINAL 0   */
        DO WHILE( CPOS > 1 AND C( CPOS ) = 0 );
           CPOS = CPOS - 1;
        END;
        C( 0 ) = CPOS;
     END MULTIPLY$ELEMENT ;
     /* THE RESULT WILL BE COMPUTED IN MRESULT, ALLOWING A OR B TO BE C    */
     DO RPOS = 1 TO LAST( MRESULT ); MRESULT( RPOS ) = 0; END;
     /* MULTIPLY BY EACH DIGIT AND ADD TO THE RESULT                       */
     DO RPOS = 1 TO A( 0 );
        IF A( RPOS ) <> 0 THEN DO;
           CALL MULTIPLY$ELEMENT( A( RPOS ), B$PTR, .MRESULT, RPOS );
        END;
     END;
     /* RETURN THE RESULT IN C                                             */
     DO RPOS = 0 TO MRESULT( 0 ); C( RPOS ) = MRESULT( RPOS ); END;
  END;

  /* INTEGER SUARE ROOT: BASED ON THE ONE IN THE PL/M FOR FROBENIUS NUMBERS */
  SQRT: PROCEDURE( N )ADDRESS;
     DECLARE ( N, X0, X1 ) ADDRESS;
     IF N <= 3 THEN DO;
         IF N = 0 THEN X0 = 0; ELSE X0 = 1;
         END;
     ELSE DO;
        X0 = SHR( N, 1 );
        DO WHILE( ( X1 := SHR( X0 + ( N / X0 ), 1 ) ) < X0 );
           X0 = X1;
        END;
     END;
     RETURN X0;
  END SQRT;

  /* FIND THE PYTHAGORIAN TRIPLET                                          */
  DECLARE ( A, B, C, M, N, M2, N2, SQRT$1000 ) ADDRESS;
  DECLARE ( LA, LB, LC, ABC ) LONG$INTEGER;
  SQRT$1000 = SQRT( 1000 );
  DO N = 1 TO SQRT$1000;             /* M AND N MUST HAD DIFFERENT PARITY, */
     DO M = N + 1 TO SQRT$1000 BY 2;             /* I.E. ONE ODD, ONE EVEN */
        /* NOTE: A = M2 - N2, B = 2MN, C = M2 + N2                         */
        /* A + B + C = M2 - N2 + 2MN + M2 + N2 = 2( M2 + MN ) = 2M( M + N )*/ 
        IF ( M * ( M + N ) ) = 500 THEN DO;
           M2 = M * M;
           N2 = N * N;
           CALL SET$LONG$INTEGER( .A, M2 - N2 );
           CALL SET$LONG$INTEGER( .B, 2 * M * N );
           CALL SET$LONG$INTEGER( .C, M2 + N2 );
           CALL LONG$MULTIPLY( .A,   .B, .ABC );
           CALL LONG$MULTIPLY( .ABC, .C, .ABC );
           CALL PRINT$LONG$INTEGER( .ABC );
           CALL PRINT$NL;
        END;
     END;
  END;

EOF</lang>

Output:
31875000

Python

Python 3.8.8 (default, Apr 13 2021, 15:08:03)
Type "help", "copyright", "credits" or "license" for more information.
>>> [(a, b, c) for a in range(1, 1000) for b in range(a, 1000) for c in range(1000) if a + b + c == 1000 and a*a + b*b == c*c]
[(200, 375, 425)]

Raku

<lang perl6>hyper for 1..998 -> $a { 

   my $a2 = $a²;
   for $a + 1 .. 999 -> $b {
       my $c = 1000 - $a - $b;
       last if $c < $b;
       say "$a² + $b² = $c²\n$a  + $b  + $c = {$a+$b+$c}\n$a  × $b  × $c = {$a×$b×$c}"
           and exit if $a2 + $b² == $c²
   }

}</lang>

Output:
200² + 375² = 425²
200  + 375  + 425 = 1000
200  × 375  × 425 = 31875000

REXX

Some optimizations were done such as pre-computing the squares of all possible A's, B's, and C's.

Also, there were multiple shortcuts to limit an otherwise exhaustive search;   Once a sum or a square was too big,
the next integer was used. <lang rexx>/*REXX pgm computes integers A, B, C that solve: 0<A<B<C; A+B+C = 1000; A^2+B^2 = C^2 */ parse arg s hi n . /*obtain optional argument from the CL.*/ if s== | s=="," then s= 1000 /*Not specified? Then use the default.*/ if hi== | hi=="," then hi= 1000 /* " " " " " " */ if n== | n=="," then n= 1 /* " " " " " " */ hi2= hi-2 /*N: number of solutions to find/show.*/ 

                    do j=1  for hi;   @.j= j*j  /*pre─compute squares ──► HI, inclusive*/
                    end  /*j*/
  1. = 0; pad= left(, 9) /*#: the number of solutions found. */
      do       a=2 by 2 for hi2%2;    aa= @.a   /*go hunting for solutions to equations*/
         do    b=a+1  for hi2-a;      ab= a + b /*calculate sum of two  (A,B)  squares.*/
         if ab>hi            then iterate a     /*Sum of A+B>HI?    Then stop with B's */
         aabb= aa + @.b                         /*compute the sum of:   A^2 + B^2      */
            do c=b+1  for hi2-b                 /*test integers that satisfy equations.*/
            if @.c   > aabb  then iterate b     /*Square> A^2+B^2?  Then stop with C's.*/
            if @.c \== aabb  then iterate       /*   "  \=A^2+B^2?  Then keep searching*/
            abc= ab + c                         /*compute the sum of  A  +  B  +  C    */
            if abc == s      then call show     /*Does A+B+C = S?   Then solution found*/
            if abc >  s      then iterate b     /*Is     "   > S?   Then stop with C's.*/
            end   /*c*/
         end      /*b*/
      end         /*a*/

done: say pad pad pad # ' solutions found.' exit 0 /*stick a fork in it, we're all done. */ /*──────────────────────────────────────────────────────────────────────────────────────*/ show: #= #+1; say pad 'a=' a pad "b=" b pad 'c=' c; if #>=n then signal done; return</lang>

output   when using the default inputs:
          a= 200           b= 375           c= 425
                              1  solutions found.

Ring

<lang ring> load "stdlib.ring" see "working..." + nl

time1 = clock() for a = 1 to 1000 

    for b = a+1 to 1000   
        for c = b+1 to 1000
            if (pow(a,2)+pow(b,2)=pow(c,2))
                if a+b+c = 1000
                   see "a = " + a + " b = " + b + " c = " + c + nl
                   exit 3
                ok
            ok
         next
    next

next

time2 = clock() time3 = time2/1000 - time1/1000 see "Elapsed time = " + time3 + " s" + nl see "done..." + nl </lang>

Output:
working...
a = 200 b = 375 c = 425
Elapsed time = 497.61 s
done...

Wren

Very simple approach, only takes 0.013 seconds even in Wren. <lang ecmascript>var a = 3 while (true) { 

   var b = a + 1
   while (true) {
       var c = 1000 - a - b
       if (c <= b) break
       if (a*a + b*b == c*c) {
           System.print("a = %(a), b = %(b), c = %(c)")
           System.print("a + b + c = %(a + b + c)")
           System.print("a * b * c = %(a * b * c)")
           return
       }
       b = b + 1
   }
   a = a + 1

}</lang>

Output:
a = 200, b = 375, c = 425
a + b + c = 1000
a * b * c = 31875000


Incidentally, even though we are told there is only one solution, it is almost as quick to verify this by observing that, since a < b < c, the maximum value of a must be such that 3a + 2 = 1000 or max(a) = 332. The following version ran in 0.015 seconds and, of course, produced the same output: <lang ecmascript>for (a in 3..332) { 

   var b = a + 1
   while (true) {
       var c = 1000 - a - b
       if (c <= b) break
       if (a*a + b*b == c*c) {
           System.print("a = %(a), b = %(b), c = %(c)")
           System.print("a + b + c = %(a + b + c)")
           System.print("a * b * c = %(a * b * c)")
       }
       b = b + 1
   }

}</lang>

XPL0

<lang XPL0>int N, M, A, B, C; for N:= 1 to sqrt(1000) do 

   for M:= N+1 to sqrt(1000) do
       [A:= M*M - N*N; \Euclid's formula
        B:= 2*M*N;
        C:= M*M + N*N;
        if A+B+C = 1000 then
           IntOut(0, A*B*C);
       ]</lang>
Output:
31875000
Retrieved from "https://rosettacode.org/wiki/Special_pythagorean_triplet?oldid=315522"