Qucs-core
0.0.19
|
00001 /* 00002 * differentiate.cpp - the Qucs equation differentiator implementations 00003 * 00004 * Copyright (C) 2007, 2008 Stefan Jahn <stefan@lkcc.org> 00005 * 00006 * This is free software; you can redistribute it and/or modify 00007 * it under the terms of the GNU General Public License as published by 00008 * the Free Software Foundation; either version 2, or (at your option) 00009 * any later version. 00010 * 00011 * This software is distributed in the hope that it will be useful, 00012 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00013 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00014 * GNU General Public License for more details. 00015 * 00016 * You should have received a copy of the GNU General Public License 00017 * along with this package; see the file COPYING. If not, write to 00018 * the Free Software Foundation, Inc., 51 Franklin Street - Fifth Floor, 00019 * Boston, MA 02110-1301, USA. 00020 * 00021 * $Id$ 00022 * 00023 */ 00024 00025 #if HAVE_CONFIG_H 00026 # include <config.h> 00027 #endif 00028 00029 #include <stdio.h> 00030 #include <stdlib.h> 00031 #include <string.h> 00032 #include <ctype.h> 00033 #include <cmath> 00034 00035 #include "logging.h" 00036 #include "complex.h" 00037 #include "object.h" 00038 #include "consts.h" 00039 #include "equation.h" 00040 #include "differentiate.h" 00041 00042 using namespace qucs; 00043 using namespace qucs::eqn; 00044 00045 // Short helper macros. 00046 #define C(con) ((constant *) (con)) 00047 #define A(con) ((application *) (con)) 00048 #define R(con) ((reference *) (con)) 00049 #define D(con) (C(con)->d) 00050 00051 #define isConst(n) ((n)->getTag()==CONSTANT && C(n)->getType()==TAG_DOUBLE) 00052 #define isRef(r,v) ((r)->getTag()==REFERENCE && !strcmp(R(r)->n,v)) 00053 #define isZero(n) (isConst(n) && D(n) == 0.0) 00054 #define isOne(n) (isConst(n) && D(n) == 1.0) 00055 #define isNeg(n) (isConst(n) && D(n) == -1.0) 00056 #define isEuler(n) ((isConst(n) && D(n) == euler) || isRef(n,"e")) 00057 #define isval(n,v) (isConst(n) && D(n) == v) 00058 00059 #define isVar(v) ((v)->getTag()==REFERENCE) 00060 #define isApp(v) ((v)->getTag()==APPLICATION) 00061 #define isMul(v) (isApp(v) && !strcmp(A(v)->n,"*")) 00062 #define isSqr(v) (isApp(v) && !strcmp(A(v)->n,"sqr")) 00063 00064 #define retCon(val) \ 00065 constant * res = new constant (TAG_DOUBLE); res->d = val; return res; 00066 #define defCon(def,val) \ 00067 constant * def = new constant (TAG_DOUBLE); def->d = val; 00068 #define defRef(def,var) \ 00069 reference * def = new reference (); def->n = strdup (var); 00070 #define retApp1(op,f0) \ 00071 application * res = new application (); res->n = strdup (op); \ 00072 res->nargs = 1; res->args = f0; res->args->setNext (NULL); return res; 00073 #define defApp1(def,op,f0) \ 00074 application * def = new application (); def->n = strdup (op); \ 00075 def->nargs = 1; def->args = f0; def->args->setNext (NULL); 00076 #define defApp2(def,op,f0,f1) \ 00077 application * def = new application (); def->n = strdup (op); \ 00078 def->nargs = 2; def->args = f0; def->args->append (f1); 00079 #define retApp2(op,f0,f1) \ 00080 application * res = new application (); res->n = strdup (op); \ 00081 res->nargs = 2; res->args = f0; res->args->append (f1); return res; 00082 #define retApp3(op,f0,f1,f2) \ 00083 application * res = new application (); res->n = strdup (op); \ 00084 res->nargs = 3; res->args = f0; res->args->append (f1); \ 00085 res->args->append (f2); return res; 00086 #define defApp3(def,op,f0,f1,f2) \ 00087 application * def = new application (); def->n = strdup (op); \ 00088 def->nargs = 3; def->args = f0; def->args->append (f1); \ 00089 def->args->append (f2); 00090 00091 #define _A(idx) app->args->get(idx) 00092 #define _A0 _A(0) 00093 #define _A1 _A(1) 00094 #define _A2 _A(2) 00095 #define _D0 _A(0)->differentiate (derivative) 00096 #define _D1 _A(1)->differentiate (derivative) 00097 #define _D2 _A(2)->differentiate (derivative) 00098 00099 #define _AF0(var) node * var = _A0; 00100 #define _AF1(var) node * var = _A1; 00101 #define _AF2(var) node * var = _A2; 00102 #define _AD0(var) node * var = _D0; 00103 #define _AD1(var) node * var = _D1; 00104 #define _AD2(var) node * var = _D2; 00105 00106 #define _AA(a,idx) A(a)->args->get(idx) 00107 #define _AA0(a) _AA(a,0) 00108 #define _AA1(a) _AA(a,1) 00109 00110 #define _AAF0(a,var) node * var = _AA0(a); 00111 #define _AAF1(a,var) node * var = _AA1(a); 00112 00113 node * differentiate::plus_binary (application * app, char * derivative) { 00114 _AD0 (d0); 00115 _AD1 (d1); 00116 if (!isConst (d0) && !isConst (d1)) { 00117 retApp2 ("+", d0, d1); 00118 } 00119 return plus_reduce (d0, d1); 00120 } 00121 00122 node * differentiate::plus_unary (application * app, char * derivative) { 00123 _AD0 (d0); 00124 return d0; 00125 } 00126 00127 node * differentiate::plus_reduce (node * f0, node * f1) { 00128 if (isZero (f0) && isZero (f1)) { 00129 delete f0; delete f1; 00130 retCon (0); 00131 } else if (isZero (f0)) { 00132 delete f0; 00133 return f1; 00134 } else if (isZero (f1)) { 00135 delete f1; 00136 return f0; 00137 } else if (isConst (f0) && isConst (f1)) { 00138 nr_double_t t = D(f0) + D(f1); 00139 delete f0; delete f1; 00140 retCon (t); 00141 } else { 00142 retApp2 ("+", f0, f1); 00143 } 00144 } 00145 00146 node * differentiate::minus_binary (application * app, char * derivative) { 00147 _AD0 (d0); 00148 _AD1 (d1); 00149 if (!isConst (d0) && !isConst (d1)) { 00150 retApp2 ("-", d0, d1); 00151 } 00152 return minus_reduce (d0, d1); 00153 } 00154 00155 node * differentiate::minus_unary (application * app, char * derivative) { 00156 _AD0 (d0); 00157 return minus_reduce (d0); 00158 } 00159 00160 node * differentiate::minus_reduce (node * f0) { 00161 if (isZero (f0)) { 00162 delete f0; 00163 retCon (0); 00164 } else if (isConst (f0)) { 00165 nr_double_t t = -D(f0); 00166 delete f0; 00167 retCon (t); 00168 } 00169 retApp1 ("-", f0); 00170 } 00171 00172 node * differentiate::minus_reduce (node * f0, node * f1) { 00173 if (isZero (f0) && isZero (f1)) { 00174 delete f0; delete f1; 00175 retCon (0); 00176 } else if (isZero (f0)) { 00177 delete f0; 00178 return minus_reduce (f1); 00179 } else if (isZero (f1)) { 00180 delete f1; 00181 return f0; 00182 } else if (isConst (f0) && isConst (f1)) { 00183 nr_double_t t = D(f0) - D(f1); 00184 delete f0; delete f1; 00185 retCon (t); 00186 } else { 00187 retApp2 ("-", f0, f1); 00188 } 00189 } 00190 00191 node * differentiate::times (application * app, char * derivative) { 00192 _AF0 (f0); 00193 _AF1 (f1); 00194 if (isConst (f0) && isConst (f1)) { 00195 retCon (0); 00196 } 00197 _AD0 (d0); 00198 _AD1 (d1); 00199 node * t1 = times_reduce (f0->recreate(), d1); 00200 node * t2 = times_reduce (f1->recreate(), d0); 00201 return plus_reduce (t1, t2); 00202 } 00203 00204 node * differentiate::times_reduce (node * f0, node * f1) { 00205 if (isZero (f0) || isZero (f1)) { 00206 delete f0; delete f1; 00207 retCon (0); 00208 } else if (isOne (f0)) { 00209 delete f0; 00210 return f1; 00211 } else if (isNeg (f0)) { 00212 delete f0; 00213 return minus_reduce (f1); 00214 } else if (isOne (f1)) { 00215 delete f1; 00216 return f0; 00217 } else if (isNeg (f1)) { 00218 delete f1; 00219 return minus_reduce (f0); 00220 } else if (isConst (f0) && isConst (f1)) { 00221 nr_double_t t = D(f0) * D(f1); 00222 delete f0; delete f1; 00223 retCon (t); 00224 } else { 00225 retApp2 ("*", f0, f1); 00226 } 00227 } 00228 00229 node * differentiate::over (application * app, char * derivative) { 00230 _AF0 (f0); 00231 _AF1 (f1); 00232 if (isConst (f0) && isConst (f1)) { 00233 retCon (0); 00234 } 00235 _AD0 (d0); 00236 _AD1 (d1); 00237 node * t1 = times_reduce (f0->recreate(), d1); 00238 node * t2 = times_reduce (f1->recreate(), d0); 00239 node * t3 = minus_reduce (t2, t1); 00240 node * t4 = sqr_reduce (f1->recreate()); 00241 return over_reduce (t3, t4); 00242 } 00243 00244 node * differentiate::over_reduce (node * f0, node * f1) { 00245 if (isOne (f0) && isOne (f1)) { 00246 delete f0; delete f1; 00247 retCon (1); 00248 } else if (isZero (f0)) { 00249 delete f0; delete f1; 00250 retCon (0); 00251 } else if (isConst (f0) && isConst (f1)) { 00252 if (isZero (f1)) { 00253 retApp2 ("/", f0, f1); 00254 } 00255 nr_double_t t = D(f0) / D(f1); 00256 delete f0; delete f1; 00257 retCon (t); 00258 } else if (isOne (f1)) { 00259 delete f1; 00260 return f0; 00261 } else if (isNeg (f1)) { 00262 delete f1; 00263 return minus_reduce (f0); 00264 } else { 00265 over_reduce_adv (f0, f1); 00266 retApp2 ("/", f0, f1); 00267 } 00268 } 00269 00270 void differentiate::over_reduce_adv (node * &f0, node * &f1) { 00271 if (isVar (f0)) { 00272 if (isSqr (f1)) { 00273 _AAF0 (f1,g1); 00274 if (isVar (g1)) { 00275 if (!strcmp (R(f0)->n, R(g1)->n)) { 00276 defCon (one, 1); 00277 reference * var = new reference (*R(g1)); 00278 delete f0; 00279 delete f1; 00280 f0 = one; 00281 f1 = var; 00282 } 00283 } 00284 } 00285 } 00286 } 00287 00288 node * differentiate::power (application * app, char * derivative) { 00289 _AF0 (f0); 00290 _AF1 (f1); 00291 if (isConst (f0) && isConst (f1)) { 00292 retCon (0); 00293 } 00294 _AD0 (d0); 00295 _AD1 (d1); 00296 if (isZero (d1)) { 00297 defCon (one, 1); 00298 node * t1 = minus_reduce (f1->recreate(), one); 00299 node * t2 = power_reduce (f0->recreate(), t1); 00300 node * t3 = times_reduce (f1->recreate(), t2); 00301 return times_reduce (t3, d0); 00302 } 00303 else { 00304 node * t1 = power_reduce (f0->recreate(), f1->recreate()); 00305 node * ln = ln_reduce (f0->recreate()); 00306 node * t2 = times_reduce (d1, ln); 00307 node * t3 = times_reduce (f1->recreate(), d0); 00308 node * t4 = over_reduce (t3, f0->recreate()); 00309 node * t5 = plus_reduce (t2, t4); 00310 return times_reduce (t1, t5); 00311 } 00312 } 00313 00314 node * differentiate::power_reduce (node * f0, node * f1) { 00315 if (isOne (f0)) { 00316 delete f0; delete f1; 00317 retCon (1); 00318 } else if (isZero (f0)) { 00319 delete f0; delete f1; 00320 retCon (0); 00321 } else if (isConst (f0) && isConst (f1)) { 00322 if (isZero (f1)) { 00323 delete f0; delete f1; 00324 retCon (1); 00325 } 00326 nr_double_t t = std::pow (D(f0), D(f1)); 00327 delete f0; delete f1; 00328 retCon (t); 00329 } else if (isOne (f1)) { 00330 delete f1; 00331 return f0; 00332 } else { 00333 retApp2 ("^", f0, f1); 00334 } 00335 } 00336 00337 node * differentiate::ln (application * app, char * derivative) { 00338 _AF0 (f0); 00339 _AD0 (d0); 00340 return over_reduce (d0, f0->recreate ()); 00341 } 00342 00343 node * differentiate::ln_reduce (node * f0) { 00344 if (isOne (f0)) { 00345 delete f0; 00346 retCon (0); 00347 } else if (isEuler (f0)) { 00348 delete f0; 00349 retCon (1); 00350 } 00351 retApp1 ("ln", f0); 00352 } 00353 00354 node * differentiate::log10 (application * app, char * derivative) { 00355 _AF0 (f0); 00356 _AD0 (d0); 00357 node * t1 = over_reduce (d0, f0->recreate ()); 00358 defCon (ten, 10); 00359 return over_reduce (t1, ln_reduce (ten)); 00360 } 00361 00362 node * differentiate::log2 (application * app, char * derivative) { 00363 _AF0 (f0); 00364 _AD0 (d0); 00365 node * t1 = over_reduce (d0, f0->recreate ()); 00366 defCon (two, 2); 00367 return over_reduce (t1, ln_reduce (two)); 00368 } 00369 00370 node * differentiate::sqrt (application * app, char * derivative) { 00371 _AF0 (f0); 00372 _AD0 (d0); 00373 defCon (half, 0.5); 00374 node * t1 = times_reduce (half, d0); 00375 node * t2 = sqrt_reduce (f0->recreate()); 00376 return over_reduce (t1, t2); 00377 } 00378 00379 node * differentiate::sqrt_reduce (node * f0) { 00380 if (isOne (f0)) { 00381 delete f0; 00382 retCon (1); 00383 } else if (isZero (f0)) { 00384 delete f0; 00385 retCon (0); 00386 } 00387 retApp1 ("sqrt", f0); 00388 } 00389 00390 node * differentiate::app_reduce (const char * func, node * d0, node * f0) { 00391 if (isOne (d0)) { 00392 delete d0; 00393 retApp1 (func, f0); 00394 } else if (isZero (d0)) { 00395 delete d0; delete f0; 00396 retCon (0); 00397 } 00398 defApp1 (app, func, f0); 00399 return times_reduce (d0, app); 00400 } 00401 00402 node * differentiate::exp (application * app, char * derivative) { 00403 _AF0 (f0); 00404 _AD0 (d0); 00405 return app_reduce ("exp", d0, f0->recreate()); 00406 } 00407 00408 node * differentiate::limexp (application * app, char * derivative) { 00409 _AF0 (f0); 00410 _AD0 (d0); 00411 defCon (lexp, ::exp (limitexp)); 00412 defCon (lcon, limitexp); 00413 defApp2 (ask, "<", f0->recreate(), lcon); 00414 defApp1 (exp, "exp", f0->recreate()); 00415 defApp3 (ite, "?:", ask, exp, lexp); 00416 return times_reduce (d0, ite); 00417 } 00418 00419 node * differentiate::sin (application * app, char * derivative) { 00420 _AF0 (f0); 00421 _AD0 (d0); 00422 return app_reduce ("cos", d0, f0->recreate()); 00423 } 00424 00425 node * differentiate::cos (application * app, char * derivative) { 00426 _AF0 (f0); 00427 _AD0 (d0); 00428 node * t1 = minus_reduce (d0); 00429 return app_reduce ("sin", t1, f0->recreate()); 00430 } 00431 00432 node * differentiate::tan (application * app, char * derivative) { 00433 _AF0 (f0); 00434 _AD0 (d0); 00435 defApp1 (sec, "sec", f0->recreate()); 00436 defCon (two, 2); 00437 node * t1 = power_reduce (sec, two); 00438 return times_reduce (d0, t1); 00439 } 00440 00441 node * differentiate::sec (application * app, char * derivative) { 00442 _AF0 (f0); 00443 _AD0 (d0); 00444 defApp1 (sin, "sin", f0->recreate()); 00445 defApp1 (cos, "cos", f0->recreate()); 00446 defCon (two, 2); 00447 node * t1 = power_reduce (cos, two); 00448 node * t2 = over_reduce (sin, t1); 00449 return times_reduce (d0, t2); 00450 } 00451 00452 node * differentiate::cot (application * app, char * derivative) { 00453 _AF0 (f0); 00454 _AD0 (d0); 00455 defApp1 (cosec, "cosec", f0->recreate()); 00456 defCon (two, 2); 00457 node * t1 = minus_reduce (d0); 00458 node * t2 = power_reduce (cosec, two); 00459 return times_reduce (t1, t2); 00460 } 00461 00462 node * differentiate::cosec (application * app, char * derivative) { 00463 _AF0 (f0); 00464 _AD0 (d0); 00465 defApp1 (sin, "sin", f0->recreate()); 00466 defApp1 (cos, "cos", f0->recreate()); 00467 defCon (two, 2); 00468 node * t1 = minus_reduce (d0); 00469 node * t2 = power_reduce (sin, two); 00470 node * t3 = over_reduce (cos, t2); 00471 return times_reduce (t1, t3); 00472 } 00473 00474 node * differentiate::abs (application * app, char * derivative) { 00475 _AF0 (f0); 00476 _AD0 (d0); 00477 return app_reduce ("sign", d0, f0->recreate()); 00478 } 00479 00480 node * differentiate::step (application *, char *) { 00481 retCon (0); 00482 } 00483 00484 node * differentiate::sign (application *, char *) { 00485 retCon (0); 00486 } 00487 00488 node * differentiate::arcsin (application * app, char * derivative) { 00489 _AF0 (f0); 00490 _AD0 (d0); 00491 node * sqr = sqr_reduce (f0->recreate()); 00492 defCon (one, 1); 00493 node * t1 = minus_reduce (one, sqr); 00494 node * t2 = sqrt_reduce (t1); 00495 return over_reduce (d0, t2); 00496 } 00497 00498 node * differentiate::square (application * app, char * derivative) { 00499 _AF0 (f0); 00500 _AD0 (d0); 00501 defCon (two, 2); 00502 node * t1 = times_reduce (two, d0); 00503 return times_reduce (t1, f0->recreate()); 00504 } 00505 00506 node * differentiate::sqr_reduce (node * f0) { 00507 if (isOne (f0)) { 00508 delete f0; 00509 retCon (1); 00510 } else if (isZero (f0)) { 00511 delete f0; 00512 retCon (0); 00513 } else if (isConst (f0)) { 00514 nr_double_t t = D(f0) * D(f0); 00515 delete f0; 00516 retCon (t); 00517 } else { 00518 retApp1 ("sqr", f0); 00519 } 00520 } 00521 00522 node * differentiate::arccos (application * app, char * derivative) { 00523 _AF0 (f0); 00524 _AD0 (d0); 00525 node * sqr = sqr_reduce (f0->recreate()); 00526 defCon (one, 1); 00527 node * t1 = minus_reduce (one, sqr); 00528 node * t2 = sqrt_reduce (t1); 00529 node * t3 = minus_reduce (d0); 00530 return over_reduce (t3, t2); 00531 } 00532 00533 node * differentiate::arctan (application * app, char * derivative) { 00534 _AF0 (f0); 00535 _AD0 (d0); 00536 node * sqr = sqr_reduce (f0->recreate()); 00537 defCon (one, 1); 00538 node * t1 = plus_reduce (one, sqr); 00539 return over_reduce (d0, t1); 00540 } 00541 00542 node * differentiate::arccot (application * app, char * derivative) { 00543 _AF0 (f0); 00544 _AD0 (d0); 00545 node * sqr = sqr_reduce (f0->recreate()); 00546 defCon (one, 1); 00547 node * t1 = plus_reduce (one, sqr); 00548 node * t2 = minus_reduce (d0); 00549 return over_reduce (t2, t1); 00550 } 00551 00552 node * differentiate::arcsec (application * app, char * derivative) { 00553 _AF0 (f0); 00554 _AD0 (d0); 00555 node * sqr = sqr_reduce (f0->recreate()); 00556 defCon (one, 1); 00557 node * t1 = minus_reduce (sqr, one); 00558 node * t2 = sqrt_reduce (t1); 00559 node * t3 = times_reduce (f0->recreate(), t2); 00560 return over_reduce (d0, t3); 00561 } 00562 00563 node * differentiate::arccosec (application * app, char * derivative) { 00564 _AF0 (f0); 00565 _AD0 (d0); 00566 node * sqr = sqr_reduce (f0->recreate()); 00567 defCon (one, 1); 00568 node * t1 = minus_reduce (sqr, one); 00569 node * t2 = sqrt_reduce (t1); 00570 node * t3 = times_reduce (f0->recreate(), t2); 00571 node * t4 = minus_reduce (d0); 00572 return over_reduce (t4, t3); 00573 } 00574 00575 node * differentiate::sinh (application * app, char * derivative) { 00576 _AF0 (f0); 00577 _AD0 (d0); 00578 return app_reduce ("cosh", d0, f0->recreate()); 00579 } 00580 00581 node * differentiate::cosh (application * app, char * derivative) { 00582 _AF0 (f0); 00583 _AD0 (d0); 00584 return app_reduce ("sinh", d0, f0->recreate()); 00585 } 00586 00587 node * differentiate::tanh (application * app, char * derivative) { 00588 _AF0 (f0); 00589 _AD0 (d0); 00590 defApp1 (cosh, "cosh", f0->recreate()); 00591 defCon (two, 2); 00592 node * t1 = power_reduce (cosh, two); 00593 return over_reduce (d0, t1); 00594 } 00595 00596 node * differentiate::coth (application * app, char * derivative) { 00597 _AF0 (f0); 00598 _AD0 (d0); 00599 defApp1 (sinh, "sinh", f0->recreate()); 00600 defCon (two, 2); 00601 node * t1 = power_reduce (sinh, two); 00602 node * t2 = minus_reduce (d0); 00603 return over_reduce (t2, t1); 00604 } 00605 00606 node * differentiate::artanh (application * app, char * derivative) { 00607 _AF0 (f0); 00608 _AD0 (d0); 00609 node * sqr = sqr_reduce (f0->recreate()); 00610 defCon (one, 1); 00611 node * t1 = minus_reduce (one, sqr); 00612 return over_reduce (d0, t1); 00613 } 00614 00615 node * differentiate::arcoth (application * app, char * derivative) { 00616 _AF0 (f0); 00617 _AD0 (d0); 00618 node * sqr = sqr_reduce (f0->recreate()); 00619 defCon (one, 1); 00620 node * t1 = minus_reduce (sqr, one); 00621 node * t2 = minus_reduce (d0); 00622 return over_reduce (t2, t1); 00623 } 00624 00625 node * differentiate::arcosh (application * app, char * derivative) { 00626 _AF0 (f0); 00627 _AD0 (d0); 00628 node * sqr = sqr_reduce (f0->recreate()); 00629 defCon (one, 1); 00630 node * t1 = minus_reduce (sqr, one); 00631 node * t2 = sqrt_reduce (t1); 00632 return over_reduce (d0, t2); 00633 } 00634 00635 node * differentiate::arsinh (application * app, char * derivative) { 00636 _AF0 (f0); 00637 _AD0 (d0); 00638 node * sqr = sqr_reduce (f0->recreate()); 00639 defCon (one, 1); 00640 node * t1 = plus_reduce (sqr, one); 00641 node * t2 = sqrt_reduce (t1); 00642 return over_reduce (d0, t2); 00643 } 00644 00645 node * differentiate::arsech (application * app, char * derivative) { 00646 _AF0 (f0); 00647 _AD0 (d0); 00648 node * sqr = sqr_reduce (f0->recreate()); 00649 defCon (one, 1); 00650 node * t1 = minus_reduce (one, sqr); 00651 node * t2 = sqrt_reduce (t1); 00652 node * t3 = times_reduce (f0->recreate(), t2); 00653 node * t4 = minus_reduce (d0); 00654 return over_reduce (t4, t3); 00655 } 00656 00657 node * differentiate::arcosech (application * app, char * derivative) { 00658 _AF0 (f0); 00659 _AD0 (d0); 00660 node * sqr = sqr_reduce (f0->recreate()); 00661 defCon (one, 1); 00662 node * t1 = plus_reduce (one, sqr); 00663 node * t2 = sqrt_reduce (t1); 00664 node * t3 = times_reduce (f0->recreate(), t2); 00665 node * t4 = minus_reduce (d0); 00666 return over_reduce (t4, t3); 00667 } 00668 00669 node * differentiate::ifthenelse (application * app, char * derivative) { 00670 _AF0 (f0); 00671 _AD1 (d1); 00672 _AD2 (d2); 00673 if (isConst (d1) && isConst (d2)) { 00674 if (D(d1) == D(d2)) { 00675 nr_double_t t = D(d1); 00676 delete d1; delete d2; 00677 retCon (t); 00678 } 00679 } 00680 retApp3 ("?:", f0->recreate(), d1, d2); 00681 } 00682 00683 node * differentiate::sinc (application * app, char * derivative) { 00684 _AF0 (f0); 00685 _AD0 (d0); 00686 defApp1 (sinc, "sinc", f0->recreate()); 00687 defApp1 (cos, "cos", f0->recreate()); 00688 node * t1 = minus_reduce (cos, sinc); 00689 node * t2 = over_reduce (t1, f0->recreate()); 00690 return times_reduce (d0, t2); 00691 } 00692 00693 node * differentiate::norm (application * app, char * derivative) { 00694 _AF0 (f0); 00695 _AD0 (d0); 00696 defCon (two, 2); 00697 node * t1 = times_reduce (d0, two); 00698 return times_reduce (t1, f0->recreate()); 00699 } 00700 00701 node * differentiate::xhypot (application * app, char * derivative) { 00702 _AF0 (f0); 00703 _AF1 (f1); 00704 _AD0 (d0); 00705 _AD1 (d1); 00706 node * t1 = hypot_reduce (f0->recreate(), f1->recreate()); 00707 node * t2 = times_reduce (d0, f0->recreate()); 00708 node * t3 = times_reduce (d1, f1->recreate()); 00709 node * t4 = plus_reduce (t2, t3); 00710 return over_reduce (t4, t1); 00711 } 00712 00713 node * differentiate::hypot_reduce (node * f0, node * f1) { 00714 if (isZero (f0) && isZero (f1)) { 00715 delete f0; delete f1; 00716 retCon (0); 00717 } else if (isZero (f0)) { 00718 delete f0; 00719 return sqrt_reduce (sqr_reduce (f1)); 00720 } else if (isZero (f1)) { 00721 delete f1; 00722 return sqrt_reduce (sqr_reduce (f0)); 00723 } else if (isConst (f0) && isConst (f1)) { 00724 nr_double_t t = ::xhypot (D(f0), D(f1)); 00725 delete f0; delete f1; 00726 retCon (t); 00727 } else { 00728 retApp2 ("hypot", f0, f1); 00729 } 00730 } 00731 00732 #include "constants.h" 00733 00734 node * differentiate::vt (application * app, char * derivative) { 00735 _AD0 (d0); 00736 defCon (con, kBoverQ); 00737 return times_reduce (d0, con); 00738 } 00739 00740 // List of differentiators. 00741 struct differentiation_t eqn::differentiations[] = { 00742 { "+", differentiate::plus_binary, 2 }, 00743 { "+", differentiate::plus_unary, 1 }, 00744 { "-", differentiate::minus_binary, 2 }, 00745 { "-", differentiate::minus_unary, 1 }, 00746 { "*", differentiate::times, 2 }, 00747 { "/", differentiate::over, 2 }, 00748 { "^", differentiate::power, 2 }, 00749 00750 { "?:", differentiate::ifthenelse, 3 }, 00751 00752 { "ln", differentiate::ln, 1 }, 00753 { "log10", differentiate::log10, 1 }, 00754 { "log2", differentiate::log2, 1 }, 00755 { "sqrt", differentiate::sqrt, 1 }, 00756 { "exp", differentiate::exp, 1 }, 00757 { "sinc", differentiate::sinc, 1 }, 00758 { "norm", differentiate::norm, 1 }, 00759 { "sin", differentiate::sin, 1 }, 00760 { "cos", differentiate::cos, 1 }, 00761 { "tan", differentiate::tan, 1 }, 00762 { "sec", differentiate::sec, 1 }, 00763 { "cot", differentiate::cot, 1 }, 00764 { "cosec", differentiate::cosec, 1 }, 00765 { "abs", differentiate::abs, 1 }, 00766 { "step", differentiate::step, 1 }, 00767 { "sign", differentiate::sign, 1 }, 00768 { "arcsin", differentiate::arcsin, 1 }, 00769 { "arccos", differentiate::arccos, 1 }, 00770 { "arctan", differentiate::arctan, 1 }, 00771 { "arccot", differentiate::arccot, 1 }, 00772 { "arcsec", differentiate::arcsec, 1 }, 00773 { "arccosec", differentiate::arccosec, 1 }, 00774 { "sqr", differentiate::square, 1 }, 00775 { "sinh", differentiate::sinh, 1 }, 00776 { "cosh", differentiate::cosh, 1 }, 00777 { "tanh", differentiate::tanh, 1 }, 00778 { "coth", differentiate::coth, 1 }, 00779 { "arsinh", differentiate::arsinh, 1 }, 00780 { "arcosh", differentiate::arcosh, 1 }, 00781 { "artanh", differentiate::artanh, 1 }, 00782 { "arcoth", differentiate::arcoth, 1 }, 00783 { "arsech", differentiate::arsech, 1 }, 00784 { "arcosech", differentiate::arcosech, 1 }, 00785 { "hypot", differentiate::xhypot, 2 }, 00786 { "limexp", differentiate::limexp, 1 }, 00787 { "vt", differentiate::vt, 1 }, 00788 00789 { NULL, NULL, 0 /* end of list */ } 00790 };