You mean printf, not fprint. Here is code that prints using printf a minimal-length representation of each integer satisfying your requirements. Note how cat and $ are used in the printf to construct a format string with a variable number of format codes (the %ds).

restart:
Squares:= [$1..12]^~2;
for n from 129 to 129+13^2 do
    found:= false;
    for m to nops(Squares) do
        for c in combinat:-choose(Squares,m) do
            if add(c)=n then
                printf(cat("\n%d = %d^2", " + %d^2"$(m-1)), n, sqrt~(c)[]);
                found:= true
            fi
        until found
    until found;
    if not found then printf("\n%d: no solution", n) fi
od:

129 = 2^2 + 5^2 + 10^2
130 = 3^2 + 11^2
131 = 1^2 + 3^2 + 11^2
132 = 1^2 + 3^2 + 4^2 + 5^2 + 9^2
133 = 4^2 + 6^2 + 9^2
134 = 2^2 + 3^2 + 11^2
135 = 1^2 + 2^2 + 3^2 + 11^2
136 = 6^2 + 10^2
137 = 4^2 + 11^2
138 = 1^2 + 4^2 + 11^2
139 = 3^2 + 7^2 + 9^2
140 = 2^2 + 6^2 + 10^2
141 = 2^2 + 4^2 + 11^2
142 = 5^2 + 6^2 + 9^2
143 = 1^2 + 5^2 + 6^2 + 9^2
144 = 12^2
145 = 1^2 + 12^2
146 = 5^2 + 11^2
147 = 1^2 + 5^2 + 11^2
148 = 2^2 + 12^2
149 = 7^2 + 10^2
150 = 1^2 + 7^2 + 10^2
151 = 1^2 + 2^2 + 5^2 + 11^2
152 = 4^2 + 6^2 + 10^2
153 = 3^2 + 12^2
154 = 1^2 + 3^2 + 12^2
155 = 3^2 + 5^2 + 11^2
156 = 1^2 + 3^2 + 5^2 + 11^2
157 = 6^2 + 11^2
158 = 1^2 + 6^2 + 11^2
159 = 1^2 + 3^2 + 7^2 + 10^2
160 = 4^2 + 12^2
161 = 1^2 + 4^2 + 12^2
162 = 4^2 + 5^2 + 11^2
163 = 1^2 + 4^2 + 5^2 + 11^2
164 = 8^2 + 10^2
165 = 1^2 + 8^2 + 10^2
166 = 3^2 + 6^2 + 11^2
167 = 1^2 + 3^2 + 6^2 + 11^2
168 = 2^2 + 8^2 + 10^2
169 = 5^2 + 12^2
170 = 7^2 + 11^2
171 = 1^2 + 7^2 + 11^2
172 = 1^2 + 3^2 + 4^2 + 5^2 + 11^2
173 = 2^2 + 5^2 + 12^2
174 = 2^2 + 7^2 + 11^2
175 = 1^2 + 2^2 + 7^2 + 11^2
176 = 1^2 + 3^2 + 6^2 + 7^2 + 9^2
177 = 2^2 + 3^2 + 8^2 + 10^2
178 = 3^2 + 5^2 + 12^2
179 = 3^2 + 7^2 + 11^2
180 = 6^2 + 12^2
181 = 9^2 + 10^2
182 = 1^2 + 9^2 + 10^2
183 = 1^2 + 5^2 + 6^2 + 11^2
184 = 2^2 + 6^2 + 12^2
185 = 8^2 + 11^2
186 = 1^2 + 8^2 + 11^2
187 = 1^2 + 4^2 + 7^2 + 11^2
188 = 1^2 + 2^2 + 3^2 + 5^2 + 7^2 + 10^2
189 = 2^2 + 8^2 + 11^2
190 = 3^2 + 9^2 + 10^2
191 = 1^2 + 3^2 + 9^2 + 10^2
192 = 1^2 + 3^2 + 5^2 + 6^2 + 11^2
193 = 7^2 + 12^2
194 = 1^2 + 7^2 + 12^2
195 = 5^2 + 7^2 + 11^2
196 = 4^2 + 6^2 + 12^2
197 = 2^2 + 7^2 + 12^2
198 = 1^2 + 2^2 + 7^2 + 12^2
199 = 2^2 + 5^2 + 7^2 + 11^2
200 = 6^2 + 8^2 + 10^2
201 = 4^2 + 8^2 + 11^2
202 = 9^2 + 11^2
203 = 1^2 + 9^2 + 11^2
204 = 2^2 + 6^2 + 8^2 + 10^2
205 = 5^2 + 6^2 + 12^2
206 = 2^2 + 9^2 + 11^2
207 = 1^2 + 2^2 + 9^2 + 11^2
208 = 8^2 + 12^2
209 = 1^2 + 8^2 + 12^2
210 = 5^2 + 8^2 + 11^2
211 = 3^2 + 9^2 + 11^2
212 = 2^2 + 8^2 + 12^2
213 = 7^2 + 8^2 + 10^2
214 = 1^2 + 7^2 + 8^2 + 10^2
215 = 2^2 + 3^2 + 9^2 + 11^2
216 = 4^2 + 6^2 + 8^2 + 10^2
217 = 3^2 + 8^2 + 12^2
218 = 4^2 + 9^2 + 11^2
219 = 1^2 + 4^2 + 9^2 + 11^2
220 = 1^2 + 3^2 + 5^2 + 8^2 + 11^2
221 = 10^2 + 11^2
222 = 1^2 + 10^2 + 11^2
223 = 1^2 + 2^2 + 4^2 + 9^2 + 11^2
224 = 4^2 + 8^2 + 12^2
225 = 9^2 + 12^2
226 = 1^2 + 9^2 + 12^2
227 = 5^2 + 9^2 + 11^2
228 = 1^2 + 5^2 + 9^2 + 11^2
229 = 2^2 + 9^2 + 12^2
230 = 3^2 + 10^2 + 11^2
231 = 1^2 + 3^2 + 10^2 + 11^2
232 = 1^2 + 2^2 + 5^2 + 9^2 + 11^2
233 = 5^2 + 8^2 + 12^2
234 = 3^2 + 9^2 + 12^2
235 = 1^2 + 3^2 + 9^2 + 12^2
236 = 3^2 + 5^2 + 9^2 + 11^2
237 = 4^2 + 10^2 + 11^2
238 = 6^2 + 9^2 + 11^2
239 = 1^2 + 6^2 + 9^2 + 11^2
240 = 1^2 + 3^2 + 7^2 + 9^2 + 10^2
241 = 4^2 + 9^2 + 12^2
242 = 1^2 + 4^2 + 9^2 + 12^2
243 = 3^2 + 7^2 + 8^2 + 11^2
244 = 10^2 + 12^2
245 = 1^2 + 10^2 + 12^2
246 = 5^2 + 10^2 + 11^2
247 = 1^2 + 5^2 + 10^2 + 11^2
248 = 2^2 + 10^2 + 12^2
249 = 1^2 + 2^2 + 10^2 + 12^2
250 = 5^2 + 9^2 + 12^2
251 = 7^2 + 9^2 + 11^2
252 = 1^2 + 7^2 + 9^2 + 11^2
253 = 3^2 + 10^2 + 12^2
254 = 1^2 + 3^2 + 10^2 + 12^2
255 = 2^2 + 7^2 + 9^2 + 11^2
256 = 1^2 + 2^2 + 7^2 + 9^2 + 11^2
257 = 6^2 + 10^2 + 11^2
258 = 1^2 + 6^2 + 10^2 + 11^2
259 = 3^2 + 5^2 + 9^2 + 12^2
260 = 4^2 + 10^2 + 12^2
261 = 6^2 + 9^2 + 12^2
262 = 1^2 + 6^2 + 9^2 + 12^2
263 = 5^2 + 6^2 + 9^2 + 11^2
264 = 2^2 + 4^2 + 10^2 + 12^2
265 = 11^2 + 12^2
266 = 1^2 + 11^2 + 12^2
267 = 1^2 + 8^2 + 9^2 + 11^2
268 = 1^2 + 4^2 + 7^2 + 9^2 + 11^2
269 = 2^2 + 11^2 + 12^2
270 = 7^2 + 10^2 + 11^2
271 = 1^2 + 7^2 + 10^2 + 11^2
272 = 3^2 + 5^2 + 6^2 + 9^2 + 11^2
273 = 2^2 + 5^2 + 10^2 + 12^2
274 = 3^2 + 11^2 + 12^2
275 = 1^2 + 3^2 + 11^2 + 12^2
276 = 5^2 + 7^2 + 9^2 + 11^2
277 = 4^2 + 6^2 + 9^2 + 12^2
278 = 2^2 + 3^2 + 11^2 + 12^2
279 = 3^2 + 7^2 + 10^2 + 11^2
280 = 6^2 + 10^2 + 12^2
281 = 4^2 + 11^2 + 12^2
282 = 1^2 + 4^2 + 11^2 + 12^2
283 = 3^2 + 7^2 + 9^2 + 12^2
284 = 2^2 + 6^2 + 10^2 + 12^2
285 = 8^2 + 10^2 + 11^2
286 = 1^2 + 8^2 + 10^2 + 11^2
287 = 6^2 + 7^2 + 9^2 + 11^2
288 = 1^2 + 6^2 + 7^2 + 9^2 + 11^2
289 = 8^2 + 9^2 + 12^2
290 = 5^2 + 11^2 + 12^2
291 = 1^2 + 5^2 + 11^2 + 12^2
292 = 1^2 + 5^2 + 8^2 + 9^2 + 11^2
293 = 7^2 + 10^2 + 12^2
294 = 1^2 + 7^2 + 10^2 + 12^2
295 = 5^2 + 7^2 + 10^2 + 11^2
296 = 4^2 + 6^2 + 10^2 + 12^2
297 = 2^2 + 7^2 + 10^2 + 12^2
298 = 3^2 + 8^2 + 9^2 + 12^2

