1 /*
2 * Find and display the minimum sum of the squares for integers.
3 *
4 * Copyright (C) 2007 George Gesslein II.
5
6 This library is free software; you can redistribute it and/or
7 modify it under the terms of the GNU Lesser General Public
8 License as published by the Free Software Foundation; either
9 version 2.1 of the License, or (at your option) any later version.
10
11 This library is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 Lesser General Public License for more details.
15
16 You should have received a copy of the GNU Lesser General Public
17 License along with this library; if not, write to the Free Software
18 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
19
20 The chief copyright holder can be contacted at gesslein@mathomatic.org, or
21 George Gesslein II, P.O. Box 224, Lansing, NY 14882-0224 USA.
22
23 */
24
25 /*
26 * Usage: matho-sumsq [numbers]
27 *
28 * This program finds the minimum number of positive integers squared
29 * and summed to equal a given positive integer. If the number of squared
30 * integers is 2 (default "multi"), all combinations with 2 squares are displayed,
31 * otherwise only the first solution found is displayed.
32 *
33 * If nothing is specified on the command line, the program gets its numbers from standard input.
34 *
35 * This file is mostly useful for testing the long integer square root function lsqrt() in file "lsqrt.c".
36 */
37
38 #include <stdio.h>
39 #include <stdlib.h>
40 #include <math.h>
41 #include <errno.h>
42 #include <assert.h>
43
44 #define true 1
45 #define false 0
46
47 #define squared(a) ((a) * (a))
48
49 int multi = 2; /* display all combinations for this number of squares or less */
50
51 long squares[4];
52
53 long lsqrt(long y);
54
55 /*
56 * Find the sum of "n" squares that equal "d1" and display if found.
57 *
58 * Return false if not found.
59 */
60 int
sumsq(long d1,int n)61 sumsq(long d1, int n)
62 {
63 int i = 0;
64 long d2;
65 long save = 0;
66 int found_one = false;
67
68 d2 = d1;
69 try_again:
70 for (;;) {
71 if (i == 2) {
72 save = d2;
73 }
74 squares[i] = lsqrt(d2);
75 assert(squares[i] >= 0);
76 d2 -= squared(squares[i]);
77 i++;
78 if (i >= n) {
79 break;
80 }
81 }
82 if (d2 == 0) {
83 for (i = 0; i < n; i++) {
84 d2 += squared(squares[i]);
85 }
86 if (d2 != d1) {
87 fprintf(stderr, "Error: Result doesn't compare identical to original number!\n");
88 exit(EXIT_FAILURE);
89 }
90 for (i = 0; i < (n - 1); i++) {
91 if (squares[i] < squares[i+1]) {
92 goto skip_print;
93 }
94 }
95 found_one = true;
96 printf("%ld = %ld^2", d1, squares[0]);
97 for (i = 1; i < n; i++) {
98 if (squares[i] != 0)
99 printf(" + %ld^2", squares[i]);
100 }
101 printf("\n");
102 skip_print:
103 assert(found_one);
104 if (n < 2 || n > multi) {
105 return found_one;
106 }
107 }
108 switch (n) {
109 case 4:
110 if (squares[2] > squares[n-1]) {
111 squares[2] -= 1;
112 d2 = save - squared(squares[2]);
113 i = 3;
114 goto try_again;
115 }
116 case 3:
117 if (squares[1] > squares[n-1]) {
118 squares[1] -= 1;
119 d2 = d1 - squared(squares[0]) - squared(squares[1]);
120 i = 2;
121 goto try_again;
122 }
123 case 2:
124 if (squares[0] > squares[n-1]) {
125 squares[0] -= 1;
126 d2 = d1 - squared(squares[0]);
127 i = 1;
128 goto try_again;
129 }
130 }
131 return found_one;
132 }
133
134 /*
135 * Display the shortest sums of squares for long integer "d1".
136 *
137 * Return the minimum number of summed squares required to represent "d1".
138 */
139 int
findsq(long d1)140 findsq(long d1)
141 {
142 if (sumsq(d1, 1))
143 return 1;
144 if (sumsq(d1, 2))
145 return 2;
146 if (sumsq(d1, 3))
147 return 3;
148 if (sumsq(d1, 4))
149 return 4;
150 fprintf(stderr, "Whoops! Can't find the sum of four squares that equal %ld.\n", d1);
151 exit(EXIT_FAILURE);
152 }
153
154 int
main(int argc,char * argv[])155 main(int argc, char *argv[])
156 {
157 int i;
158 long d1 = 0;
159 char *cp, buf[1000];
160
161 if (argc > 1) {
162 for (i = 1; i < argc; i++) {
163 errno = 0;
164 d1 = strtol(argv[i], &cp, 10);
165 if (errno) {
166 perror(argv[i]);
167 exit(EXIT_FAILURE);
168 }
169 if (d1 < 0) {
170 fprintf(stderr, "Invalid command-line argument: \"%s\", positive integer required.\n", argv[i]);
171 exit(EXIT_FAILURE);
172 }
173 if (*cp == '+') {
174 for (;; d1 += 1) {
175 findsq(d1);
176 }
177 exit(EXIT_SUCCESS);
178 }
179 if (argv[i][0] == '\0' || *cp) {
180 fprintf(stderr, "Invalid number: \"%s\".\n", argv[i]);
181 exit(EXIT_FAILURE);
182 }
183 findsq(d1);
184 }
185 } else {
186 for (fflush(NULL); fgets(buf, sizeof(buf), stdin); fflush(NULL)) {
187 errno = 0;
188 d1 = strtol(buf, &cp, 10);
189 if (errno) {
190 perror(NULL);
191 continue;
192 }
193 if (cp != buf && d1 == 0) {
194 exit(EXIT_SUCCESS);
195 }
196 if (d1 <= 0) {
197 fprintf(stderr, "Positive integer required; 0 to quit.\n");
198 continue;
199 }
200 findsq(d1);
201 }
202 }
203 exit(EXIT_SUCCESS);
204 }
205