My Project  0.0.16
QUCS Mapping
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
differentiate.cpp
Go to the documentation of this file.
1 /*
2  * differentiate.cpp - the Qucs equation differentiator implementations
3  *
4  * Copyright (C) 2007, 2008 Stefan Jahn <stefan@lkcc.org>
5  *
6  * This is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 2, or (at your option)
9  * any later version.
10  *
11  * This software 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
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with this package; see the file COPYING. If not, write to
18  * the Free Software Foundation, Inc., 51 Franklin Street - Fifth Floor,
19  * Boston, MA 02110-1301, USA.
20  *
21  * $Id: differentiate.cpp 1825 2011-03-11 20:42:14Z ela $
22  *
23  */
24 
25 #if HAVE_CONFIG_H
26 # include <config.h>
27 #endif
28 
29 #include <stdio.h>
30 #include <stdlib.h>
31 #include <string.h>
32 #include <ctype.h>
33 #include <math.h>
34 
35 #include "logging.h"
36 #include "complex.h"
37 #include "object.h"
38 #include "consts.h"
39 #include "equation.h"
40 #include "differentiate.h"
41 
42 using namespace eqn;
43 
44 // Short helper macros.
45 #define C(con) ((constant *) (con))
46 #define A(con) ((application *) (con))
47 #define R(con) ((reference *) (con))
48 #define D(con) (C(con)->d)
49 
50 #define isConst(n) ((n)->getTag()==CONSTANT && C(n)->getType()==TAG_DOUBLE)
51 #define isRef(r,v) ((r)->getTag()==REFERENCE && !strcmp(R(r)->n,v))
52 #define isZero(n) (isConst(n) && D(n) == 0.0)
53 #define isOne(n) (isConst(n) && D(n) == 1.0)
54 #define isNeg(n) (isConst(n) && D(n) == -1.0)
55 #define isEuler(n) ((isConst(n) && D(n) == M_E) || isRef(n,"e"))
56 #define isval(n,v) (isConst(n) && D(n) == v)
57 
58 #define isVar(v) ((v)->getTag()==REFERENCE)
59 #define isApp(v) ((v)->getTag()==APPLICATION)
60 #define isMul(v) (isApp(v) && !strcmp(A(v)->n,"*"))
61 #define isSqr(v) (isApp(v) && !strcmp(A(v)->n,"sqr"))
62 
63 #define retCon(val) \
64  constant * res = new constant (TAG_DOUBLE); res->d = val; return res;
65 #define defCon(def,val) \
66  constant * def = new constant (TAG_DOUBLE); def->d = val;
67 #define defRef(def,var) \
68  reference * def = new reference (); def->n = strdup (var);
69 #define retApp1(op,f0) \
70  application * res = new application (); res->n = strdup (op); \
71  res->nargs = 1; res->args = f0; res->args->setNext (NULL); return res;
72 #define defApp1(def,op,f0) \
73  application * def = new application (); def->n = strdup (op); \
74  def->nargs = 1; def->args = f0; def->args->setNext (NULL);
75 #define defApp2(def,op,f0,f1) \
76  application * def = new application (); def->n = strdup (op); \
77  def->nargs = 2; def->args = f0; def->args->append (f1);
78 #define retApp2(op,f0,f1) \
79  application * res = new application (); res->n = strdup (op); \
80  res->nargs = 2; res->args = f0; res->args->append (f1); return res;
81 #define retApp3(op,f0,f1,f2) \
82  application * res = new application (); res->n = strdup (op); \
83  res->nargs = 3; res->args = f0; res->args->append (f1); \
84  res->args->append (f2); return res;
85 #define defApp3(def,op,f0,f1,f2) \
86  application * def = new application (); def->n = strdup (op); \
87  def->nargs = 3; def->args = f0; def->args->append (f1); \
88  def->args->append (f2);
89 
90 #define _A(idx) app->args->get(idx)
91 #define _A0 _A(0)
92 #define _A1 _A(1)
93 #define _A2 _A(2)
94 #define _D0 _A(0)->differentiate (derivative)
95 #define _D1 _A(1)->differentiate (derivative)
96 #define _D2 _A(2)->differentiate (derivative)
97 
98 #define _AF0(var) node * var = _A0;
99 #define _AF1(var) node * var = _A1;
100 #define _AF2(var) node * var = _A2;
101 #define _AD0(var) node * var = _D0;
102 #define _AD1(var) node * var = _D1;
103 #define _AD2(var) node * var = _D2;
104 
105 #define _AA(a,idx) A(a)->args->get(idx)
106 #define _AA0(a) _AA(a,0)
107 #define _AA1(a) _AA(a,1)
108 
109 #define _AAF0(a,var) node * var = _AA0(a);
110 #define _AAF1(a,var) node * var = _AA1(a);
111 
112 node * differentiate::plus_binary (application * app, char * derivative) {
113  _AD0 (d0);
114  _AD1 (d1);
115  if (!isConst (d0) && !isConst (d1)) {
116  retApp2 ("+", d0, d1);
117  }
118  return plus_reduce (d0, d1);
119 }
120 
121 node * differentiate::plus_unary (application * app, char * derivative) {
122  _AD0 (d0);
123  return d0;
124 }
125 
126 node * differentiate::plus_reduce (node * f0, node * f1) {
127  if (isZero (f0) && isZero (f1)) {
128  delete f0; delete f1;
129  retCon (0);
130  } else if (isZero (f0)) {
131  delete f0;
132  return f1;
133  } else if (isZero (f1)) {
134  delete f1;
135  return f0;
136  } else if (isConst (f0) && isConst (f1)) {
137  nr_double_t t = D(f0) + D(f1);
138  delete f0; delete f1;
139  retCon (t);
140  } else {
141  retApp2 ("+", f0, f1);
142  }
143 }
144 
145 node * differentiate::minus_binary (application * app, char * derivative) {
146  _AD0 (d0);
147  _AD1 (d1);
148  if (!isConst (d0) && !isConst (d1)) {
149  retApp2 ("-", d0, d1);
150  }
151  return minus_reduce (d0, d1);
152 }
153 
154 node * differentiate::minus_unary (application * app, char * derivative) {
155  _AD0 (d0);
156  return minus_reduce (d0);
157 }
158 
159 node * differentiate::minus_reduce (node * f0) {
160  if (isZero (f0)) {
161  delete f0;
162  retCon (0);
163  } else if (isConst (f0)) {
164  nr_double_t t = -D(f0);
165  delete f0;
166  retCon (t);
167  }
168  retApp1 ("-", f0);
169 }
170 
171 node * differentiate::minus_reduce (node * f0, node * f1) {
172  if (isZero (f0) && isZero (f1)) {
173  delete f0; delete f1;
174  retCon (0);
175  } else if (isZero (f0)) {
176  delete f0;
177  return minus_reduce (f1);
178  } else if (isZero (f1)) {
179  delete f1;
180  return f0;
181  } else if (isConst (f0) && isConst (f1)) {
182  nr_double_t t = D(f0) - D(f1);
183  delete f0; delete f1;
184  retCon (t);
185  } else {
186  retApp2 ("-", f0, f1);
187  }
188 }
189 
190 node * differentiate::times (application * app, char * derivative) {
191  _AF0 (f0);
192  _AF1 (f1);
193  if (isConst (f0) && isConst (f1)) {
194  retCon (0);
195  }
196  _AD0 (d0);
197  _AD1 (d1);
198  node * t1 = times_reduce (f0->recreate(), d1);
199  node * t2 = times_reduce (f1->recreate(), d0);
200  return plus_reduce (t1, t2);
201 }
202 
203 node * differentiate::times_reduce (node * f0, node * f1) {
204  if (isZero (f0) || isZero (f1)) {
205  delete f0; delete f1;
206  retCon (0);
207  } else if (isOne (f0)) {
208  delete f0;
209  return f1;
210  } else if (isNeg (f0)) {
211  delete f0;
212  return minus_reduce (f1);
213  } else if (isOne (f1)) {
214  delete f1;
215  return f0;
216  } else if (isNeg (f1)) {
217  delete f1;
218  return minus_reduce (f0);
219  } else if (isConst (f0) && isConst (f1)) {
220  nr_double_t t = D(f0) * D(f1);
221  delete f0; delete f1;
222  retCon (t);
223  } else {
224  retApp2 ("*", f0, f1);
225  }
226 }
227 
228 node * differentiate::over (application * app, char * derivative) {
229  _AF0 (f0);
230  _AF1 (f1);
231  if (isConst (f0) && isConst (f1)) {
232  retCon (0);
233  }
234  _AD0 (d0);
235  _AD1 (d1);
236  node * t1 = times_reduce (f0->recreate(), d1);
237  node * t2 = times_reduce (f1->recreate(), d0);
238  node * t3 = minus_reduce (t2, t1);
239  node * t4 = sqr_reduce (f1->recreate());
240  return over_reduce (t3, t4);
241 }
242 
243 node * differentiate::over_reduce (node * f0, node * f1) {
244  if (isOne (f0) && isOne (f1)) {
245  delete f0; delete f1;
246  retCon (1);
247  } else if (isZero (f0)) {
248  delete f0; delete f1;
249  retCon (0);
250  } else if (isConst (f0) && isConst (f1)) {
251  if (isZero (f1)) {
252  retApp2 ("/", f0, f1);
253  }
254  nr_double_t t = D(f0) / D(f1);
255  delete f0; delete f1;
256  retCon (t);
257  } else if (isOne (f1)) {
258  delete f1;
259  return f0;
260  } else if (isNeg (f1)) {
261  delete f1;
262  return minus_reduce (f0);
263  } else {
264  over_reduce_adv (f0, f1);
265  retApp2 ("/", f0, f1);
266  }
267 }
268 
269 void differentiate::over_reduce_adv (node * &f0, node * &f1) {
270  if (isVar (f0)) {
271  if (isSqr (f1)) {
272  _AAF0 (f1,g1);
273  if (isVar (g1)) {
274  if (!strcmp (R(f0)->n, R(g1)->n)) {
275  defCon (one, 1);
276  reference * var = new reference (*R(g1));
277  delete f0;
278  delete f1;
279  f0 = one;
280  f1 = var;
281  }
282  }
283  }
284  }
285 }
286 
287 node * differentiate::power (application * app, char * derivative) {
288  _AF0 (f0);
289  _AF1 (f1);
290  if (isConst (f0) && isConst (f1)) {
291  retCon (0);
292  }
293  _AD0 (d0);
294  _AD1 (d1);
295  if (isZero (d1)) {
296  defCon (one, 1);
297  node * t1 = minus_reduce (f1->recreate(), one);
298  node * t2 = power_reduce (f0->recreate(), t1);
299  node * t3 = times_reduce (f1->recreate(), t2);
300  return times_reduce (t3, d0);
301  }
302  else {
303  node * t1 = power_reduce (f0->recreate(), f1->recreate());
304  node * ln = ln_reduce (f0->recreate());
305  node * t2 = times_reduce (d1, ln);
306  node * t3 = times_reduce (f1->recreate(), d0);
307  node * t4 = over_reduce (t3, f0->recreate());
308  node * t5 = plus_reduce (t2, t4);
309  return times_reduce (t1, t5);
310  }
311 }
312 
313 node * differentiate::power_reduce (node * f0, node * f1) {
314  if (isOne (f0)) {
315  delete f0; delete f1;
316  retCon (1);
317  } else if (isZero (f0)) {
318  delete f0; delete f1;
319  retCon (0);
320  } else if (isConst (f0) && isConst (f1)) {
321  if (isZero (f1)) {
322  delete f0; delete f1;
323  retCon (1);
324  }
325  nr_double_t t = pow (D(f0), D(f1));
326  delete f0; delete f1;
327  retCon (t);
328  } else if (isOne (f1)) {
329  delete f1;
330  return f0;
331  } else {
332  retApp2 ("^", f0, f1);
333  }
334 }
335 
336 node * differentiate::ln (application * app, char * derivative) {
337  _AF0 (f0);
338  _AD0 (d0);
339  return over_reduce (d0, f0->recreate ());
340 }
341 
342 node * differentiate::ln_reduce (node * f0) {
343  if (isOne (f0)) {
344  delete f0;
345  retCon (0);
346  } else if (isEuler (f0)) {
347  delete f0;
348  retCon (1);
349  }
350  retApp1 ("ln", f0);
351 }
352 
353 node * differentiate::log10 (application * app, char * derivative) {
354  _AF0 (f0);
355  _AD0 (d0);
356  node * t1 = over_reduce (d0, f0->recreate ());
357  defCon (ten, 10);
358  return over_reduce (t1, ln_reduce (ten));
359 }
360 
361 node * differentiate::log2 (application * app, char * derivative) {
362  _AF0 (f0);
363  _AD0 (d0);
364  node * t1 = over_reduce (d0, f0->recreate ());
365  defCon (two, 2);
366  return over_reduce (t1, ln_reduce (two));
367 }
368 
369 node * differentiate::sqrt (application * app, char * derivative) {
370  _AF0 (f0);
371  _AD0 (d0);
372  defCon (half, 0.5);
373  node * t1 = times_reduce (half, d0);
374  node * t2 = sqrt_reduce (f0->recreate());
375  return over_reduce (t1, t2);
376 }
377 
378 node * differentiate::sqrt_reduce (node * f0) {
379  if (isOne (f0)) {
380  delete f0;
381  retCon (1);
382  } else if (isZero (f0)) {
383  delete f0;
384  retCon (0);
385  }
386  retApp1 ("sqrt", f0);
387 }
388 
389 node * differentiate::app_reduce (const char * func, node * d0, node * f0) {
390  if (isOne (d0)) {
391  delete d0;
392  retApp1 (func, f0);
393  } else if (isZero (d0)) {
394  delete d0; delete f0;
395  retCon (0);
396  }
397  defApp1 (app, func, f0);
398  return times_reduce (d0, app);
399 }
400 
401 node * differentiate::exp (application * app, char * derivative) {
402  _AF0 (f0);
403  _AD0 (d0);
404  return app_reduce ("exp", d0, f0->recreate());
405 }
406 
407 node * differentiate::limexp (application * app, char * derivative) {
408  _AF0 (f0);
409  _AD0 (d0);
410  defCon (lexp, ::exp (M_LIMEXP));
411  defCon (lcon, M_LIMEXP);
412  defApp2 (ask, "<", f0->recreate(), lcon);
413  defApp1 (exp, "exp", f0->recreate());
414  defApp3 (ite, "?:", ask, exp, lexp);
415  return times_reduce (d0, ite);
416 }
417 
418 node * differentiate::sin (application * app, char * derivative) {
419  _AF0 (f0);
420  _AD0 (d0);
421  return app_reduce ("cos", d0, f0->recreate());
422 }
423 
424 node * differentiate::cos (application * app, char * derivative) {
425  _AF0 (f0);
426  _AD0 (d0);
427  node * t1 = minus_reduce (d0);
428  return app_reduce ("sin", t1, f0->recreate());
429 }
430 
431 node * differentiate::tan (application * app, char * derivative) {
432  _AF0 (f0);
433  _AD0 (d0);
434  defApp1 (sec, "sec", f0->recreate());
435  defCon (two, 2);
436  node * t1 = power_reduce (sec, two);
437  return times_reduce (d0, t1);
438 }
439 
440 node * differentiate::sec (application * app, char * derivative) {
441  _AF0 (f0);
442  _AD0 (d0);
443  defApp1 (sin, "sin", f0->recreate());
444  defApp1 (cos, "cos", f0->recreate());
445  defCon (two, 2);
446  node * t1 = power_reduce (cos, two);
447  node * t2 = over_reduce (sin, t1);
448  return times_reduce (d0, t2);
449 }
450 
451 node * differentiate::cot (application * app, char * derivative) {
452  _AF0 (f0);
453  _AD0 (d0);
454  defApp1 (cosec, "cosec", f0->recreate());
455  defCon (two, 2);
456  node * t1 = minus_reduce (d0);
457  node * t2 = power_reduce (cosec, two);
458  return times_reduce (t1, t2);
459 }
460 
461 node * differentiate::cosec (application * app, char * derivative) {
462  _AF0 (f0);
463  _AD0 (d0);
464  defApp1 (sin, "sin", f0->recreate());
465  defApp1 (cos, "cos", f0->recreate());
466  defCon (two, 2);
467  node * t1 = minus_reduce (d0);
468  node * t2 = power_reduce (sin, two);
469  node * t3 = over_reduce (cos, t2);
470  return times_reduce (t1, t3);
471 }
472 
473 node * differentiate::abs (application * app, char * derivative) {
474  _AF0 (f0);
475  _AD0 (d0);
476  return app_reduce ("sign", d0, f0->recreate());
477 }
478 
480  retCon (0);
481 }
482 
484  retCon (0);
485 }
486 
487 node * differentiate::arcsin (application * app, char * derivative) {
488  _AF0 (f0);
489  _AD0 (d0);
490  node * sqr = sqr_reduce (f0->recreate());
491  defCon (one, 1);
492  node * t1 = minus_reduce (one, sqr);
493  node * t2 = sqrt_reduce (t1);
494  return over_reduce (d0, t2);
495 }
496 
497 node * differentiate::square (application * app, char * derivative) {
498  _AF0 (f0);
499  _AD0 (d0);
500  defCon (two, 2);
501  node * t1 = times_reduce (two, d0);
502  return times_reduce (t1, f0->recreate());
503 }
504 
505 node * differentiate::sqr_reduce (node * f0) {
506  if (isOne (f0)) {
507  delete f0;
508  retCon (1);
509  } else if (isZero (f0)) {
510  delete f0;
511  retCon (0);
512  } else if (isConst (f0)) {
513  nr_double_t t = D(f0) * D(f0);
514  delete f0;
515  retCon (t);
516  } else {
517  retApp1 ("sqr", f0);
518  }
519 }
520 
521 node * differentiate::arccos (application * app, char * derivative) {
522  _AF0 (f0);
523  _AD0 (d0);
524  node * sqr = sqr_reduce (f0->recreate());
525  defCon (one, 1);
526  node * t1 = minus_reduce (one, sqr);
527  node * t2 = sqrt_reduce (t1);
528  node * t3 = minus_reduce (d0);
529  return over_reduce (t3, t2);
530 }
531 
532 node * differentiate::arctan (application * app, char * derivative) {
533  _AF0 (f0);
534  _AD0 (d0);
535  node * sqr = sqr_reduce (f0->recreate());
536  defCon (one, 1);
537  node * t1 = plus_reduce (one, sqr);
538  return over_reduce (d0, t1);
539 }
540 
541 node * differentiate::arccot (application * app, char * derivative) {
542  _AF0 (f0);
543  _AD0 (d0);
544  node * sqr = sqr_reduce (f0->recreate());
545  defCon (one, 1);
546  node * t1 = plus_reduce (one, sqr);
547  node * t2 = minus_reduce (d0);
548  return over_reduce (t2, t1);
549 }
550 
551 node * differentiate::arcsec (application * app, char * derivative) {
552  _AF0 (f0);
553  _AD0 (d0);
554  node * sqr = sqr_reduce (f0->recreate());
555  defCon (one, 1);
556  node * t1 = minus_reduce (sqr, one);
557  node * t2 = sqrt_reduce (t1);
558  node * t3 = times_reduce (f0->recreate(), t2);
559  return over_reduce (d0, t3);
560 }
561 
562 node * differentiate::arccosec (application * app, char * derivative) {
563  _AF0 (f0);
564  _AD0 (d0);
565  node * sqr = sqr_reduce (f0->recreate());
566  defCon (one, 1);
567  node * t1 = minus_reduce (sqr, one);
568  node * t2 = sqrt_reduce (t1);
569  node * t3 = times_reduce (f0->recreate(), t2);
570  node * t4 = minus_reduce (d0);
571  return over_reduce (t4, t3);
572 }
573 
574 node * differentiate::sinh (application * app, char * derivative) {
575  _AF0 (f0);
576  _AD0 (d0);
577  return app_reduce ("cosh", d0, f0->recreate());
578 }
579 
580 node * differentiate::cosh (application * app, char * derivative) {
581  _AF0 (f0);
582  _AD0 (d0);
583  return app_reduce ("sinh", d0, f0->recreate());
584 }
585 
586 node * differentiate::tanh (application * app, char * derivative) {
587  _AF0 (f0);
588  _AD0 (d0);
589  defApp1 (cosh, "cosh", f0->recreate());
590  defCon (two, 2);
591  node * t1 = power_reduce (cosh, two);
592  return over_reduce (d0, t1);
593 }
594 
595 node * differentiate::coth (application * app, char * derivative) {
596  _AF0 (f0);
597  _AD0 (d0);
598  defApp1 (sinh, "sinh", f0->recreate());
599  defCon (two, 2);
600  node * t1 = power_reduce (sinh, two);
601  node * t2 = minus_reduce (d0);
602  return over_reduce (t2, t1);
603 }
604 
605 node * differentiate::artanh (application * app, char * derivative) {
606  _AF0 (f0);
607  _AD0 (d0);
608  node * sqr = sqr_reduce (f0->recreate());
609  defCon (one, 1);
610  node * t1 = minus_reduce (one, sqr);
611  return over_reduce (d0, t1);
612 }
613 
614 node * differentiate::arcoth (application * app, char * derivative) {
615  _AF0 (f0);
616  _AD0 (d0);
617  node * sqr = sqr_reduce (f0->recreate());
618  defCon (one, 1);
619  node * t1 = minus_reduce (sqr, one);
620  node * t2 = minus_reduce (d0);
621  return over_reduce (t2, t1);
622 }
623 
624 node * differentiate::arcosh (application * app, char * derivative) {
625  _AF0 (f0);
626  _AD0 (d0);
627  node * sqr = sqr_reduce (f0->recreate());
628  defCon (one, 1);
629  node * t1 = minus_reduce (sqr, one);
630  node * t2 = sqrt_reduce (t1);
631  return over_reduce (d0, t2);
632 }
633 
634 node * differentiate::arsinh (application * app, char * derivative) {
635  _AF0 (f0);
636  _AD0 (d0);
637  node * sqr = sqr_reduce (f0->recreate());
638  defCon (one, 1);
639  node * t1 = plus_reduce (sqr, one);
640  node * t2 = sqrt_reduce (t1);
641  return over_reduce (d0, t2);
642 }
643 
644 node * differentiate::arsech (application * app, char * derivative) {
645  _AF0 (f0);
646  _AD0 (d0);
647  node * sqr = sqr_reduce (f0->recreate());
648  defCon (one, 1);
649  node * t1 = minus_reduce (one, sqr);
650  node * t2 = sqrt_reduce (t1);
651  node * t3 = times_reduce (f0->recreate(), t2);
652  node * t4 = minus_reduce (d0);
653  return over_reduce (t4, t3);
654 }
655 
656 node * differentiate::arcosech (application * app, char * derivative) {
657  _AF0 (f0);
658  _AD0 (d0);
659  node * sqr = sqr_reduce (f0->recreate());
660  defCon (one, 1);
661  node * t1 = plus_reduce (one, sqr);
662  node * t2 = sqrt_reduce (t1);
663  node * t3 = times_reduce (f0->recreate(), t2);
664  node * t4 = minus_reduce (d0);
665  return over_reduce (t4, t3);
666 }
667 
668 node * differentiate::ifthenelse (application * app, char * derivative) {
669  _AF0 (f0);
670  _AD1 (d1);
671  _AD2 (d2);
672  if (isConst (d1) && isConst (d2)) {
673  if (D(d1) == D(d2)) {
674  nr_double_t t = D(d1);
675  delete d1; delete d2;
676  retCon (t);
677  }
678  }
679  retApp3 ("?:", f0->recreate(), d1, d2);
680 }
681 
682 node * differentiate::sinc (application * app, char * derivative) {
683  _AF0 (f0);
684  _AD0 (d0);
685  defApp1 (sinc, "sinc", f0->recreate());
686  defApp1 (cos, "cos", f0->recreate());
687  node * t1 = minus_reduce (cos, sinc);
688  node * t2 = over_reduce (t1, f0->recreate());
689  return times_reduce (d0, t2);
690 }
691 
692 node * differentiate::norm (application * app, char * derivative) {
693  _AF0 (f0);
694  _AD0 (d0);
695  defCon (two, 2);
696  node * t1 = times_reduce (d0, two);
697  return times_reduce (t1, f0->recreate());
698 }
699 
700 node * differentiate::xhypot (application * app, char * derivative) {
701  _AF0 (f0);
702  _AF1 (f1);
703  _AD0 (d0);
704  _AD1 (d1);
705  node * t1 = hypot_reduce (f0->recreate(), f1->recreate());
706  node * t2 = times_reduce (d0, f0->recreate());
707  node * t3 = times_reduce (d1, f1->recreate());
708  node * t4 = plus_reduce (t2, t3);
709  return over_reduce (t4, t1);
710 }
711 
712 node * differentiate::hypot_reduce (node * f0, node * f1) {
713  if (isZero (f0) && isZero (f1)) {
714  delete f0; delete f1;
715  retCon (0);
716  } else if (isZero (f0)) {
717  delete f0;
718  return sqrt_reduce (sqr_reduce (f1));
719  } else if (isZero (f1)) {
720  delete f1;
721  return sqrt_reduce (sqr_reduce (f0));
722  } else if (isConst (f0) && isConst (f1)) {
723  nr_double_t t = ::xhypot (D(f0), D(f1));
724  delete f0; delete f1;
725  retCon (t);
726  } else {
727  retApp2 ("hypot", f0, f1);
728  }
729 }
730 
731 #include "constants.h"
732 
733 node * differentiate::vt (application * app, char * derivative) {
734  _AD0 (d0);
735  defCon (con, kBoverQ);
736  return times_reduce (d0, con);
737 }
738 
739 // List of differentiators.
741  { "+", differentiate::plus_binary, 2 },
742  { "+", differentiate::plus_unary, 1 },
743  { "-", differentiate::minus_binary, 2 },
744  { "-", differentiate::minus_unary, 1 },
745  { "*", differentiate::times, 2 },
746  { "/", differentiate::over, 2 },
747  { "^", differentiate::power, 2 },
748 
749  { "?:", differentiate::ifthenelse, 3 },
750 
751  { "ln", differentiate::ln, 1 },
752  { "log10", differentiate::log10, 1 },
753  { "log2", differentiate::log2, 1 },
754  { "sqrt", differentiate::sqrt, 1 },
755  { "exp", differentiate::exp, 1 },
756  { "sinc", differentiate::sinc, 1 },
757  { "norm", differentiate::norm, 1 },
758  { "sin", differentiate::sin, 1 },
759  { "cos", differentiate::cos, 1 },
760  { "tan", differentiate::tan, 1 },
761  { "sec", differentiate::sec, 1 },
762  { "cot", differentiate::cot, 1 },
763  { "cosec", differentiate::cosec, 1 },
764  { "abs", differentiate::abs, 1 },
765  { "step", differentiate::step, 1 },
766  { "sign", differentiate::sign, 1 },
767  { "arcsin", differentiate::arcsin, 1 },
768  { "arccos", differentiate::arccos, 1 },
769  { "arctan", differentiate::arctan, 1 },
770  { "arccot", differentiate::arccot, 1 },
771  { "arcsec", differentiate::arcsec, 1 },
772  { "arccosec", differentiate::arccosec, 1 },
773  { "sqr", differentiate::square, 1 },
774  { "sinh", differentiate::sinh, 1 },
775  { "cosh", differentiate::cosh, 1 },
776  { "tanh", differentiate::tanh, 1 },
777  { "coth", differentiate::coth, 1 },
778  { "arsinh", differentiate::arsinh, 1 },
779  { "arcosh", differentiate::arcosh, 1 },
780  { "artanh", differentiate::artanh, 1 },
781  { "arcoth", differentiate::arcoth, 1 },
782  { "arsech", differentiate::arsech, 1 },
783  { "arcosech", differentiate::arcosech, 1 },
784  { "hypot", differentiate::xhypot, 2 },
785  { "limexp", differentiate::limexp, 1 },
786  { "vt", differentiate::vt, 1 },
787 
788  { NULL, NULL, 0 /* end of list */ }
789 };