- introduced info_flags::cinteger, info_flags::crational,
[ginac.git] / ginac / inifcns_trans.cpp
1 /** @file inifcns_trans.cpp
2  *
3  *  Implementation of transcendental (and trigonometric and hyperbolic)
4  *  functions. */
5
6 /*
7  *  GiNaC Copyright (C) 1999 Johannes Gutenberg University Mainz, Germany
8  *
9  *  This program is free software; you can redistribute it and/or modify
10  *  it under the terms of the GNU General Public License as published by
11  *  the Free Software Foundation; either version 2 of the License, or
12  *  (at your option) any later version.
13  *
14  *  This program is distributed in the hope that it will be useful,
15  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17  *  GNU General Public License for more details.
18  *
19  *  You should have received a copy of the GNU General Public License
20  *  along with this program; if not, write to the Free Software
21  *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
22  */
23
24 #include <vector>
25 #include <stdexcept>
26
27 #include "inifcns.h"
28 #include "ex.h"
29 #include "constant.h"
30 #include "numeric.h"
31 #include "power.h"
32 #include "relational.h"
33 #include "symbol.h"
34 #include "utils.h"
35
36 #ifndef NO_GINAC_NAMESPACE
37 namespace GiNaC {
38 #endif // ndef NO_GINAC_NAMESPACE
39
40 //////////
41 // exponential function
42 //////////
43
44 static ex exp_evalf(ex const & x)
45 {
46     BEGIN_TYPECHECK
47         TYPECHECK(x,numeric)
48     END_TYPECHECK(exp(x))
49     
50     return exp(ex_to_numeric(x)); // -> numeric exp(numeric)
51 }
52
53 static ex exp_eval(ex const & x)
54 {
55     // exp(0) -> 1
56     if (x.is_zero()) {
57         return exONE();
58     }
59     // exp(n*Pi*I/2) -> {+1|+I|-1|-I}
60     ex TwoExOverPiI=(2*x)/(Pi*I);
61     if (TwoExOverPiI.info(info_flags::integer)) {
62         numeric z=mod(ex_to_numeric(TwoExOverPiI),numeric(4));
63         if (z.is_equal(numZERO()))
64             return exONE();
65         if (z.is_equal(numONE()))
66             return ex(I);
67         if (z.is_equal(numTWO()))
68             return exMINUSONE();
69         if (z.is_equal(numTHREE()))
70             return ex(-I);
71     }
72     // exp(log(x)) -> x
73     if (is_ex_the_function(x, log))
74         return x.op(0);
75     
76     // exp(float)
77     if (x.info(info_flags::numeric) && !x.info(info_flags::crational))
78         return exp_evalf(x);
79     
80     return exp(x).hold();
81 }
82
83 static ex exp_diff(ex const & x, unsigned diff_param)
84 {
85     GINAC_ASSERT(diff_param==0);
86
87     return exp(x);
88 }
89
90 REGISTER_FUNCTION(exp, exp_eval, exp_evalf, exp_diff, NULL);
91
92 //////////
93 // natural logarithm
94 //////////
95
96 static ex log_evalf(ex const & x)
97 {
98     BEGIN_TYPECHECK
99         TYPECHECK(x,numeric)
100     END_TYPECHECK(log(x))
101     
102     return log(ex_to_numeric(x)); // -> numeric log(numeric)
103 }
104
105 static ex log_eval(ex const & x)
106 {
107     if (x.info(info_flags::numeric)) {
108         // log(1) -> 0
109         if (x.is_equal(exONE()))
110             return exZERO();
111         // log(-1) -> I*Pi
112         if (x.is_equal(exMINUSONE()))
113             return (I*Pi);
114         // log(I) -> Pi*I/2
115         if (x.is_equal(I))
116             return (I*Pi*numeric(1,2));
117         // log(-I) -> -Pi*I/2
118         if (x.is_equal(-I))
119             return (I*Pi*numeric(-1,2));
120         // log(0) -> throw singularity
121         if (x.is_equal(exZERO()))
122             throw(std::domain_error("log_eval(): log(0)"));
123         // log(float)
124         if (!x.info(info_flags::crational))
125             return log_evalf(x);
126     }
127     // log(exp(t)) -> t (for real-valued t):
128     if (is_ex_the_function(x, exp)) {
129         ex t=x.op(0);
130         if (t.info(info_flags::real))
131             return t;
132     }
133     
134     return log(x).hold();
135 }
136
137 static ex log_diff(ex const & x, unsigned diff_param)
138 {
139     GINAC_ASSERT(diff_param==0);
140
141     return power(x, -1);
142 }
143
144 REGISTER_FUNCTION(log, log_eval, log_evalf, log_diff, NULL);
145
146 //////////
147 // sine (trigonometric function)
148 //////////
149
150 static ex sin_evalf(ex const & x)
151 {
152     BEGIN_TYPECHECK
153        TYPECHECK(x,numeric)
154     END_TYPECHECK(sin(x))
155     
156     return sin(ex_to_numeric(x)); // -> numeric sin(numeric)
157 }
158
159 static ex sin_eval(ex const & x)
160 {
161     ex xOverPi=x/Pi;
162     if (xOverPi.info(info_flags::numeric)) {
163         // sin(n*Pi) -> 0
164         if (xOverPi.info(info_flags::integer))
165             return exZERO();
166         
167         // sin((2n+1)*Pi/2) -> {+|-}1
168         ex xOverPiMinusHalf=xOverPi-exHALF();
169         if (xOverPiMinusHalf.info(info_flags::even))
170             return exONE();
171         else if (xOverPiMinusHalf.info(info_flags::odd))
172             return exMINUSONE();
173     }
174     
175     if (is_ex_exactly_of_type(x, function)) {
176         ex t=x.op(0);
177         // sin(asin(x)) -> x
178         if (is_ex_the_function(x, asin))
179             return t;
180         // sin(acos(x)) -> (1-x^2)^(1/2)
181         if (is_ex_the_function(x, acos))
182             return power(exONE()-power(t,exTWO()),exHALF());
183         // sin(atan(x)) -> x*(1+x^2)^(-1/2)
184         if (is_ex_the_function(x, atan))
185             return t*power(exONE()+power(t,exTWO()),exMINUSHALF());
186     }
187     
188     // sin(float) -> float
189     if (x.info(info_flags::numeric) && !x.info(info_flags::crational))
190         return sin_evalf(x);
191     
192     return sin(x).hold();
193 }
194
195 static ex sin_diff(ex const & x, unsigned diff_param)
196 {
197     GINAC_ASSERT(diff_param==0);
198     
199     return cos(x);
200 }
201
202 REGISTER_FUNCTION(sin, sin_eval, sin_evalf, sin_diff, NULL);
203
204 //////////
205 // cosine (trigonometric function)
206 //////////
207
208 static ex cos_evalf(ex const & x)
209 {
210     BEGIN_TYPECHECK
211         TYPECHECK(x,numeric)
212     END_TYPECHECK(cos(x))
213     
214     return cos(ex_to_numeric(x)); // -> numeric cos(numeric)
215 }
216
217 static ex cos_eval(ex const & x)
218 {
219     ex xOverPi=x/Pi;
220     if (xOverPi.info(info_flags::numeric)) {
221         // cos(n*Pi) -> {+|-}1
222         if (xOverPi.info(info_flags::even))
223             return exONE();
224         else if (xOverPi.info(info_flags::odd))
225             return exMINUSONE();
226         
227         // cos((2n+1)*Pi/2) -> 0
228         ex xOverPiMinusHalf=xOverPi-exHALF();
229         if (xOverPiMinusHalf.info(info_flags::integer))
230             return exZERO();
231     }
232     
233     if (is_ex_exactly_of_type(x, function)) {
234         ex t=x.op(0);
235         // cos(acos(x)) -> x
236         if (is_ex_the_function(x, acos))
237             return t;
238         // cos(asin(x)) -> (1-x^2)^(1/2)
239         if (is_ex_the_function(x, asin))
240             return power(exONE()-power(t,exTWO()),exHALF());
241         // cos(atan(x)) -> (1+x^2)^(-1/2)
242         if (is_ex_the_function(x, atan))
243             return power(exONE()+power(t,exTWO()),exMINUSHALF());
244     }
245     
246     // cos(float) -> float
247     if (x.info(info_flags::numeric) && !x.info(info_flags::crational))
248         return cos_evalf(x);
249     
250     return cos(x).hold();
251 }
252
253 static ex cos_diff(ex const & x, unsigned diff_param)
254 {
255     GINAC_ASSERT(diff_param==0);
256
257     return numMINUSONE()*sin(x);
258 }
259
260 REGISTER_FUNCTION(cos, cos_eval, cos_evalf, cos_diff, NULL);
261
262 //////////
263 // tangent (trigonometric function)
264 //////////
265
266 static ex tan_evalf(ex const & x)
267 {
268     BEGIN_TYPECHECK
269        TYPECHECK(x,numeric)
270     END_TYPECHECK(tan(x)) // -> numeric tan(numeric)
271     
272     return tan(ex_to_numeric(x));
273 }
274
275 static ex tan_eval(ex const & x)
276 {
277     // tan(n*Pi/3) -> {0|3^(1/2)|-(3^(1/2))}
278     ex ThreeExOverPi=numTHREE()*x/Pi;
279     if (ThreeExOverPi.info(info_flags::integer)) {
280         numeric z=mod(ex_to_numeric(ThreeExOverPi),numeric(3));
281         if (z.is_equal(numZERO()))
282             return exZERO();
283         if (z.is_equal(numONE()))
284             return power(exTHREE(),exHALF());
285         if (z.is_equal(numTWO()))
286             return -power(exTHREE(),exHALF());
287     }
288     
289     // tan((2n+1)*Pi/2) -> throw
290     ex ExOverPiMinusHalf=x/Pi-exHALF();
291     if (ExOverPiMinusHalf.info(info_flags::integer))
292         throw (std::domain_error("tan_eval(): infinity"));
293     
294     if (is_ex_exactly_of_type(x, function)) {
295         ex t=x.op(0);
296         // tan(atan(x)) -> x
297         if (is_ex_the_function(x, atan))
298             return t;
299         // tan(asin(x)) -> x*(1+x^2)^(-1/2)
300         if (is_ex_the_function(x, asin))
301             return t*power(exONE()-power(t,exTWO()),exMINUSHALF());
302         // tan(acos(x)) -> (1-x^2)^(1/2)/x
303         if (is_ex_the_function(x, acos))
304             return power(t,exMINUSONE())*power(exONE()-power(t,exTWO()),exHALF());
305     }
306     
307     // tan(float) -> float
308     if (x.info(info_flags::numeric) && !x.info(info_flags::crational)) {
309         return tan_evalf(x);
310     }
311     
312     return tan(x).hold();
313 }
314
315 static ex tan_diff(ex const & x, unsigned diff_param)
316 {
317     GINAC_ASSERT(diff_param==0);
318     
319     return (1+power(tan(x),exTWO()));
320 }
321
322 static ex tan_series(ex const & x, symbol const & s, ex const & point, int order)
323 {
324     // method:
325     // Taylor series where there is no pole falls back to tan_diff.
326     // On a pole simply expand sin(x)/cos(x).
327     ex xpoint = x.subs(s==point);
328     if (!(2*xpoint/Pi).info(info_flags::odd))
329         throw do_taylor();
330     // if we got here we have to care for a simple pole
331     return (sin(x)/cos(x)).series(s, point, order+2);
332 }
333
334 REGISTER_FUNCTION(tan, tan_eval, tan_evalf, tan_diff, tan_series);
335
336 //////////
337 // inverse sine (arc sine)
338 //////////
339
340 static ex asin_evalf(ex const & x)
341 {
342     BEGIN_TYPECHECK
343        TYPECHECK(x,numeric)
344     END_TYPECHECK(asin(x))
345     
346     return asin(ex_to_numeric(x)); // -> numeric asin(numeric)
347 }
348
349 static ex asin_eval(ex const & x)
350 {
351     if (x.info(info_flags::numeric)) {
352         // asin(0) -> 0
353         if (x.is_zero())
354             return x;
355         // asin(1/2) -> Pi/6
356         if (x.is_equal(exHALF()))
357             return numeric(1,6)*Pi;
358         // asin(1) -> Pi/2
359         if (x.is_equal(exONE()))
360             return numeric(1,2)*Pi;
361         // asin(-1/2) -> -Pi/6
362         if (x.is_equal(exMINUSHALF()))
363             return numeric(-1,6)*Pi;
364         // asin(-1) -> -Pi/2
365         if (x.is_equal(exMINUSONE()))
366             return numeric(-1,2)*Pi;
367         // asin(float) -> float
368         if (!x.info(info_flags::crational))
369             return asin_evalf(x);
370     }
371     
372     return asin(x).hold();
373 }
374
375 static ex asin_diff(ex const & x, unsigned diff_param)
376 {
377     GINAC_ASSERT(diff_param==0);
378     
379     return power(1-power(x,exTWO()),exMINUSHALF());
380 }
381
382 REGISTER_FUNCTION(asin, asin_eval, asin_evalf, asin_diff, NULL);
383
384 //////////
385 // inverse cosine (arc cosine)
386 //////////
387
388 static ex acos_evalf(ex const & x)
389 {
390     BEGIN_TYPECHECK
391        TYPECHECK(x,numeric)
392     END_TYPECHECK(acos(x))
393     
394     return acos(ex_to_numeric(x)); // -> numeric acos(numeric)
395 }
396
397 static ex acos_eval(ex const & x)
398 {
399     if (x.info(info_flags::numeric)) {
400         // acos(1) -> 0
401         if (x.is_equal(exONE()))
402             return exZERO();
403         // acos(1/2) -> Pi/3
404         if (x.is_equal(exHALF()))
405             return numeric(1,3)*Pi;
406         // acos(0) -> Pi/2
407         if (x.is_zero())
408             return numeric(1,2)*Pi;
409         // acos(-1/2) -> 2/3*Pi
410         if (x.is_equal(exMINUSHALF()))
411             return numeric(2,3)*Pi;
412         // acos(-1) -> Pi
413         if (x.is_equal(exMINUSONE()))
414             return Pi;
415         // acos(float) -> float
416         if (!x.info(info_flags::crational))
417             return acos_evalf(x);
418     }
419     
420     return acos(x).hold();
421 }
422
423 static ex acos_diff(ex const & x, unsigned diff_param)
424 {
425     GINAC_ASSERT(diff_param==0);
426     
427     return numMINUSONE()*power(1-power(x,exTWO()),exMINUSHALF());
428 }
429
430 REGISTER_FUNCTION(acos, acos_eval, acos_evalf, acos_diff, NULL);
431
432 //////////
433 // inverse tangent (arc tangent)
434 //////////
435
436 static ex atan_evalf(ex const & x)
437 {
438     BEGIN_TYPECHECK
439         TYPECHECK(x,numeric)
440     END_TYPECHECK(atan(x))
441     
442     return atan(ex_to_numeric(x)); // -> numeric atan(numeric)
443 }
444
445 static ex atan_eval(ex const & x)
446 {
447     if (x.info(info_flags::numeric)) {
448         // atan(0) -> 0
449         if (x.is_equal(exZERO()))
450             return exZERO();
451         // atan(float) -> float
452         if (!x.info(info_flags::crational))
453             return atan_evalf(x);
454     }
455     
456     return atan(x).hold();
457 }    
458
459 static ex atan_diff(ex const & x, unsigned diff_param)
460 {
461     GINAC_ASSERT(diff_param==0);
462
463     return power(1+x*x, -1);
464 }
465
466 REGISTER_FUNCTION(atan, atan_eval, atan_evalf, atan_diff, NULL);
467
468 //////////
469 // inverse tangent (atan2(y,x))
470 //////////
471
472 static ex atan2_evalf(ex const & y, ex const & x)
473 {
474     BEGIN_TYPECHECK
475         TYPECHECK(y,numeric)
476         TYPECHECK(x,numeric)
477     END_TYPECHECK(atan2(y,x))
478     
479     return atan(ex_to_numeric(y),ex_to_numeric(x)); // -> numeric atan(numeric)
480 }
481
482 static ex atan2_eval(ex const & y, ex const & x)
483 {
484     if (y.info(info_flags::numeric) && !y.info(info_flags::crational) &&
485         x.info(info_flags::numeric) && !x.info(info_flags::crational)) {
486         return atan2_evalf(y,x);
487     }
488     
489     return atan2(y,x).hold();
490 }    
491
492 static ex atan2_diff(ex const & y, ex const & x, unsigned diff_param)
493 {
494     GINAC_ASSERT(diff_param<2);
495     
496     if (diff_param==0) {
497         // d/dy atan(y,x)
498         return x*pow(pow(x,2)+pow(y,2),-1);
499     }
500     // d/dx atan(y,x)
501     return -y*pow(pow(x,2)+pow(y,2),-1);
502 }
503
504 REGISTER_FUNCTION(atan2, atan2_eval, atan2_evalf, atan2_diff, NULL);
505
506 //////////
507 // hyperbolic sine (trigonometric function)
508 //////////
509
510 static ex sinh_evalf(ex const & x)
511 {
512     BEGIN_TYPECHECK
513        TYPECHECK(x,numeric)
514     END_TYPECHECK(sinh(x))
515     
516     return sinh(ex_to_numeric(x)); // -> numeric sinh(numeric)
517 }
518
519 static ex sinh_eval(ex const & x)
520 {
521     if (x.info(info_flags::numeric)) {
522         // sinh(0) -> 0
523         if (x.is_zero())
524             return exZERO();
525         // sinh(float) -> float
526         if (!x.info(info_flags::crational))
527             return sinh_evalf(x);
528     }
529     
530     ex xOverPiI=x/Pi/I;
531     if (xOverPiI.info(info_flags::numeric)) {
532         // sinh(n*Pi*I) -> 0
533         if (xOverPiI.info(info_flags::integer))
534             return exZERO();
535         
536         // sinh((2n+1)*Pi*I/2) -> {+|-}I
537         ex xOverPiIMinusHalf=xOverPiI-exHALF();
538         if (xOverPiIMinusHalf.info(info_flags::even))
539             return I;
540         else if (xOverPiIMinusHalf.info(info_flags::odd))
541             return -I;
542     }
543     
544     if (is_ex_exactly_of_type(x, function)) {
545         ex t=x.op(0);
546         // sinh(asinh(x)) -> x
547         if (is_ex_the_function(x, asinh))
548             return t;
549         // sinh(acosh(x)) -> (x-1)^(1/2) * (x+1)^(1/2)
550         if (is_ex_the_function(x, acosh))
551             return power(t-exONE(),exHALF())*power(t+exONE(),exHALF());
552         // sinh(atanh(x)) -> x*(1-x^2)^(-1/2)
553         if (is_ex_the_function(x, atanh))
554             return t*power(exONE()-power(t,exTWO()),exMINUSHALF());
555     }
556     
557     return sinh(x).hold();
558 }
559
560 static ex sinh_diff(ex const & x, unsigned diff_param)
561 {
562     GINAC_ASSERT(diff_param==0);
563     
564     return cosh(x);
565 }
566
567 REGISTER_FUNCTION(sinh, sinh_eval, sinh_evalf, sinh_diff, NULL);
568
569 //////////
570 // hyperbolic cosine (trigonometric function)
571 //////////
572
573 static ex cosh_evalf(ex const & x)
574 {
575     BEGIN_TYPECHECK
576        TYPECHECK(x,numeric)
577     END_TYPECHECK(cosh(x))
578     
579     return cosh(ex_to_numeric(x)); // -> numeric cosh(numeric)
580 }
581
582 static ex cosh_eval(ex const & x)
583 {
584     if (x.info(info_flags::numeric)) {
585         // cosh(0) -> 1
586         if (x.is_zero())
587             return exONE();
588         // cosh(float) -> float
589         if (!x.info(info_flags::crational))
590             return cosh_evalf(x);
591     }
592     
593     ex xOverPiI=x/Pi/I;
594     if (xOverPiI.info(info_flags::numeric)) {
595         // cosh(n*Pi*I) -> {+|-}1
596         if (xOverPiI.info(info_flags::even))
597             return exONE();
598         else if (xOverPiI.info(info_flags::odd))
599             return exMINUSONE();
600         
601         // cosh((2n+1)*Pi*I/2) -> 0
602         ex xOverPiIMinusHalf=xOverPiI-exHALF();
603         if (xOverPiIMinusHalf.info(info_flags::integer))
604             return exZERO();
605     }
606     
607     if (is_ex_exactly_of_type(x, function)) {
608         ex t=x.op(0);
609         // cosh(acosh(x)) -> x
610         if (is_ex_the_function(x, acosh))
611             return t;
612         // cosh(asinh(x)) -> (1+x^2)^(1/2)
613         if (is_ex_the_function(x, asinh))
614             return power(exONE()+power(t,exTWO()),exHALF());
615         // cosh(atanh(x)) -> (1-x^2)^(-1/2)
616         if (is_ex_the_function(x, atanh))
617             return power(exONE()-power(t,exTWO()),exMINUSHALF());
618     }
619     
620     return cosh(x).hold();
621 }
622
623 static ex cosh_diff(ex const & x, unsigned diff_param)
624 {
625     GINAC_ASSERT(diff_param==0);
626     
627     return sinh(x);
628 }
629
630 REGISTER_FUNCTION(cosh, cosh_eval, cosh_evalf, cosh_diff, NULL);
631
632 //////////
633 // hyperbolic tangent (trigonometric function)
634 //////////
635
636 static ex tanh_evalf(ex const & x)
637 {
638     BEGIN_TYPECHECK
639        TYPECHECK(x,numeric)
640     END_TYPECHECK(tanh(x))
641     
642     return tanh(ex_to_numeric(x)); // -> numeric tanh(numeric)
643 }
644
645 static ex tanh_eval(ex const & x)
646 {
647     if (x.info(info_flags::numeric)) {
648         // tanh(0) -> 0
649         if (x.is_zero())
650             return exZERO();
651         // tanh(float) -> float
652         if (!x.info(info_flags::crational))
653             return tanh_evalf(x);
654     }
655     
656     if (is_ex_exactly_of_type(x, function)) {
657         ex t=x.op(0);
658         // tanh(atanh(x)) -> x
659         if (is_ex_the_function(x, atanh))
660             return t;
661         // tanh(asinh(x)) -> x*(1+x^2)^(-1/2)
662         if (is_ex_the_function(x, asinh))
663             return t*power(exONE()+power(t,exTWO()),exMINUSHALF());
664         // tanh(acosh(x)) -> (x-1)^(1/2)*(x+1)^(1/2)/x
665         if (is_ex_the_function(x, acosh))
666             return power(t-exONE(),exHALF())*power(t+exONE(),exHALF())*power(t,exMINUSONE());
667     }
668     
669     return tanh(x).hold();
670 }
671
672 static ex tanh_diff(ex const & x, unsigned diff_param)
673 {
674     GINAC_ASSERT(diff_param==0);
675     
676     return exONE()-power(tanh(x),exTWO());
677 }
678
679 static ex tanh_series(ex const & x, symbol const & s, ex const & point, int order)
680 {
681     // method:
682     // Taylor series where there is no pole falls back to tanh_diff.
683     // On a pole simply expand sinh(x)/cosh(x).
684     ex xpoint = x.subs(s==point);
685     if (!(2*I*xpoint/Pi).info(info_flags::odd))
686         throw do_taylor();
687     // if we got here we have to care for a simple pole
688     return (sinh(x)/cosh(x)).series(s, point, order+2);
689 }
690
691 REGISTER_FUNCTION(tanh, tanh_eval, tanh_evalf, tanh_diff, tanh_series);
692
693 //////////
694 // inverse hyperbolic sine (trigonometric function)
695 //////////
696
697 static ex asinh_evalf(ex const & x)
698 {
699     BEGIN_TYPECHECK
700        TYPECHECK(x,numeric)
701     END_TYPECHECK(asinh(x))
702     
703     return asinh(ex_to_numeric(x)); // -> numeric asinh(numeric)
704 }
705
706 static ex asinh_eval(ex const & x)
707 {
708     if (x.info(info_flags::numeric)) {
709         // asinh(0) -> 0
710         if (x.is_zero())
711             return exZERO();
712         // asinh(float) -> float
713         if (!x.info(info_flags::crational))
714             return asinh_evalf(x);
715     }
716     
717     return asinh(x).hold();
718 }
719
720 static ex asinh_diff(ex const & x, unsigned diff_param)
721 {
722     GINAC_ASSERT(diff_param==0);
723     
724     return power(1+power(x,exTWO()),exMINUSHALF());
725 }
726
727 REGISTER_FUNCTION(asinh, asinh_eval, asinh_evalf, asinh_diff, NULL);
728
729 //////////
730 // inverse hyperbolic cosine (trigonometric function)
731 //////////
732
733 static ex acosh_evalf(ex const & x)
734 {
735     BEGIN_TYPECHECK
736        TYPECHECK(x,numeric)
737     END_TYPECHECK(acosh(x))
738     
739     return acosh(ex_to_numeric(x)); // -> numeric acosh(numeric)
740 }
741
742 static ex acosh_eval(ex const & x)
743 {
744     if (x.info(info_flags::numeric)) {
745         // acosh(0) -> Pi*I/2
746         if (x.is_zero())
747             return Pi*I*numeric(1,2);
748         // acosh(1) -> 0
749         if (x.is_equal(exONE()))
750             return exZERO();
751         // acosh(-1) -> Pi*I
752         if (x.is_equal(exMINUSONE()))
753             return Pi*I;
754         // acosh(float) -> float
755         if (!x.info(info_flags::crational))
756             return acosh_evalf(x);
757     }
758     
759     return acosh(x).hold();
760 }
761
762 static ex acosh_diff(ex const & x, unsigned diff_param)
763 {
764     GINAC_ASSERT(diff_param==0);
765     
766     return power(x-1,exMINUSHALF())*power(x+1,exMINUSHALF());
767 }
768
769 REGISTER_FUNCTION(acosh, acosh_eval, acosh_evalf, acosh_diff, NULL);
770
771 //////////
772 // inverse hyperbolic tangent (trigonometric function)
773 //////////
774
775 static ex atanh_evalf(ex const & x)
776 {
777     BEGIN_TYPECHECK
778        TYPECHECK(x,numeric)
779     END_TYPECHECK(atanh(x))
780     
781     return atanh(ex_to_numeric(x)); // -> numeric atanh(numeric)
782 }
783
784 static ex atanh_eval(ex const & x)
785 {
786     if (x.info(info_flags::numeric)) {
787         // atanh(0) -> 0
788         if (x.is_zero())
789             return exZERO();
790         // atanh({+|-}1) -> throw
791         if (x.is_equal(exONE()) || x.is_equal(exONE()))
792             throw (std::domain_error("atanh_eval(): infinity"));
793         // atanh(float) -> float
794         if (!x.info(info_flags::crational))
795             return atanh_evalf(x);
796     }
797     
798     return atanh(x).hold();
799 }
800
801 static ex atanh_diff(ex const & x, unsigned diff_param)
802 {
803     GINAC_ASSERT(diff_param==0);
804     
805     return power(exONE()-power(x,exTWO()),exMINUSONE());
806 }
807
808 REGISTER_FUNCTION(atanh, atanh_eval, atanh_evalf, atanh_diff, NULL);
809
810 #ifndef NO_GINAC_NAMESPACE
811 } // namespace GiNaC
812 #endif // ndef NO_GINAC_NAMESPACE