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. */
7 * GiNaC Copyright (C) 1999 Johannes Gutenberg University Mainz, Germany
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.
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.
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
24 #include <ginac/ginac.h>
26 #ifndef NO_GINAC_NAMESPACE
27 using namespace GiNaC;
28 #endif // ndef NO_GINAC_NAMESPACE
30 static unsigned powerlaws1(void)
38 ex e1=power(power(x,a),b);
39 if (!(is_ex_exactly_of_type(e1,power) &&
40 is_ex_exactly_of_type(e1.op(0),power) &&
41 is_ex_exactly_of_type(e1.op(0).op(0),symbol) &&
42 is_ex_exactly_of_type(e1.op(0).op(1),symbol) &&
43 is_ex_exactly_of_type(e1.op(1),symbol) &&
44 e1.is_equal(power(power(x,a),b)) )) {
45 clog << "(x^a)^b, x,a,b symbolic wrong" << endl;
46 clog << "returned: " << e1 << endl;
51 if (!(is_ex_exactly_of_type(e2,power) &&
52 is_ex_exactly_of_type(e2.op(0),symbol) &&
53 is_ex_exactly_of_type(e2.op(1),symbol) &&
54 e2.is_equal(power(x,b)) )) {
55 clog << "(x^a)^b, x,b symbolic, a==1 wrong" << endl;
56 clog << "returned: " << e2 << endl;
61 if (!(is_ex_exactly_of_type(e3,power) &&
62 is_ex_exactly_of_type(e3.op(0),power) &&
63 is_ex_exactly_of_type(e3.op(0).op(0),symbol) &&
64 is_ex_exactly_of_type(e3.op(0).op(1),numeric) &&
65 is_ex_exactly_of_type(e3.op(1),symbol) &&
66 e3.is_equal(power(power(x,-1),b)) )) {
67 clog << "(x^a)^b, x,b symbolic, a==-1 wrong" << endl;
68 clog << "returned: " << e3 << endl;
72 ex e4=e1.subs(lst(a==-1,b==2.5));
73 if (!(is_ex_exactly_of_type(e4,power) &&
74 is_ex_exactly_of_type(e4.op(0),power) &&
75 is_ex_exactly_of_type(e4.op(0).op(0),symbol) &&
76 is_ex_exactly_of_type(e4.op(0).op(1),numeric) &&
77 is_ex_exactly_of_type(e4.op(1),numeric) &&
78 e4.is_equal(power(power(x,-1),2.5)) )) {
79 clog << "(x^a)^b, x symbolic, a==-1, b==2.5 wrong" << endl;
80 clog << "returned: " << e4 << endl;
84 ex e5=e1.subs(lst(a==-0.9,b==2.5));
85 if (!(is_ex_exactly_of_type(e5,power) &&
86 is_ex_exactly_of_type(e5.op(0),symbol) &&
87 is_ex_exactly_of_type(e5.op(1),numeric) &&
88 e5.is_equal(power(x,numeric(-0.9)*numeric(2.5))) )) {
89 clog << "(x^a)^b, x symbolic, a==-0.9, b==2.5 wrong" << endl;
90 clog << "returned: " << e5 << endl;
94 ex e6=e1.subs(lst(a==numeric(3)+numeric(5.3)*I,b==-5));
95 if (!(is_ex_exactly_of_type(e6,power) &&
96 is_ex_exactly_of_type(e6.op(0),symbol) &&
97 is_ex_exactly_of_type(e6.op(1),numeric) &&
98 e6.is_equal(power(x,numeric(-15)+numeric(5.3)*numeric(-5)*I)) )) {
99 clog << "(x^a)^b, x symbolic, a==3+5.3*I, b==-5 wrong" << endl;
100 clog << "returned: " << e6 << endl;
106 static unsigned powerlaws2(void)
108 // (a*x)^b = a^b * x^b
115 if (!(is_ex_exactly_of_type(e1,power) &&
116 is_ex_exactly_of_type(e1.op(0),mul) &&
117 (e1.op(0).nops()==2) &&
118 is_ex_exactly_of_type(e1.op(0).op(0),symbol) &&
119 is_ex_exactly_of_type(e1.op(0).op(1),symbol) &&
120 is_ex_exactly_of_type(e1.op(1),symbol) &&
121 e1.is_equal(power(a*x,b)) )) {
122 clog << "(a*x)^b, x,a,b symbolic wrong" << endl;
123 clog << "returned: " << e1 << endl;
128 if (!(is_ex_exactly_of_type(e2,power) &&
129 is_ex_exactly_of_type(e2.op(0),mul) &&
130 (e2.op(0).nops()==2) &&
131 is_ex_exactly_of_type(e2.op(0).op(0),symbol) &&
132 is_ex_exactly_of_type(e2.op(0).op(1),numeric) &&
133 is_ex_exactly_of_type(e2.op(1),symbol) &&
134 e2.is_equal(power(3*x,b)) )) {
135 clog << "(a*x)^b, x,b symbolic, a==3 wrong" << endl;
136 clog << "returned: " << e2 << endl;
140 ex e3=e1.subs(b==-3);
141 if (!(is_ex_exactly_of_type(e3,mul) &&
143 is_ex_exactly_of_type(e3.op(0),power) &&
144 is_ex_exactly_of_type(e3.op(1),power) &&
145 e3.is_equal(power(a,-3)*power(x,-3)) )) {
146 clog << "(a*x)^b, x,a symbolic, b==-3 wrong" << endl;
147 clog << "returned: " << e3 << endl;
151 ex e4=e1.subs(b==4.5);
152 if (!(is_ex_exactly_of_type(e4,power) &&
153 is_ex_exactly_of_type(e4.op(0),mul) &&
154 (e4.op(0).nops()==2) &&
155 is_ex_exactly_of_type(e4.op(0).op(0),symbol) &&
156 is_ex_exactly_of_type(e4.op(0).op(1),symbol) &&
157 is_ex_exactly_of_type(e4.op(1),numeric) &&
158 e4.is_equal(power(a*x,4.5)) )) {
159 clog << "(a*x)^b, x,a symbolic, b==4.5 wrong" << endl;
160 clog << "returned: " << e4 << endl;
164 ex e5=e1.subs(lst(a==3.2,b==3+numeric(5)*I));
165 if (!(is_ex_exactly_of_type(e5,mul) &&
167 is_ex_exactly_of_type(e5.op(0),power) &&
168 is_ex_exactly_of_type(e5.op(1),numeric) &&
169 e5.is_equal(power(x,3+numeric(5)*I)*
170 power(numeric(3.2),3+numeric(5)*I)) )) {
171 clog << "(a*x)^b, x symbolic, a==3.2, b==3+5*I wrong" << endl;
172 clog << "returned: " << e5 << endl;
176 ex e6=e1.subs(lst(a==-3.2,b==3+numeric(5)*I));
177 if (!(is_ex_exactly_of_type(e6,mul) &&
179 is_ex_exactly_of_type(e6.op(0),power) &&
180 is_ex_exactly_of_type(e6.op(1),numeric) &&
181 e6.is_equal(power(-x,3+numeric(5)*I)*
182 power(numeric(3.2),3+numeric(5)*I)) )) {
183 clog << "(a*x)^b, x symbolic, a==-3.2, b==3+5*I wrong" << endl;
184 clog << "returned: " << e6 << endl;
188 ex e7=e1.subs(lst(a==3+numeric(5)*I,b==3.2));
189 if (!(is_ex_exactly_of_type(e7,power) &&
190 is_ex_exactly_of_type(e7.op(0),mul) &&
191 (e7.op(0).nops()==2) &&
192 is_ex_exactly_of_type(e7.op(0).op(0),symbol) &&
193 is_ex_exactly_of_type(e7.op(0).op(1),numeric) &&
194 is_ex_exactly_of_type(e7.op(1),numeric) &&
195 e7.is_equal(power((3+numeric(5)*I)*x,3.2)) )) {
196 clog << "(a*x)^b, x symbolic, a==3+5*I, b==3.2 wrong" << endl;
197 clog << "returned: " << e7 << endl;
204 static unsigned powerlaws3(void)
206 // numeric evaluation
208 ex e1=power(numeric(4),numeric(1)/numeric(2));
210 clog << "4^(1/2) wrongly returned " << e1 << endl;
214 ex e2=power(numeric(27),numeric(2)/numeric(3));
216 clog << "27^(2/3) wrongly returned " << e2 << endl;
220 ex e3=power(numeric(5),numeric(1)/numeric(2));
221 if (!(is_ex_exactly_of_type(e3,power) &&
222 e3.op(0).is_equal(numeric(5)) &&
223 e3.op(1).is_equal(numeric(1)/numeric(2)) )) {
224 clog << "5^(1/2) wrongly returned " << e3 << endl;
228 ex e4=power(numeric(5),evalf(numeric(1)/numeric(2)));
229 if (!(is_ex_exactly_of_type(e4,numeric))) {
230 clog << "5^(0.5) wrongly returned " << e4 << endl;
234 ex e5=power(evalf(numeric(5)),numeric(1)/numeric(2));
235 if (!(is_ex_exactly_of_type(e5,numeric))) {
236 clog << "5.0^(1/2) wrongly returned " << e5 << endl;
243 static unsigned powerlaws4(void)
245 // test for mul::eval()
251 ex f1=power(a*b,ex(1)/ex(2));
252 ex f2=power(a*b,ex(3)/ex(2));
261 clog << "(a*b)^(1/2)*(a*b)^(3/2)*c wrongly returned " << e1 << endl;
267 unsigned powerlaws(void)
271 cout << "checking power laws..." << flush;
272 clog << "---------power laws:" << endl;
274 result += powerlaws1();
275 result += powerlaws2();
276 result += powerlaws3();
277 result += powerlaws4();
281 clog << "(no output)" << endl;