check/exam_powerlaws.cpp

   1 /** @file exam_powerlaws.cpp
   2  *
   3  *  Tests for power laws.  You shouldn't try to draw much inspiration from
   4  *  this code, it is a sanity check rather deeply rooted in GiNaC's classes. */
   5
   6 /*
   7  *  GiNaC Copyright (C) 1999-2000 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 "exams.h"
  25
  26 static unsigned exam_powerlaws1(void)
  27 {
  28         // (x^a)^b = x^(a*b)
  29
  30         symbol x("x");
  31         symbol a("a");
  32         symbol b("b");
  33
  34         ex e1=power(power(x,a),b);
  35         if (!(is_ex_exactly_of_type(e1,power) &&
  36                   is_ex_exactly_of_type(e1.op(0),power) &&
  37                   is_ex_exactly_of_type(e1.op(0).op(0),symbol) &&
  38                   is_ex_exactly_of_type(e1.op(0).op(1),symbol) &&
  39                   is_ex_exactly_of_type(e1.op(1),symbol) &&
  40                   e1.is_equal(power(power(x,a),b)) )) {
  41                 clog << "(x^a)^b, x,a,b symbolic wrong" << endl;
  42                 clog << "returned: " << e1 << endl;
  43                 return 1;
  44         }
  45
  46         ex e2=e1.subs(a==1);
  47         if (!(is_ex_exactly_of_type(e2,power) &&
  48                   is_ex_exactly_of_type(e2.op(0),symbol) &&
  49                   is_ex_exactly_of_type(e2.op(1),symbol) &&
  50                   e2.is_equal(power(x,b)) )) {
  51                 clog << "(x^a)^b, x,b symbolic, a==1 wrong" << endl;
  52                 clog << "returned: " << e2 << endl;
  53                 return 1;
  54         }
  55
  56         ex e3=e1.subs(a==-1);
  57         if (!(is_ex_exactly_of_type(e3,power) &&
  58                   is_ex_exactly_of_type(e3.op(0),power) &&
  59                   is_ex_exactly_of_type(e3.op(0).op(0),symbol) &&
  60                   is_ex_exactly_of_type(e3.op(0).op(1),numeric) &&
  61                   is_ex_exactly_of_type(e3.op(1),symbol) &&
  62                   e3.is_equal(power(power(x,-1),b)) )) {
  63                 clog << "(x^a)^b, x,b symbolic, a==-1 wrong" << endl;
  64                 clog << "returned: " << e3 << endl;
  65                 return 1;
  66         }
  67
  68         ex e4=e1.subs(lst(a==-1,b==2.5));
  69         if (!(is_ex_exactly_of_type(e4,power) &&
  70                   is_ex_exactly_of_type(e4.op(0),power) &&
  71                   is_ex_exactly_of_type(e4.op(0).op(0),symbol) &&
  72                   is_ex_exactly_of_type(e4.op(0).op(1),numeric) &&
  73                   is_ex_exactly_of_type(e4.op(1),numeric) &&
  74                   e4.is_equal(power(power(x,-1),2.5)) )) {
  75                 clog << "(x^a)^b, x symbolic, a==-1, b==2.5 wrong" << endl;
  76                 clog << "returned: " << e4 << endl;
  77                 return 1;
  78         }
  79
  80         ex e5=e1.subs(lst(a==-0.9,b==2.5));
  81         if (!(is_ex_exactly_of_type(e5,power) &&
  82                   is_ex_exactly_of_type(e5.op(0),symbol) &&
  83                   is_ex_exactly_of_type(e5.op(1),numeric) &&
  84                   e5.is_equal(power(x,numeric(-0.9)*numeric(2.5))) )) {
  85                 clog << "(x^a)^b, x symbolic, a==-0.9, b==2.5 wrong" << endl;
  86                 clog << "returned: " << e5 << endl;
  87                 return 1;
  88         }
  89
  90         ex e6=e1.subs(lst(a==numeric(3)+numeric(5.3)*I,b==-5));
  91         if (!(is_ex_exactly_of_type(e6,power) &&
  92                   is_ex_exactly_of_type(e6.op(0),symbol) &&
  93                   is_ex_exactly_of_type(e6.op(1),numeric) &&
  94                   e6.is_equal(power(x,numeric(-15)+numeric(5.3)*numeric(-5)*I)) )) {
  95                 clog << "(x^a)^b, x symbolic, a==3+5.3*I, b==-5 wrong" << endl;
  96                 clog << "returned: " << e6 << endl;
  97                 return 1;
  98         }
  99
 100         return 0;
 101 }
 102
 103 static unsigned exam_powerlaws2(void)
 104 {
 105         // (a*x)^b = a^b * x^b
 106
 107         symbol x("x");
 108         symbol a("a");
 109         symbol b("b");
 110
 111         ex e1=power(a*x,b);
 112         if (!(is_ex_exactly_of_type(e1,power) &&
 113                   is_ex_exactly_of_type(e1.op(0),mul) &&
 114                   (e1.op(0).nops()==2) &&
 115                   is_ex_exactly_of_type(e1.op(0).op(0),symbol) &&
 116                   is_ex_exactly_of_type(e1.op(0).op(1),symbol) &&
 117                   is_ex_exactly_of_type(e1.op(1),symbol) &&
 118                   e1.is_equal(power(a*x,b)) )) {
 119                 clog << "(a*x)^b, x,a,b symbolic wrong" << endl;
 120                 clog << "returned: " << e1 << endl;
 121                 return 1;
 122         }
 123
 124         ex e2=e1.subs(a==3);
 125         if (!(is_ex_exactly_of_type(e2,power) &&
 126                   is_ex_exactly_of_type(e2.op(0),mul) &&
 127                   (e2.op(0).nops()==2) &&
 128                   is_ex_exactly_of_type(e2.op(0).op(0),symbol) &&
 129                   is_ex_exactly_of_type(e2.op(0).op(1),numeric) &&
 130                   is_ex_exactly_of_type(e2.op(1),symbol) &&
 131                   e2.is_equal(power(3*x,b)) )) {
 132                 clog << "(a*x)^b, x,b symbolic, a==3 wrong" << endl;
 133                 clog << "returned: " << e2 << endl;
 134                 return 1;
 135         }
 136
 137         ex e3=e1.subs(b==-3);
 138         if (!(is_ex_exactly_of_type(e3,mul) &&
 139                   (e3.nops()==2) &&
 140                   is_ex_exactly_of_type(e3.op(0),power) &&
 141                   is_ex_exactly_of_type(e3.op(1),power) &&
 142                   e3.is_equal(power(a,-3)*power(x,-3)) )) {
 143                 clog << "(a*x)^b, x,a symbolic, b==-3 wrong" << endl;
 144                 clog << "returned: " << e3 << endl;
 145                 return 1;
 146         }
 147
 148         ex e4=e1.subs(b==4.5);
 149         if (!(is_ex_exactly_of_type(e4,power) &&
 150                   is_ex_exactly_of_type(e4.op(0),mul) &&
 151                   (e4.op(0).nops()==2) &&
 152                   is_ex_exactly_of_type(e4.op(0).op(0),symbol) &&
 153                   is_ex_exactly_of_type(e4.op(0).op(1),symbol) &&
 154                   is_ex_exactly_of_type(e4.op(1),numeric) &&
 155                   e4.is_equal(power(a*x,4.5)) )) {
 156                 clog << "(a*x)^b, x,a symbolic, b==4.5 wrong" << endl;
 157                 clog << "returned: " << e4 << endl;
 158                 return 1;
 159         }
 160
 161         ex e5=e1.subs(lst(a==3.2,b==3+numeric(5)*I));
 162         if (!(is_ex_exactly_of_type(e5,mul) &&
 163                   (e5.nops()==2) &&
 164                   is_ex_exactly_of_type(e5.op(0),power) &&
 165                   is_ex_exactly_of_type(e5.op(1),numeric) &&
 166                   e5.is_equal(power(x,3+numeric(5)*I)*
 167                                           power(numeric(3.2),3+numeric(5)*I)) )) {
 168                 clog << "(a*x)^b, x symbolic, a==3.2, b==3+5*I wrong" << endl;
 169                 clog << "returned: " << e5 << endl;
 170                 return 1;
 171         }
 172
 173         ex e6=e1.subs(lst(a==-3.2,b==3+numeric(5)*I));
 174         if (!(is_ex_exactly_of_type(e6,mul) &&
 175                   (e6.nops()==2) &&
 176                   is_ex_exactly_of_type(e6.op(0),power) &&
 177                   is_ex_exactly_of_type(e6.op(1),numeric) &&
 178                   e6.is_equal(power(-x,3+numeric(5)*I)*
 179                                           power(numeric(3.2),3+numeric(5)*I)) )) {
 180                 clog << "(a*x)^b, x symbolic, a==-3.2, b==3+5*I wrong" << endl;
 181                 clog << "returned: " << e6 << endl;
 182                 return 1;
 183         }
 184
 185         ex e7=e1.subs(lst(a==3+numeric(5)*I,b==3.2));
 186         if (!(is_ex_exactly_of_type(e7,power) &&
 187                   is_ex_exactly_of_type(e7.op(0),mul) &&
 188                   (e7.op(0).nops()==2) &&
 189                   is_ex_exactly_of_type(e7.op(0).op(0),symbol) &&
 190                   is_ex_exactly_of_type(e7.op(0).op(1),numeric) &&
 191                   is_ex_exactly_of_type(e7.op(1),numeric) &&
 192                   e7.is_equal(power((3+numeric(5)*I)*x,3.2)) )) {
 193                 clog << "(a*x)^b, x symbolic, a==3+5*I, b==3.2 wrong" << endl;
 194                 clog << "returned: " << e7 << endl;
 195                 return 1;
 196         }
 197
 198         return 0;
 199 }
 200
 201 static unsigned exam_powerlaws3(void)
 202 {
 203         // numeric evaluation
 204
 205         ex e1 = power(numeric(4),numeric(1,2));
 206         if (e1 != 2) {
 207                 clog << "4^(1/2) wrongly returned " << e1 << endl;
 208                 return 1;
 209         }
 210
 211         ex e2 = power(numeric(27),numeric(2,3));
 212         if (e2 != 9) {
 213                 clog << "27^(2/3) wrongly returned " << e2 << endl;
 214                 return 1;
 215         }
 216
 217         ex e3 = power(numeric(5),numeric(1,2));
 218         if (!(is_ex_exactly_of_type(e3,power) &&
 219                   e3.op(0).is_equal(numeric(5)) &&
 220                   e3.op(1).is_equal(numeric(1,2)))) {
 221                 clog << "5^(1/2) wrongly returned " << e3 << endl;
 222                 return 1;
 223         }
 224
 225         ex e4 = power(numeric(5),evalf(numeric(1,2)));
 226         if (!(is_ex_exactly_of_type(e4,numeric))) {
 227                 clog << "5^(0.5) wrongly returned " << e4 << endl;
 228                 return 1;
 229         }
 230
 231         ex e5 = power(evalf(numeric(5)),numeric(1,2));
 232         if (!(is_ex_exactly_of_type(e5,numeric))) {
 233                 clog << "5.0^(1/2) wrongly returned " << e5 << endl;
 234                 return 1;
 235         }
 236
 237         return 0;
 238 }
 239
 240 static unsigned exam_powerlaws4(void)
 241 {
 242         // test for mul::eval()
 243
 244         symbol a("a");
 245         symbol b("b");
 246         symbol c("c");
 247
 248         ex f1 = power(a*b,ex(1)/ex(2));
 249         ex f2 = power(a*b,ex(3)/ex(2));
 250         ex f3 = c;
 251
 252         exvector v;
 253         v.push_back(f1);
 254         v.push_back(f2);
 255         v.push_back(f3);
 256         ex e1 = mul(v);
 257         if (e1!=a*a*b*b*c) {
 258                 clog << "(a*b)^(1/2)*(a*b)^(3/2)*c wrongly returned " << e1 << endl;
 259                 return 1;
 260         }
 261
 262         return 0;
 263 }
 264
 265 static unsigned exam_powerlaws5(void)
 266 {
 267         // cabinet of slightly pathological cases
 268
 269         symbol a("a");
 270
 271         ex e1 = pow(1,a);
 272         if (e1 != 1) {
 273                 clog << "1^a wrongly returned " << e1 << endl;
 274                 return 1;
 275         }
 276
 277         ex e2 = pow(0,a);
 278         if (!(is_ex_exactly_of_type(e2,power))) {
 279                 clog << "0^a was evaluated to " << e2
 280                          << " though nothing is known about a." << endl;
 281                 return 1;
 282         }
 283
 284         return 0;
 285 }
 286
 287 unsigned exam_powerlaws(void)
 288 {
 289         unsigned result = 0;
 290
 291         cout << "examining power laws" << flush;
 292         clog << "----------power laws:" << endl;
 293
 294         result += exam_powerlaws1();  cout << '.' << flush;
 295         result += exam_powerlaws2();  cout << '.' << flush;
 296         result += exam_powerlaws3();  cout << '.' << flush;
 297         result += exam_powerlaws4();  cout << '.' << flush;
 298         result += exam_powerlaws5();  cout << '.' << flush;
 299
 300         if (!result) {
 301                 cout << " passed " << endl;
 302                 clog << "(no output)" << endl;
 303         } else {
 304                 cout << " failed " << endl;
 305         }
 306
 307         return result;
 308 }