Text 150, 232 rader
Skriven 2004-08-27 09:17:06 av DAVID WILLIAMS (1:250/514)
Ärende: Pythagorean Triples
===========================
I changed my priorities. Now, I want a program that will go to the 
highest numbers possible, regardless of speed. So I changed all the 
integer operations ("\", MOD, Booleans) to floating-point ones, and 
made all the numbers double-precision floating-point ones. The result 
is this version. Previously, the program bombed at about A=46000, 
because long-integers overflowed. This one will handle values of A of 
more than 100 million! 
  
Of course, it would take almost forever to print out all the triples up 
to there. This version is *much* slower than the old one. But you can 
start it wherever you want. 
  
                                 dow 
  
---------------------------------------------------- 
  
' PYTRIP-L.BAS - Pythagorean Triples 
  
' David O. Williams. 2004 
' david.williams@ablelink.org 
  
' Calculates and prints integer triples, A, B, C, such that 
' A < B < C,  they have no common factors (>1), and A^2 + B^2 = C^2. 
' The list is printed in order of increasing A, then of B. 
' A counter of triples is also shown. 
  
' A more detailed explanation is at the end of the main module. 
  
DECLARE FUNCTION NoComFacs# (X#, Y#) 
DECLARE SUB PrintOut (A#, B#, C#) 
  
DEFDBL A-Z 
  
CLS 
  
DO 
  PRINT "Both values must be >= 0, and Maximum >= Minimum" 
  PRINT 
  INPUT "Minimum value of smallest number in triple"; Mn 
  INPUT "Maximum value of smallest number in triple"; Mx 
  PRINT 
  Mx = INT(Mx) 
  Mn = INT(Mn) 
  IF Mn >= 0 AND Mx >= 0 AND Mx >= Mn THEN EXIT DO 
  BEEP 
  PRINT "Illegal value/s!" 
  PRINT 
LOOP 
  
T = Mx > Mn OR (Mn <> 4 AND Mn - 4 * INT(Mn / 4) <> 2) 
IF Mx > 2 AND T THEN 
  PRINT "Count"; TAB(21); "A"; TAB(35); "B"; TAB(58); "C" 
  PRINT 
ELSE 
  PRINT "There are no triples in this range!" 
  END 
END IF 
  
Z = SQR(2) + 1 
  
H = SQR(Mn / Z) 
T = 2 * INT(H / 2) 
IF H - T > 1 THEN Q = T + 1 ELSE Q = T - 1 
  ' Q is initial loop limit for odd cases 
J = Q + 2 
K = INT(J * J * Z) 
  
R = INT(SQR(Mn / (Z + Z))) ' initial loop limit for even cases 
L = R + 1 
M = INT(L * L * Z) 
  
' start of main loop 
FOR A = Mn TO Mx 
  SELECT CASE A - 4 * INT(A / 4) 
    CASE 0 
      S = A / 2 
      IF S > M THEN 
        R = L 
        L = L + 1 
        M = INT(L * L * Z) 
      END IF 
      FOR F = R TO 1 STEP -1 
        G = S / F 
        IF G = INT(G) THEN 
          IF F + G - 2 * (INT(F / 2) + INT(G / 2)) = 1 THEN 
            IF NoComFacs(F, G) THEN 
              X = G * G 
              Y = F * F 
              B = X - Y 
              C = X + Y 
              PrintOut A, B, C 
            END IF 
          END IF 
        END IF 
      NEXT 
    CASE 1, 3 
      IF A > K THEN 
        Q = J 
        J = J + 2 
        K = INT(J * J * Z) 
      END IF 
      FOR D = Q TO 1 STEP -2 
        E = A / D 
        IF E = INT(E) THEN 
          IF NoComFacs(D, E) THEN 
            T = D * D 
            B = (E * E - T) / 2 
            C = B + T 
            PrintOut A, B, C 
          END IF 
        END IF 
      NEXT 
  END SELECT 
NEXT 
  
END 
  
'---------------------------------------------------------- 
  
' Brief explanation: 
  
' Pythagorean triples, if they are in lowest terms with no common 
' factors, can be written as: Odd#^2 + Even#^2 = Big#^2. Big# is the 
' largest integer (corresponding to the hypotenuse of the right-angled 
' triangle). However, Odd# may be smaller or larger than Even#. 
  
' If the three above numbers have no common factors, it is possible 
' to define two odd positive integers, D and E, with E > D, such that: 
  
' Odd# = D * E 
' Even# = (E^2 - D^2) / 2 
' Big# = (E^2 + D^2) / 2 
  
' These definitions satisfy Pythagoras's Theorem. It is possible to 
' prove that all valid triples, in lowest terms, can be written this 
' way, with odd-integer values of D and E. (They must both be odd, so 
' that their product is Odd#.) 
  
' If a triple is written as A, B, C, with A < B < C, C must correspond 
' to Big#, but A can be either Odd# or Even#, and B will be the other. 
  
' The program treats these two possibilities separately. If A is odd, 
' so it must be Odd#, the program simply searches for two odd integers 
' whose product is A. After confirming that they have no common 
' factors, which would mean that the triple is not in lowest terms, 
' the program calculates B and C from them, and prints them out. 
  
' If A is even, it is useful to define two further numbers, F and G, 
' with G > F, such that: 
  
' D = G - F 
' E = F + G 
  
' This means that: 
  
' F = (E - D) / 2 
' G = (D + E) / 2 
  
' Since D and E are both odd, their sum and difference are both even, 
' so F and G are integers. However, the sum and difference of F and G 
' are E and D, which are odd, which implies that the parities of F 
' and G must be opposites, so one is odd and the other even. 
  
' Since, in this case, A is the same as Even#, it is given by: 
  
' A = (E^2 - D^2) / 2 
  
' Writing G - F for D and G + F for E, and simplifying, this gives: 
  
' A / 2 = F * G 
  
' Since F and G are of opposite parity, this means that A / 2 must be 
' even, so A must be a multiple of 4. There are no valid triples when 
' A MOD 4 = 2. 
  
' When A is a multiple of 4, the program looks for factors of A / 2. 
' It checks that they are of opposite parity, and have no common 
' factors. It then uses these values of F and G to calculate B and 
' C, using the easily-proved formulae: 
  
' B = G^2 - F^2 
' C = G^2 + F^2 
  
' (Strictly, the parity check could be omitted. Since they are factors 
' of an even number, at least one of the found values of F and G must 
' be even. The possibility that they are both even would be detected 
' by the common-factor test. However, the parity check is much simpler 
' and faster than the common-factor routine, so having it in the 
' program saves some time.) 
  
' The FOR ... NEXT loops in the odd and even CASEs search for the 
' factors D and E, or F and G, respectively, The loops run from 
' high to low values of the counting variables, since this makes the 
' triples appear in the desired order. Also, the ranges of the loops 
' are limited so that only triples in which B > A appear. The 
' variables Q and R govern these ranges. They are initialized near 
' the start of the program, and are incremented in the main loop as 
' the value of A increases. The numbers K and M are pre-calculated 
' limits. When A (or S, which is just A / 2 when A is even) passes 
' the relevant limit, Q or R is incremented appropriately, and a new 
' value of K or M is calculated. In the vast majority of iterations 
' around the loops, it is not necessary to increment Q or R, so all 
' that has to be done is a simple comparison, e.g IF A > K THEN, 
' which turns out not to be true. This saves a lot of time. Also, 
' this method does not use the slow operation SQR inside any loops. 
' A full mathematical treatment of the situation uses the number 
' (SQR(2) + 1) several times. This number is therefore treated as a 
' constant in the program, named Z. 
  
'                       = end = 
  
  
FUNCTION NoComFacs (X, Y) ' non-zero if X and Y have no common factors 
  U = Y 
  V = X 
  DO WHILE V > 1 
    W = U - V * INT(U / V) 
    U = V 
    V = W 
  LOOP 
  NoComFacs = V 
END FUNCTION 
  
SUB PrintOut (A, B, C) 
  STATIC N 
  N = N + 1 
  PRINT N; TAB(20); A; TAB(34); B; TAB(57); C 
END SUB 
  
---------------------------------------------------- 
--- Platinum Xpress/Win/WINServer v3.0pr5
 * Origin: The Bayman BBS,Toronto, (416)698-6573 - 1:250/514 (1:250/514)