In your expression sin(4*pi)+(1/3)*cos(6*pi), pi is treated as the independent variable, not the well-known constant. The constant in Maple is Pi, with a capital P. So perhaps you mean to plot sin(4*x)+(1/3)*cos(6*x) or sin(4*Pi*x)+(1/3)*cos(6*Pi*x)? The plotting domain (the range of the independent variable) is specified as the second argument to plot:

plot(sin(4*x)+(1/3)*cos(6*x), x= 0..2*Pi)

If you don't specify a plotting domain, the command picks one for you, usually -10..10.

I concur with the VV and Axel: The limit (over the real numbers) exists. If it makes you more comfortable, translate the problem to an unbounded domain:

g:= x-> sin(sin(x)/x)/(sin(x)/x):
limit(g(x), x= infinity);
                           1
plot(g(x), x= 1..50);

Now clearly for any epsilon > 0, there exist uncountably many x such that |g(x) - 1| < epsilon, and 1 is the only value for which this is true.

Or, using a more traditional definition: Let D be the domain of g. For any epsilon > 0, there exists M such that for all x in (D intersect [M, infinity)), |g(x) - 1| < epsilon. (This formulation automatically implies the uniqueness of 1.)

Is this what you have in mind?

plots:-animate(
    plots:-display, 
    ['''plot3d''']~(
        [[sqrt(8*y),t,y], [y^2,t,y]], 
        y= 0..Y, t= 0..2*Pi, coords= cylindrical, color=~["Khaki", "LightBlue"],
        labels= ['x', 'z', 'y'], axes= normal, scaling= constrained,
        view= [-4..4, -4..4, 0..4]
    ),   
    Y= 0..2, frames= 60, paraminfo= false
);

