| blob | c434b6f0fbe5de7ad7cb0ea3280cef515b8e6530 |
1 ;; Maxima functions for finding the maximum or minimum
2 ;; Copyright (C) 2005, 2007, 2021 Barton Willis
4 ;; Barton Willis
5 ;; Department of Mathematics
6 ;; University of Nebraska at Kearney
7 ;; Kearney NE 68847
8 ;; willisb@unk.edu
10 ;; This source code is licensed under the terms of the Lisp Lesser
11 ;; GNU Public License (LLGPL). The LLGPL consists of a preamble, published
12 ;; by Franz Inc. (http://opensource.franz.com/preamble.html), and the GNU
13 ;; Library General Public License (LGPL), version 2, or (at your option)
14 ;; any later version. When the preamble conflicts with the LGPL,
15 ;; the preamble takes precedence.
17 ;; This library is distributed in the hope that it will be useful,
18 ;; but WITHOUT ANY WARRANTY; without even the implied warranty of
19 ;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 ;; Library General Public License for details.
22 ;; You should have received a copy of the GNU Library General Public
23 ;; License along with this library; if not, write to the
24 ;; Free Software Foundation, Inc., 51 Franklin St, Fifth Floor,
25 ;; Boston, MA 02110-1301, USA.
32 ;; Let's remove built-in symbols from list for user-defined properties.
47 (merror (intl:gettext "The value of maxmin_effort must be a nonnegative integer; found ~M ~%") y))))
51 ;; Return true if there is pi and pj in the CL list p such that
52 ;; x is between pi and pj. This means that either pi <= x <= pj or
53 ;; pj <= x <= pi. For example, 2x is between x and 3x.
55 ;; Strangely, sign((a-b)*(b-a)) --> pnz but sign(expand((a-b)*(b-a))) --> nz.
56 ;; To workaround this weirdness, we could call expand instead of factor on
57 ;; (mul (sub x pk) (sub qk x))), but let's not do that.
59 ;; The betweenp simplification is done last; this has some interesting effects:
60 ;; max(x^2,x^4,x^6,x^2+1) (standard simplification) --> max(x^4,x^6,x^2+1)
61 ;; (betweenp) --> max(x^4,x^6,x^2+1). If the betweenp simplification were done
62 ;; first, we'd have max(x^2,x^4,x^6,x^2+1) --> max(x^2,x^6,x^2+1) --> max(x^6,x^2+1).
64 ;; The function $factor refuses to factor an expression that has an exponent that
65 ;; exceeds factor_max_degree. I think this is a better heuristic than the conssize
66 ;; method used by factor-if-small (defined in compar.lisp). So let's use $factor, not
67 ;; factor-if-small. We could locally set the value of factor_max_degree, but let's not.
69 ;; Removing factor from (csign ($factor (mul (sub x pk) (sub qk x))) causes max
70 ;; to miss the simplification max(x^2,x^4,x^6) --> max(x^2, x^4). Arguably, csign
71 ;; should be more semantically neutral--until it is, let's keep factor in here.
84 ;; Define a simplim%function to handle a limit of $max.
92 ;; True iff e is a GRE expression of the form max(...)
96 ;; True iff e is a GRE expression of the form min(...)
100 ;; Undone: max(1-x,1+x) - max(x,-x) --> 1. Doing this is possible--we could
101 ;; simplify max(1-x,1+x) --> 1 + max(-x,x).
103 ;; Return either the maximum of the arguments of l or a simplified max nounform.
104 ;; When max(a1,a2, ..., am) simplifies to a nounform, say max(b1, b2, ... , bn), each
105 ;; bk is *explicitly* one of a1 ... am. Much the same when max(a1,a2, ..., am)
106 ;; simplifies to a non-nounform; thus, for example,
107 ;;
108 ;; max(cos(x)^2+sin(x)^2,-107) --> cos(x)^2 + sin(x)^2 (*not* 1).
109 ;;
110 ;; Since coeff is purely syntatic, we have
111 ;;
112 ;; (%i1) xxx : 42*x^(cos(a)^2+sin(a)^2) + x^(-107)$
113 ;; (%i2) coeff(xxx,x, hipow(xxx,x));
114 ;; (%o2) 42
115 ;;
116 ;; And this is, I think, what users want.
123 ;; When maxmin_effort > 0, simplify each member of l and flatten (that is, do
124 ;; max(a,max(a,b)) -> max(a,b,c)). Additionally, we accumulate the largest real
125 ;; number. Since (mnump '$%i) --> false, we don't have to worry that num-max is
126 ;; complex. The effort for this step is O(n).
135 ;; Removing minf & -inf now results in things like max(minf, %i*inf)-->%i*inf.
142 ;; For inputs such as max(max(5,a), max(7,b), max(8,c), ...), the list
143 ;; l might contain many numbers. That could make some of the following
144 ;; simplifications slow. We could flag this case and redo the loop that
145 ;; looks for the largest number, or we could remove the mnump check from the first
146 ;; loop and separately make a loop that looks for the largest number.
148 ;; Sort and remove duplicates. The effort for this step is O(n logn)).
152 ;; When maxmin_effort > 2, if e and -e are members of l, replace e & -e by
153 ;; abs(e).
166 ;; Accumulate the maximum in the list acc. For each x in l, do:
167 ;; (a) if x is > or >= every member of acc, set acc to (x),
168 ;; (b) if x is < or <= to some member of acc, do nothing,
169 ;; (c) if neither 'a' or 'b', push x into acc,
170 ;; (d) if x cannot be compared to some member of acc, set issue-warning to true.
185 ;; When issue-warning is false and maxmin_effort > 2, use the betweenp
186 ;; simplification.
196 ;; Finally, do a few clean ups:
207 ;; Return -x, but check for the special cases x = inf, minf, und, ind, and infinity.
208 ;; Also locally set negdistrib to true (this is what the function neg does). We could
209 ;; do -zeroa -> zerob and -zerob -> zeroa, I suppose.
211 ;; To catch more cases, replace this function body with ($limit (mul -1 x)).
212 ;; But ($limit (mul -1 x)) misses the case -ind --> ind & limit can be costly.
219 ;; Define a simplim%function to handle a limit of $min.
229 ;; The function simp-min piggy-backs onto simp-max. That is, we use
230 ;; min(a,b,...) = -max(-a,-b,...).
237 ;; convert min(a, min(b,c)) --> min(a,b,c)
244 ;; Is the sort needed? I think so, but need a test that requires sorting...
249 ;; Several functions (derivdegree for example) use the maximin function. Here is
250 ;; a replacement that uses simp-min or simp-max.
260 ;; Return the narrowest comparison operator op (<, <=, =, >, >=) such that
261 ;; a op b evaluates to true. Return 'unknown' when either there is no such
262 ;; operator or when Maxima's sign function isn't powerful enough to determine
263 ;; such an operator; when Maxima is able to show that either argument is not
264 ;; real valued, return 'notcomparable.'
266 ;; The subtraction can be a problem--for example, compare(0.1, 1/10)
267 ;; evaluates to "=". But for flonum floats, I believe 0.1 > 1/10.
268 ;; If you want to convert flonum and big floats to exact rational
269 ;; numbers, use $rationalize.
271 ;; I think compare(asin(x), asin(x) + 1) should evaluate to < without
272 ;; being quizzed about the sign of x. Thus the call to lenient-extended-realp.
274 ;; Maxima's sign function doesn't know that zeroa > 0 and zerob < 0
277 ;; Simplify expressions with infinities. Without these checks, we can get odd
278 ;; questions such as "Is 1 zero or nonzero?"
288 ;; We'll make compare(infinity, infinity) -> noncompareable, but
289 ;; compare(inf, inf) -> "=". Simiarly, an expression involving ind or und
290 ;; is notcompareable.
295 ;; Attempt to catch expressions that are not extended real number valued. We
296 ;; don't want to do call csign on 1+und - und, and conclude that 1+und > und.
297 ;; The check lenient-extended-realp only looks at the main operator of the
298 ;; expression. Thus lenient-extended-realp flags a<b as not real valued,
299 ;; but it fails to flag 107*(a<b).
324 ;; When it's fairly likely that the real domain of e is nonempty, return true;
325 ;; otherwise, return false. Even if z has been declared complex, the real domain
326 ;; of z is nonempty; thus (lenient-extended-realp z) --> true. When does this
327 ;; function lie? One example is acos(abs(x) + 2). The real domain of this
328 ;; expression is empty, yet lenient-extended-realp returns true for this input.
341 ;; Convert all floats and big floats in e to an exact rational representation.