source: branches/streams/abcl/src/org/armedbear/lisp/MathFunctions.java

Last change on this file was 13440, checked in by ehuelsmann, 13 years ago

Rename writeToString() to printObject() since that's what it's being used for.
Additionally, create princToString() for use in error messages, making the

required replacement where appropriate.
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 29.0 KB
Line 
1/*
2 * MathFunctions.java
3 *
4 * Copyright (C) 2004-2006 Peter Graves
5 * $Id: MathFunctions.java 13440 2011-08-05 21:25:10Z ehuelsmann $
6 *
7 * This program is free software; you can redistribute it and/or
8 * modify it under the terms of the GNU General Public License
9 * as published by the Free Software Foundation; either version 2
10 * of the License, or (at your option) any later version.
11 *
12 * This program is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
20 *
21 * As a special exception, the copyright holders of this library give you
22 * permission to link this library with independent modules to produce an
23 * executable, regardless of the license terms of these independent
24 * modules, and to copy and distribute the resulting executable under
25 * terms of your choice, provided that you also meet, for each linked
26 * independent module, the terms and conditions of the license of that
27 * module.  An independent module is a module which is not derived from
28 * or based on this library.  If you modify this library, you may extend
29 * this exception to your version of the library, but you are not
30 * obligated to do so.  If you do not wish to do so, delete this
31 * exception statement from your version.
32 */
33
34package org.armedbear.lisp;
35
36import static org.armedbear.lisp.Lisp.*;
37
38public final class MathFunctions
39{
40
41    // Implementation of section 12.1.5.3, which says:
42    // "If the result of any computation would be a complex number whose
43    //  real part is of type rational and whose imaginary part is zero,
44    //  the result is converted to the rational which is the real part."
45    private static final LispObject complexToRealFixup(LispObject result,
46                                                       LispObject arg)
47    {
48        if (result instanceof Complex
49            && ! (arg instanceof Complex)) {
50            Complex c = (Complex)result;
51            LispObject im = c.getImaginaryPart();
52            if (im.zerop())
53                return c.getRealPart();
54        }
55        return result;
56    }
57
58    // ### sin
59    private static final Primitive SIN = new Primitive("sin", "radians")
60    {
61        @Override
62        public LispObject execute(LispObject arg)
63        {
64            return sin(arg);
65        }
66    };
67
68    static LispObject sin(LispObject arg)
69    {
70        if (arg instanceof DoubleFloat)
71            return new DoubleFloat(Math.sin(((DoubleFloat)arg).value));
72        if (arg.realp())
73            return new SingleFloat((float)Math.sin(SingleFloat.coerceToFloat(arg).value));
74        if (arg instanceof Complex) {
75            LispObject n = arg.multiplyBy(Complex.getInstance(Fixnum.ZERO,
76                                                              Fixnum.ONE));
77            LispObject result = exp(n);
78            result = result.subtract(exp(n.multiplyBy(Fixnum.MINUS_ONE)));
79            return result.divideBy(Fixnum.TWO.multiplyBy(Complex.getInstance(Fixnum.ZERO,
80                                                                             Fixnum.ONE)));
81        }
82        return type_error(arg, Symbol.NUMBER);
83    }
84
85    // ### cos
86    private static final Primitive COS = new Primitive("cos", "radians")
87    {
88        @Override
89        public LispObject execute(LispObject arg)
90        {
91            return cos(arg);
92        }
93    };
94
95    static LispObject cos(LispObject arg)
96    {
97        if (arg instanceof DoubleFloat)
98            return new DoubleFloat(Math.cos(((DoubleFloat)arg).value));
99        if (arg.realp())
100            return new SingleFloat((float)Math.cos(SingleFloat.coerceToFloat(arg).value));
101        if (arg instanceof Complex) {
102            LispObject n = arg.multiplyBy(Complex.getInstance(Fixnum.ZERO,
103                                                              Fixnum.ONE));
104            LispObject result = exp(n);
105            result = result.add(exp(n.multiplyBy(Fixnum.MINUS_ONE)));
106            return result.divideBy(Fixnum.TWO);
107        }
108        return type_error(arg, Symbol.NUMBER);
109    }
110
111    // ### tan
112    private static final Primitive TAN = new Primitive("tan", "radians")
113    {
114        @Override
115        public LispObject execute(LispObject arg)
116        {
117            if (arg instanceof DoubleFloat)
118                return new DoubleFloat(Math.tan(((DoubleFloat)arg).value));
119            if (arg.realp())
120                return new SingleFloat((float)Math.tan(SingleFloat.coerceToFloat(arg).value));
121            return sin(arg).divideBy(cos(arg));
122        }
123    };
124
125    // ### asin
126    private static final Primitive ASIN = new Primitive("asin", "number")
127    {
128        @Override
129        public LispObject execute(LispObject arg)
130        {
131            return asin(arg);
132        }
133    };
134
135    static LispObject asin(LispObject arg)
136    {
137        if (arg instanceof SingleFloat) {
138            float f = ((SingleFloat)arg).value;
139            if (Math.abs(f) <= 1)
140                return new SingleFloat((float)Math.asin(f));
141        }
142        if (arg instanceof DoubleFloat) {
143            double d = ((DoubleFloat)arg).value;
144            if (Math.abs(d) <= 1)
145                return new DoubleFloat(Math.asin(d));
146        }
147        LispObject result = arg.multiplyBy(arg);
148        result = Fixnum.ONE.subtract(result);
149        result = sqrt(result);
150        LispObject n = Complex.getInstance(Fixnum.ZERO, Fixnum.ONE);
151        n = n.multiplyBy(arg);
152        result = n.add(result);
153        result = log(result);
154        result = result.multiplyBy(Complex.getInstance(Fixnum.ZERO,
155                                                       Fixnum.MINUS_ONE));
156
157        return complexToRealFixup(result, arg);
158    }
159
160    // ### acos
161    private static final Primitive ACOS = new Primitive("acos", "number")
162    {
163        @Override
164        public LispObject execute(LispObject arg)
165        {
166            return acos(arg);
167        }
168    };
169
170    static LispObject acos(LispObject arg)
171    {
172        if (arg instanceof DoubleFloat) {
173            double d = ((DoubleFloat)arg).value;
174            if (Math.abs(d) <= 1)
175                return new DoubleFloat(Math.acos(d));
176        }
177        if (arg instanceof SingleFloat) {
178            float f = ((SingleFloat)arg).value;
179            if (Math.abs(f) <= 1)
180                return new SingleFloat((float)Math.acos(f));
181        }
182        LispObject result = new DoubleFloat(Math.PI/2);
183        if (!(arg instanceof DoubleFloat)) {
184            if (arg instanceof Complex &&
185                    ((Complex)arg).getRealPart() instanceof DoubleFloat) {
186                    // do nothing; we want to keep the double float value
187            }
188            else
189                result = new SingleFloat((float)((DoubleFloat)result).value);
190        }
191        result = result.subtract(asin(arg));
192
193        return complexToRealFixup(result, arg);
194    }
195
196    // ### atan
197    private static final Primitive ATAN =
198        new Primitive("atan", "number1 &optional number2")
199    {
200        @Override
201        public LispObject execute(LispObject arg)
202        {
203            if (arg.numberp())
204                return atan(arg);
205            return type_error(arg, Symbol.NUMBER);
206        }
207
208        // "If both number1 and number2 are supplied for atan, the result is
209        // the arc tangent of number1/number2."
210
211        // y = +0     x = +0       +0
212        // y = -0     x = +0       -0
213        // y = +0     x = -0       +<PI>
214        // y = -0     x = -0       -<PI>
215        @Override
216        public LispObject execute(LispObject y, LispObject x)
217
218        {
219            if (!y.realp())
220                return type_error(y, Symbol.REAL);
221            if (!x.realp())
222                return type_error(x, Symbol.REAL);
223            double d1, d2;
224            d1 = DoubleFloat.coerceToFloat(y).value;
225            d2 = DoubleFloat.coerceToFloat(x).value;
226            double result = Math.atan2(d1, d2);
227            if (y instanceof DoubleFloat || x instanceof DoubleFloat)
228                return new DoubleFloat(result);
229            else
230                return new SingleFloat((float)result);
231        }
232    };
233
234    static LispObject atan(LispObject arg)
235    {
236        if (arg instanceof Complex) {
237            LispObject im = ((Complex)arg).imagpart;
238            if (im.zerop())
239                return Complex.getInstance(atan(((Complex)arg).realpart),
240                                           im);
241            LispObject result = arg.multiplyBy(arg);
242            result = result.add(Fixnum.ONE);
243            result = Fixnum.ONE.divideBy(result);
244            result = sqrt(result);
245            LispObject n = Complex.getInstance(Fixnum.ZERO, Fixnum.ONE);
246            n = n.multiplyBy(arg);
247            n = n.add(Fixnum.ONE);
248            result = n.multiplyBy(result);
249            result = log(result);
250            result = result.multiplyBy(Complex.getInstance(Fixnum.ZERO, Fixnum.MINUS_ONE));
251            return result;
252        }
253        if (arg instanceof DoubleFloat)
254            return new DoubleFloat(Math.atan(((DoubleFloat)arg).value));
255        return new SingleFloat((float)Math.atan(SingleFloat.coerceToFloat(arg).value));
256    }
257
258    // ### sinh
259    private static final Primitive SINH = new Primitive("sinh", "number")
260    {
261        @Override
262        public LispObject execute(LispObject arg)
263        {
264            return sinh(arg);
265        }
266    };
267
268    static LispObject sinh(LispObject arg)
269    {
270        if (arg instanceof Complex) {
271            LispObject im = ((Complex)arg).getImaginaryPart();
272            if (im.zerop())
273                return Complex.getInstance(sinh(((Complex)arg).getRealPart()),
274                                           im);
275        }
276        if (arg instanceof SingleFloat) {
277            double d = Math.sinh(((SingleFloat)arg).value);
278            return new SingleFloat((float)d);
279        } else if (arg instanceof DoubleFloat) {
280            double d = Math.sinh(((DoubleFloat)arg).value);
281            return new DoubleFloat(d);
282        }
283        LispObject result = exp(arg);
284        result = result.subtract(exp(arg.multiplyBy(Fixnum.MINUS_ONE)));
285        result = result.divideBy(Fixnum.TWO);
286
287        return complexToRealFixup(result, arg);
288    }
289
290    // ### cosh
291    private static final Primitive COSH = new Primitive("cosh", "number")
292    {
293        @Override
294        public LispObject execute(LispObject arg)
295        {
296            return cosh(arg);
297        }
298    };
299
300    static LispObject cosh(LispObject arg)
301    {
302        if (arg instanceof Complex) {
303            LispObject im = ((Complex)arg).getImaginaryPart();
304            if (im.zerop())
305                return Complex.getInstance(cosh(((Complex)arg).getRealPart()),
306                                           im);
307        }
308        if (arg instanceof SingleFloat) {
309            double d = Math.cosh(((SingleFloat)arg).value);
310            return new SingleFloat((float)d);
311        } else if (arg instanceof DoubleFloat) {
312            double d = Math.cosh(((DoubleFloat)arg).value);
313            return new DoubleFloat(d);
314        }
315        LispObject result = exp(arg);
316        result = result.add(exp(arg.multiplyBy(Fixnum.MINUS_ONE)));
317        result = result.divideBy(Fixnum.TWO);
318
319        return complexToRealFixup(result, arg);
320    }
321
322    // ### tanh
323    private static final Primitive TANH = new Primitive("tanh", "number")
324    {
325        @Override
326        public LispObject execute(LispObject arg)
327        {
328            if (arg instanceof SingleFloat) {
329                double d = Math.tanh(((SingleFloat)arg).value);
330                return new SingleFloat((float)d);
331            } else if (arg instanceof DoubleFloat) {
332                double d = Math.tanh(((DoubleFloat)arg).value);
333                return new DoubleFloat(d);
334            }
335            return sinh(arg).divideBy(cosh(arg));
336        }
337    };
338
339    // ### asinh
340    private static final Primitive ASINH = new Primitive("asinh", "number")
341    {
342        @Override
343        public LispObject execute(LispObject arg)
344        {
345            return asinh(arg);
346        }
347    };
348
349    static LispObject asinh(LispObject arg)
350    {
351        if (arg instanceof Complex) {
352            LispObject im = ((Complex)arg).getImaginaryPart();
353            if (im.zerop())
354                return Complex.getInstance(asinh(((Complex)arg).getRealPart()),
355                                           im);
356        }
357        LispObject result = arg.multiplyBy(arg);
358        result = Fixnum.ONE.add(result);
359        result = sqrt(result);
360        result = result.add(arg);
361        result = log(result);
362
363        return complexToRealFixup(result, arg);
364    }
365
366    // ### acosh
367    private static final Primitive ACOSH = new Primitive("acosh", "number")
368    {
369        @Override
370        public LispObject execute(LispObject arg)
371        {
372            return acosh(arg);
373        }
374    };
375
376    static LispObject acosh(LispObject arg)
377    {
378        if (arg instanceof Complex) {
379            LispObject im = ((Complex)arg).getImaginaryPart();
380            if (im.zerop())
381                return Complex.getInstance(acosh(((Complex)arg).getRealPart()),
382                                           im);
383        }
384        LispObject n1 = arg.add(Fixnum.ONE);
385        n1 = n1.divideBy(Fixnum.TWO);
386        n1 = sqrt(n1);
387        LispObject n2 = arg.subtract(Fixnum.ONE);
388        n2 = n2.divideBy(Fixnum.TWO);
389        n2 = sqrt(n2);
390        LispObject result = n1.add(n2);
391        result = log(result);
392        result = result.multiplyBy(Fixnum.TWO);
393
394        return complexToRealFixup(result, arg);
395    }
396
397    // ### atanh
398    private static final Primitive ATANH = new Primitive("atanh", "number")
399    {
400        @Override
401        public LispObject execute(LispObject arg)
402        {
403            return atanh(arg);
404        }
405    };
406
407    static LispObject atanh(LispObject arg)
408    {
409        if (arg instanceof Complex) {
410            LispObject im = ((Complex)arg).getImaginaryPart();
411            if (im.zerop())
412                return Complex.getInstance(atanh(((Complex)arg).getRealPart()),
413                                           im);
414        }
415        LispObject n1 = log(Fixnum.ONE.add(arg));
416        LispObject n2 = log(Fixnum.ONE.subtract(arg));
417        LispObject result = n1.subtract(n2);
418        result = result.divideBy(Fixnum.TWO);
419
420        return complexToRealFixup(result, arg);
421    }
422
423    // ### cis
424    private static final Primitive CIS = new Primitive("cis", "radians")
425    {
426        @Override
427        public LispObject execute(LispObject arg)
428        {
429            return cis(arg);
430        }
431    };
432
433    static LispObject cis(LispObject arg)
434    {
435        if (arg.realp())
436            return Complex.getInstance(cos(arg), sin(arg));
437        return type_error(arg, Symbol.REAL);
438    }
439
440    // ### exp
441    private static final Primitive EXP = new Primitive("exp", "number")
442    {
443        @Override
444        public LispObject execute(LispObject arg)
445        {
446            return exp(arg);
447        }
448    };
449
450    static LispObject exp(LispObject arg)
451    {
452        if (arg.realp()) {
453            if (arg instanceof DoubleFloat) {
454                double d = Math.pow(Math.E, ((DoubleFloat)arg).value);
455                return OverUnderFlowCheck(new DoubleFloat(d));
456            } else {
457                float f = (float) Math.pow(Math.E, SingleFloat.coerceToFloat(arg).value);
458                return OverUnderFlowCheck(new SingleFloat(f));
459            }
460        }
461        if (arg instanceof Complex) {
462            Complex c = (Complex) arg;
463            return exp(c.getRealPart()).multiplyBy(cis(c.getImaginaryPart()));
464        }
465        return type_error(arg, Symbol.NUMBER);
466    }
467
468    // ### sqrt
469    private static final Primitive SQRT = new Primitive("sqrt", "number")
470    {
471        @Override
472        public LispObject execute(LispObject arg)
473        {
474            return sqrt(arg);
475        }
476    };
477
478    static final LispObject sqrt(LispObject obj)
479    {
480        if (obj instanceof DoubleFloat) {
481            if (obj.minusp())
482                return Complex.getInstance(new DoubleFloat(0), sqrt(obj.negate()));
483            return new DoubleFloat(Math.sqrt(DoubleFloat.coerceToFloat(obj).value));
484        }
485        if (obj.realp()) {
486            if (obj.minusp())
487                return Complex.getInstance(new SingleFloat(0), sqrt(obj.negate()));
488            return new SingleFloat((float)Math.sqrt(SingleFloat.coerceToFloat(obj).value));
489        }
490        if (obj instanceof Complex) {
491            LispObject imagpart = ((Complex)obj).imagpart;
492            if (imagpart.zerop()) {
493                LispObject realpart = ((Complex)obj).realpart;
494                if (realpart.minusp())
495                    return Complex.getInstance(imagpart, sqrt(realpart.negate()));
496                else
497                    return Complex.getInstance(sqrt(realpart), imagpart);
498            }
499            return exp(log(obj).divideBy(Fixnum.TWO));
500        }
501        return type_error(obj, Symbol.NUMBER);
502    }
503
504    // ### log
505    private static final Primitive LOG =
506        new Primitive("log", "number &optional base")
507    {
508        @Override
509        public LispObject execute(LispObject arg)
510        {
511            return log(arg);
512        }
513        @Override
514        public LispObject execute(LispObject number, LispObject base)
515
516        {
517            if (number.realp() && !number.minusp()
518                && base.isEqualTo(Fixnum.getInstance(10))) {
519                double d =
520                    Math.log10(DoubleFloat.coerceToFloat(number).value);
521                if (number instanceof DoubleFloat
522                    || base instanceof DoubleFloat)
523                    return new DoubleFloat(d);
524                else
525                    return new SingleFloat((float)d);
526            }
527            return log(number).divideBy(log(base));
528        }
529    };
530
531    static final LispObject log(LispObject obj)
532    {
533        if (obj.realp() && !obj.minusp()) {
534            // Result is real.
535            if (obj instanceof Fixnum)
536                return new SingleFloat((float)Math.log(((Fixnum)obj).value));
537            if (obj instanceof Bignum)
538                return new SingleFloat((float)Math.log(((Bignum)obj).doubleValue()));
539            if (obj instanceof Ratio)
540                return new SingleFloat((float)Math.log(((Ratio)obj).doubleValue()));
541            if (obj instanceof SingleFloat)
542                return new SingleFloat((float)Math.log(((SingleFloat)obj).value));
543            if (obj instanceof DoubleFloat)
544                return new DoubleFloat(Math.log(((DoubleFloat)obj).value));
545        } else {
546            // Result is complex.
547            if (obj.realp() && obj.minusp()) {
548                if (obj instanceof DoubleFloat) {
549                    DoubleFloat re = DoubleFloat.coerceToFloat(obj);
550                    DoubleFloat abs = new DoubleFloat(Math.abs(re.value));
551                    DoubleFloat phase = new DoubleFloat(Math.PI);
552                    return Complex.getInstance(new DoubleFloat(Math.log(abs.getValue())), phase);
553                } else {
554                    SingleFloat re = SingleFloat.coerceToFloat(obj);
555                    SingleFloat abs = new SingleFloat(Math.abs(re.value));
556                    SingleFloat phase = new SingleFloat((float)Math.PI);
557                    return Complex.getInstance(new SingleFloat((float)Math.log(abs.value)), phase);
558                }
559            } else if (obj instanceof Complex) {
560                if (((Complex)obj).getRealPart() instanceof DoubleFloat) {
561                    DoubleFloat re = DoubleFloat.coerceToFloat(((Complex)obj).getRealPart());
562                    DoubleFloat im = DoubleFloat.coerceToFloat(((Complex)obj).getImaginaryPart());
563                    DoubleFloat phase =
564                        new DoubleFloat(Math.atan2(im.getValue(), re.getValue()));  // atan(y/x)
565                    DoubleFloat abs = DoubleFloat.coerceToFloat(obj.ABS());
566                    return Complex.getInstance(new DoubleFloat(Math.log(abs.getValue())), phase);
567                } else {
568                    SingleFloat re = SingleFloat.coerceToFloat(((Complex)obj).getRealPart());
569                    SingleFloat im = SingleFloat.coerceToFloat(((Complex)obj).getImaginaryPart());
570                    SingleFloat phase =
571                        new SingleFloat((float)Math.atan2(im.value, re.value));  // atan(y/x)
572                    SingleFloat abs = SingleFloat.coerceToFloat(obj.ABS());
573                    return Complex.getInstance(new SingleFloat((float)Math.log(abs.value)), phase);
574                }
575            }
576        }
577        type_error(obj, Symbol.NUMBER);
578        return NIL;
579    }
580
581    // ### expt base-number power-number => result
582    public static final Primitive EXPT =
583        new Primitive("expt", "base-number power-number")
584    {
585        @Override
586        public LispObject execute(LispObject base, LispObject power)
587
588        {
589            if (power.zerop()) {
590                if (power instanceof Fixnum) {
591                    if (base instanceof SingleFloat)
592                        return SingleFloat.ONE;
593                    if (base instanceof DoubleFloat)
594                        return DoubleFloat.ONE;
595                    if (base instanceof Complex) {
596                        if (((Complex)base).realpart instanceof SingleFloat)
597                            return Complex.getInstance(SingleFloat.ONE,
598                                                       SingleFloat.ZERO);
599                        if (((Complex)base).realpart instanceof DoubleFloat)
600                            return Complex.getInstance(DoubleFloat.ONE,
601                                                       DoubleFloat.ZERO);
602                    }
603                    return Fixnum.ONE;
604                }
605                if (power instanceof DoubleFloat)
606                    return DoubleFloat.ONE;
607                if (base instanceof DoubleFloat)
608                    return DoubleFloat.ONE;
609                return SingleFloat.ONE;
610            }
611            if (base.zerop())
612                return base;
613            if (base.isEqualTo(1))
614                return base;
615           
616            if ((power instanceof Fixnum
617                 || power instanceof Bignum)
618                 && (base.rationalp()
619                     || (base instanceof Complex
620                         && ((Complex)base).realpart.rationalp()))) {
621                // exact math version
622                return intexp(base, power);
623            }
624            // for anything not a rational or complex rational, use
625            // float approximation.
626            boolean wantDoubleFloat = false;
627            if (base instanceof DoubleFloat)
628                wantDoubleFloat = true;
629            else if (power instanceof DoubleFloat)
630                wantDoubleFloat = true;
631            else if (base instanceof Complex
632                     && (((Complex)base).getRealPart() instanceof DoubleFloat
633                         || ((Complex)base).getImaginaryPart() instanceof DoubleFloat))
634                wantDoubleFloat = true;
635            else if (power instanceof Complex
636                    && (((Complex)power).getRealPart() instanceof DoubleFloat
637                         || ((Complex)power).getImaginaryPart() instanceof DoubleFloat))
638                wantDoubleFloat = true;
639
640            if (wantDoubleFloat) {
641                if (power instanceof Complex)
642                    power = ((Complex)power).coerceToDoubleFloat();
643                else
644                    power = DoubleFloat.coerceToFloat(power);
645
646                if (base instanceof Complex)
647                    base = ((Complex)base).coerceToDoubleFloat();
648                else
649                    base = DoubleFloat.coerceToFloat(base);
650            }
651
652
653
654            if (base instanceof Complex || power instanceof Complex)
655                return exp(power.multiplyBy(log(base)));
656            final double x; // base
657            final double y; // power
658            if (base instanceof Fixnum)
659                x = ((Fixnum)base).value;
660            else if (base instanceof Bignum)
661                x = ((Bignum)base).doubleValue();
662            else if (base instanceof Ratio)
663                x = ((Ratio)base).doubleValue();
664            else if (base instanceof SingleFloat)
665                x = ((SingleFloat)base).value;
666            else if (base instanceof DoubleFloat)
667                x = ((DoubleFloat)base).value;
668            else
669                return error(new LispError("EXPT: unsupported case: base is of type " +
670                                            base.typeOf().princToString()));
671
672            if (power instanceof Fixnum)
673                y = ((Fixnum)power).value;
674            else if (power instanceof Bignum)
675                y = ((Bignum)power).doubleValue();
676            else if (power instanceof Ratio)
677                y = ((Ratio)power).doubleValue();
678            else if (power instanceof SingleFloat)
679                y = ((SingleFloat)power).value;
680            else if (power instanceof DoubleFloat)
681                y = ((DoubleFloat)power).value;
682            else
683                return error(new LispError("EXPT: unsupported case: power is of type " +
684                                            power.typeOf().princToString()));
685            double r = Math.pow(x, y);
686            if (Double.isNaN(r)) {
687                if (x < 0) {
688                    r = Math.pow(-x, y);
689                    double realPart = r * Math.cos(y * Math.PI);
690                    double imagPart = r * Math.sin(y * Math.PI);
691                    if (base instanceof DoubleFloat || power instanceof DoubleFloat)
692                        return Complex
693                            .getInstance(OverUnderFlowCheck(new DoubleFloat(realPart)),
694                                         OverUnderFlowCheck(new DoubleFloat(imagPart)));
695                    else
696                        return Complex
697                            .getInstance(OverUnderFlowCheck(new SingleFloat((float)realPart)),
698                                         OverUnderFlowCheck(new SingleFloat((float)imagPart)));
699                }
700            }
701            if (base instanceof DoubleFloat || power instanceof DoubleFloat)
702                return OverUnderFlowCheck(new DoubleFloat(r));
703            else
704                return OverUnderFlowCheck(new SingleFloat((float)r));
705        }
706    };
707
708    /** Checks number for over- or underflow values.
709     *
710     * @param number
711     * @return number or signals an appropriate error
712     */
713    final static LispObject OverUnderFlowCheck(LispObject number)
714
715    {
716        if (number instanceof Complex) {
717            OverUnderFlowCheck(((Complex)number).realpart);
718            OverUnderFlowCheck(((Complex)number).imagpart);
719            return number;
720        }
721
722        if (TRAP_OVERFLOW) {
723            if (number instanceof SingleFloat)
724                if (Float.isInfinite(((SingleFloat)number).value))
725                    return error(new FloatingPointOverflow(NIL));
726            if (number instanceof DoubleFloat)
727                if (Double.isInfinite(((DoubleFloat)number).value))
728                    return error(new FloatingPointOverflow(NIL));
729        }
730        if (TRAP_UNDERFLOW) {
731            if (number.zerop())
732                return error(new FloatingPointUnderflow(NIL));
733        }
734        return number;
735    }
736
737    /** Checks number for over- or underflow values.
738     *
739     * @param number
740     * @return number or signals an appropriate error
741     */
742    final static float OverUnderFlowCheck(float number)
743
744    {
745        if (TRAP_OVERFLOW) {
746            if (Float.isInfinite(number))
747                error(new FloatingPointOverflow(NIL));
748        }
749        if (TRAP_UNDERFLOW) {
750            if (number == 0)
751                error(new FloatingPointUnderflow(NIL));
752        }
753        return number;
754    }
755
756    /** Checks number for over- or underflow values.
757     *
758     * @param number
759     * @return number or signals an appropriate error
760     */
761    public final static double OverUnderFlowCheck(double number)
762
763    {
764        if (TRAP_OVERFLOW) {
765            if (Double.isInfinite(number))
766                error(new FloatingPointOverflow(NIL));
767        }
768        if (TRAP_UNDERFLOW) {
769            if (number == 0)
770                error(new FloatingPointUnderflow(NIL));
771        }
772        return number;
773    }
774    // Adapted from SBCL.
775    /** Return the exponent of base taken to the integer exponent power
776     *
777     * @param base A value of any type
778     * @param power An integer (fixnum or bignum) value
779     */
780    static final LispObject intexp(LispObject base, LispObject power)
781
782    {
783        if (power.isEqualTo(0))
784            return Fixnum.ONE;
785        if (base.isEqualTo(1))
786            return base;
787        if (base.isEqualTo(0))
788            return base;
789
790        if (power.minusp()) {
791            power = Fixnum.ZERO.subtract(power);
792            return Fixnum.ONE.divideBy(intexp(base, power));
793        }
794        if (base.eql(Fixnum.TWO))
795            return Fixnum.ONE.ash(power);
796
797        LispObject nextn = power.ash(Fixnum.MINUS_ONE);
798        LispObject total;
799        if (power.oddp())
800            total = base;
801        else
802            total = Fixnum.ONE;
803        while (true) {
804            if (nextn.zerop())
805                return total;
806            base = base.multiplyBy(base);
807
808            if (nextn.oddp())
809                total = base.multiplyBy(total);
810            nextn = nextn.ash(Fixnum.MINUS_ONE);
811        }
812    }
813}
Note: See TracBrowser for help on using the repository browser.