Like this:

plots:-animate(
   plot3d,
   [[[r,t,r^2/8],[r,t,sqrt(r)]], r= 0..4, t= 0..T, coords= cylindrical],
   T= 0..2*Pi, frames= 60
);
);

You can do

coeff(select(t-> degree(t,x)=1, S), x, 1)

@Mike Graber The default assignment is by value, not by address. If I execute your code, I get 7 and 9, not 9 and 9. Try doing a restart and doing your example again.

restart:
A:= 7:  B:= A:  A:= 9:  C:= A:  B, C;
                             7, 9

To refer to the address of a name rather than its value, enclose it in single forward quotes (aposthropes):

restart:
A:= 7:  B:= 'A':  A:= 9:  C:= A:  B, C
                              9, 9

Most of the rules are described on help page ?eval. There are some subtleties involved with the above when working with local variables. Occasionally, one needs to use eval to get the final value. And even less frequently, one needs to use eval with a numeric second argument to prevent evaluation to the final value.
 

I don't know if this is the source of your problem, but it is nonetheless an error: The use statement doesn't work correctly in procedures. Instead, change

use XMLTools in

to

uses XMLTools;

and remove the end use;

@Scot Gould You wrote:

• Is this a problem with just this function (and probably Maximization) or inherit in the entire Optimization package?

This will be a problem for any of the numerous numerical analysis procedures that except input functions in procedure form. This includes fsolve, dsolve, Statistics:-Fit, and numeric integration itself.

• I believe I interpret it as: “If the list [K,r] is not a list of numerical values then return the procedural call with the arguments included.”

Yes, that's correct.

• In Maple coding, is this the preferred technique for telling a calling function/procedure that the call failed?

No, it's only useful when the returned unevaluated call may be subject to subsequent evaluation. Genuine failures should be signaled with errors. Acer's technique is very common for procedures that are meant to ignore (rather than reject) nonnumeric input.

• This may require me do learn more Maple coding than I have historically thought was necessary.

