1 /** @file powerlaws.cpp
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.
6 * GiNaC Copyright (C) 1999 Johannes Gutenberg University Mainz, Germany
8 * This program is free software; you can redistribute it and/or modify
9 * it under the terms of the GNU General Public License as published by
10 * the Free Software Foundation; either version 2 of the License, or
11 * (at your option) any later version.
13 * This program is distributed in the hope that it will be useful,
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 * GNU General Public License for more details.
18 * You should have received a copy of the GNU General Public License
19 * along with this program; if not, write to the Free Software
20 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23 #include <ginac/ginac.h>
25 static unsigned powerlaws1(void)
33 ex e1=power(power(x,a),b);
34 if (!(is_ex_exactly_of_type(e1,power) &&
35 is_ex_exactly_of_type(e1.op(0),power) &&
36 is_ex_exactly_of_type(e1.op(0).op(0),symbol) &&
37 is_ex_exactly_of_type(e1.op(0).op(1),symbol) &&
38 is_ex_exactly_of_type(e1.op(1),symbol) &&
39 e1.is_equal(power(power(x,a),b)) )) {
40 clog << "(x^a)^b, x,a,b symbolic wrong" << endl;
41 clog << "returned: " << e1 << endl;
46 if (!(is_ex_exactly_of_type(e2,power) &&
47 is_ex_exactly_of_type(e2.op(0),symbol) &&
48 is_ex_exactly_of_type(e2.op(1),symbol) &&
49 e2.is_equal(power(x,b)) )) {
50 clog << "(x^a)^b, x,b symbolic, a==1 wrong" << endl;
51 clog << "returned: " << e2 << endl;
56 if (!(is_ex_exactly_of_type(e3,power) &&
57 is_ex_exactly_of_type(e3.op(0),power) &&
58 is_ex_exactly_of_type(e3.op(0).op(0),symbol) &&
59 is_ex_exactly_of_type(e3.op(0).op(1),numeric) &&
60 is_ex_exactly_of_type(e3.op(1),symbol) &&
61 e3.is_equal(power(power(x,-1),b)) )) {
62 clog << "(x^a)^b, x,b symbolic, a==-1 wrong" << endl;
63 clog << "returned: " << e3 << endl;
67 ex e4=e1.subs(lst(a==-1,b==2.5));
68 if (!(is_ex_exactly_of_type(e4,power) &&
69 is_ex_exactly_of_type(e4.op(0),power) &&
70 is_ex_exactly_of_type(e4.op(0).op(0),symbol) &&
71 is_ex_exactly_of_type(e4.op(0).op(1),numeric) &&
72 is_ex_exactly_of_type(e4.op(1),numeric) &&
73 e4.is_equal(power(power(x,-1),2.5)) )) {
74 clog << "(x^a)^b, x symbolic, a==-1, b==2.5 wrong" << endl;
75 clog << "returned: " << e4 << endl;
79 ex e5=e1.subs(lst(a==-0.9,b==2.5));
80 if (!(is_ex_exactly_of_type(e5,power) &&
81 is_ex_exactly_of_type(e5.op(0),symbol) &&
82 is_ex_exactly_of_type(e5.op(1),numeric) &&
83 e5.is_equal(power(x,numeric(-0.9)*numeric(2.5))) )) {
84 clog << "(x^a)^b, x symbolic, a==-0.9, b==2.5 wrong" << endl;
85 clog << "returned: " << e5 << endl;
89 ex e6=e1.subs(lst(a==numeric(3)+numeric(5.3)*I,b==-5));
90 if (!(is_ex_exactly_of_type(e6,power) &&
91 is_ex_exactly_of_type(e6.op(0),symbol) &&
92 is_ex_exactly_of_type(e6.op(1),numeric) &&
93 e6.is_equal(power(x,numeric(-15)+numeric(5.3)*numeric(-5)*I)) )) {
94 clog << "(x^a)^b, x symbolic, a==3+5.3*I, b==-5 wrong" << endl;
95 clog << "returned: " << e6 << endl;
101 static unsigned powerlaws2(void)
103 // (a*x)^b = a^b * x^b
110 if (!(is_ex_exactly_of_type(e1,power) &&
111 is_ex_exactly_of_type(e1.op(0),mul) &&
112 (e1.op(0).nops()==2) &&
113 is_ex_exactly_of_type(e1.op(0).op(0),symbol) &&
114 is_ex_exactly_of_type(e1.op(0).op(1),symbol) &&
115 is_ex_exactly_of_type(e1.op(1),symbol) &&
116 e1.is_equal(power(a*x,b)) )) {
117 clog << "(a*x)^b, x,a,b symbolic wrong" << endl;
118 clog << "returned: " << e1 << endl;
123 if (!(is_ex_exactly_of_type(e2,power) &&
124 is_ex_exactly_of_type(e2.op(0),mul) &&
125 (e2.op(0).nops()==2) &&
126 is_ex_exactly_of_type(e2.op(0).op(0),symbol) &&
127 is_ex_exactly_of_type(e2.op(0).op(1),numeric) &&
128 is_ex_exactly_of_type(e2.op(1),symbol) &&
129 e2.is_equal(power(3*x,b)) )) {
130 clog << "(a*x)^b, x,b symbolic, a==3 wrong" << endl;
131 clog << "returned: " << e2 << endl;
135 ex e3=e1.subs(b==-3);
136 if (!(is_ex_exactly_of_type(e3,mul) &&
138 is_ex_exactly_of_type(e3.op(0),power) &&
139 is_ex_exactly_of_type(e3.op(1),power) &&
140 e3.is_equal(power(a,-3)*power(x,-3)) )) {
141 clog << "(a*x)^b, x,a symbolic, b==-3 wrong" << endl;
142 clog << "returned: " << e3 << endl;
146 ex e4=e1.subs(b==4.5);
147 if (!(is_ex_exactly_of_type(e4,power) &&
148 is_ex_exactly_of_type(e4.op(0),mul) &&
149 (e4.op(0).nops()==2) &&
150 is_ex_exactly_of_type(e4.op(0).op(0),symbol) &&
151 is_ex_exactly_of_type(e4.op(0).op(1),symbol) &&
152 is_ex_exactly_of_type(e4.op(1),numeric) &&
153 e4.is_equal(power(a*x,4.5)) )) {
154 clog << "(a*x)^b, x,a symbolic, b==4.5 wrong" << endl;
155 clog << "returned: " << e4 << endl;
159 ex e5=e1.subs(lst(a==3.2,b==3+numeric(5)*I));
160 if (!(is_ex_exactly_of_type(e5,mul) &&
162 is_ex_exactly_of_type(e5.op(0),power) &&
163 is_ex_exactly_of_type(e5.op(1),numeric) &&
164 e5.is_equal(power(x,3+numeric(5)*I)*
165 power(numeric(3.2),3+numeric(5)*I)) )) {
166 clog << "(a*x)^b, x symbolic, a==3.2, b==3+5*I wrong" << endl;
167 clog << "returned: " << e5 << endl;
171 ex e6=e1.subs(lst(a==-3.2,b==3+numeric(5)*I));
172 if (!(is_ex_exactly_of_type(e6,mul) &&
174 is_ex_exactly_of_type(e6.op(0),power) &&
175 is_ex_exactly_of_type(e6.op(1),numeric) &&
176 e6.is_equal(power(-x,3+numeric(5)*I)*
177 power(numeric(3.2),3+numeric(5)*I)) )) {
178 clog << "(a*x)^b, x symbolic, a==-3.2, b==3+5*I wrong" << endl;
179 clog << "returned: " << e6 << endl;
183 ex e7=e1.subs(lst(a==3+numeric(5)*I,b==3.2));
184 if (!(is_ex_exactly_of_type(e7,power) &&
185 is_ex_exactly_of_type(e7.op(0),mul) &&
186 (e7.op(0).nops()==2) &&
187 is_ex_exactly_of_type(e7.op(0).op(0),symbol) &&
188 is_ex_exactly_of_type(e7.op(0).op(1),numeric) &&
189 is_ex_exactly_of_type(e7.op(1),numeric) &&
190 e7.is_equal(power((3+numeric(5)*I)*x,3.2)) )) {
191 clog << "(a*x)^b, x symbolic, a==3+5*I, b==3.2 wrong" << endl;
192 clog << "returned: " << e7 << endl;
199 static unsigned powerlaws3(void)
201 // numeric evaluation
203 ex e1=power(numeric(4),numeric(1)/numeric(2));
205 clog << "4^(1/2) wrongly returned " << e1 << endl;
209 ex e2=power(numeric(27),numeric(2)/numeric(3));
211 clog << "27^(2/3) wrongly returned " << e2 << endl;
215 ex e3=power(numeric(5),numeric(1)/numeric(2));
216 if (!(is_ex_exactly_of_type(e3,power) &&
217 e3.op(0).is_equal(numeric(5)) &&
218 e3.op(1).is_equal(numeric(1)/numeric(2)) )) {
219 clog << "5^(1/2) wrongly returned " << e3 << endl;
223 ex e4=power(numeric(5),evalf(numeric(1)/numeric(2)));
224 if (!(is_ex_exactly_of_type(e4,numeric))) {
225 clog << "5^(0.5) wrongly returned " << e4 << endl;
229 ex e5=power(evalf(numeric(5)),numeric(1)/numeric(2));
230 if (!(is_ex_exactly_of_type(e5,numeric))) {
231 clog << "5.0^(1/2) wrongly returned " << e5 << endl;
238 static unsigned powerlaws4(void)
240 // test for mul::eval()
246 ex f1=power(a*b,ex(1)/ex(2));
247 ex f2=power(a*b,ex(3)/ex(2));
256 clog << "(a*b)^(1/2)*(a*b)^(3/2)*c wrongly returned " << e1 << endl;
262 unsigned powerlaws(void)
266 cout << "checking power laws..." << flush;
267 clog << "---------power laws:" << endl;
269 result += powerlaws1();
270 result += powerlaws2();
271 result += powerlaws3();
272 result += powerlaws4();
276 clog << "(no output)" << endl;