This syntax is very simple: return, procname and args are keywords that can be used in any procedure. So, returning the unevaluated procedure call is always literally return 'procname'(args).

• In listing local variables, I failed to include the variable x, yet Maple did not flag it like it flagged f1 and f2 when I did not list them as local.

Nonlocals only get flagged when they're on the left side of a direct assignment (:=), the index variable of a for statement, or (in some versions of Maple) the index variable of s​​​​​​eq, add, or mul.

• Does not listing “x” as a local variable cause problems, or is it simply a case of sloppy coding?

It will cause problems if x has been globally assigned.

• And yes, it makes sense the K and r of the function is not the arguments K and r. I was unaware of the use of “:-K” as representing the variable of an expression.

The prefix :- on a variable name refers to the global instance of the variable. It's required when there are global of local instances of variables with the same obstensive name. The quotes, ':-K', mean "ignore the variable's assigned value, if any." Note that those are single forward quotes (aposthropes).

• In procedure proc4 execute nearly 8 times quicker when the “unapply” command is included within the procedure as opposed to placing that function call outside of the procedure, i.e., defining df as a function and not an expression?

This is definitely a "down the rabbit hole" question. What Acer is doing is directly instantiating the numeric values of K and r in the body of the procedure constructed by unapply. This makes it a purely numeric procedure, which can thus be processed much faster by the numeric integrator. To go farther down that rabbit hole, see ?evalhf. I'm not saying that evalhf is necessarily being used in this case. There are other instances (for example, externally compiled code) where a purely numeric procedure will make things faster.

The appropriate built-in command is rtable_scanblock. Although it is very fast, it has complicated syntax. See ?rtable_scanblock.

FirstPositive:= proc(A::{Vector, Array})
local ind:= (); 
    rtable_scanblock(
        A, [],
        proc(v,i) if v::realcons and is(v>0) then ind:= i[]; true fi end proc,
        'stopresult'= true
    );
    ind
end proc
:
FirstPositive(<0, 0, -1, 2, 0, 5>);
                               4

If your vectors are strictly numeric, the above can be improved. If they are also strictly real, it can be further improved.

Here is a much neater encoding of that theorem. In particular, pay attention to the angle-bracket matrix constructor < ...; ...>, which greatly simplifies the construction of block matrices. Also note the use of simultaneous multiple assignment (..., ...):= (..., ...), which eliminates the need for intermediary variables when updating multiple recurrences. (The parentheses are not required for multiple assignment; I just use them for clarity.)

Domineering:= proc(m::posint, n::posint, x::name, y::name)
local G0:= <<1>>, G1:= <<0>>, Z:= <<0>>;
    to n do
        (G1,G0,Z):= (<y*G0, Z; Z, Z>, <G0+G1, x*G0; G0, Z>, <Z, Z; Z, Z>)
    od;
    (G0^m)[1,1]
end proc
:                
Domineering(3,3,x,y); #Example call

The above uses 1D input (aka Maple Input).

In the 4th Reply to this Post, I wrote a tridiagonal solver that's much faster than anything offered by the LinearAlgebra package. As you suggest, my method uses prefactoring followed by back substitution. Thus, if the same coefficient matrix is used with different right-side vectors (as is typical), only the back substitution needs to be repeated.

Like this:

d:= x-> x^2+1: #denominator
Adj:= x-> 2*d(x)-x^3;
R:= 2 - (x^3+CurveFitting:-PolynomialInterpolation(
    [[1, Adj(1)], [2, Adj(2)], [3, Adj(3)]], x
))/d(x);

As is well documented, ListTools:-Classify is much more efficient than ListTools:-Categorize. It's also much more efficient than Acer's GL, and it can be used for your purpose via a two-line procedure:

GraphEigenvalueClassify:= (L::list([{Graph, name}, list(algebraic)]))->
    ([lhs,rhs]~@op@(g-> [map(`?[]`, g, [1])[]])~)(ListTools:-Classify(g-> g[2], L))
:

Here is an efficiency comparison:

Download grouplist_test.mw

You haven't indicated whether your lists of eigenvalues are always sorted. If they aren't, and you want to group by sorted eigenvalues, then simply change g[2] to sort(g[2]) in my procedure. This has neglible effect on the efficiency.

Most of the new syntax features introduced in Maple 2018 and Maple 2019, such as do-until loops, only work in 1D input (aka Maple Input) mode. If you change your input mode, then you'll have no syntax errors. The issue raised by Tom Leslie is completely different. It's a logic error---an infinite loop. They don't give error messages.

By the way, I like your neat coding style and your use of the newer syntax.