This file records noteworthy changes.
0.9.1 (<date>)
+* Ctors of class numeric are not explicit any more. All built-in callers for
+ pseudofunctions are now templated and default to ex arguments which relaxes
+ the need for explicit ctors.
* New functions/methods:
- find()
- remove_first(), remove_last(), sort() and unique() for lists
int main()
@{
- numeric two(2); // exact integer 2
+ numeric two = 2; // exact integer 2
numeric r(2,3); // exact fraction 2/3
numeric e(2.71828); // floating point number
- numeric p("3.1415926535897932385"); // floating point number
+ numeric p = "3.14159265358979323846"; // constructor from string
// Trott's constant in scientific notation:
numeric trott("1.0841015122311136151E-2");
@}
@end example
-Note that all those constructors are @emph{explicit} which means you are
-not allowed to write @code{numeric two=2;}. This is because the basic
-objects to be handled by GiNaC are the expressions @code{ex} and we want
-to keep things simple and wish objects like @code{pow(x,2)} to be
-handled the same way as @code{pow(x,a)}, which means that we need to
-allow a general @code{ex} as base and exponent. Therefore there is an
-implicit constructor from C-integers directly to expressions handling
-numerics at work in most of our examples. This design really becomes
-convenient when one declares own functions having more than one
-parameter but it forbids using implicit constructors because that would
-lead to compile-time ambiguities.
-
It may be tempting to construct numbers writing @code{numeric r(3/2)}.
This would, however, call C's built-in operator @code{/} for integers
first and result in a numeric holding a plain integer 1. @strong{Never
instance, all trigonometric and hyperbolic functions are implemented
(@xref{Built-in Functions}, for a complete list).
-These functions are all objects of class @code{function}. They accept
-one or more expressions as arguments and return one expression. If the
-arguments are not numerical, the evaluation of the function may be
-halted, as it does in the next example, showing how a function returns
-itself twice and finally an expression that may be really useful:
+These functions (better called @emph{pseudofunctions}) are all objects
+of class @code{function}. They accept one or more expressions as
+arguments and return one expression. If the arguments are not
+numerical, the evaluation of the function may be halted, as it does in
+the next example, showing how a function returns itself twice and
+finally an expression that may be really useful:
@cindex Gamma function
@cindex @code{subs()}
expansion and so on. Read the next chapter in order to learn more about
this.
+It must be noted that these pseudofunctions are created by inline
+functions, where the argument list is templated. This means that
+whenever you call @code{GiNaC::sin(1)} it is equivalent to
+@code{sin(ex(1))} and will therefore not result in a floating point
+numeber. Unless of course the function prototype is explicitly
+overridden -- which is the case for arguments of type @code{numeric}
+(not wrapped inside an @code{ex}). Hence, in order to obtain a floating
+point number of class @code{numeric} you should call
+@code{sin(numeric(1))}. This is almost the same as calling
+@code{sin(1).evalf()} except that the latter will return a numeric
+wrapped inside an @code{ex}.
+
@node Relations, Matrices, Mathematical functions, Basic Concepts
@c node-name, next, previous, up
// possible from then on the sequence has changed and must be compacted
bool must_copy = false;
while (itin2!=last) {
- if ((*itin1).rest.compare(itin2->rest)==0) {
- (*itin1).coeff = ex_to<numeric>(itin1->coeff).
- add_dyn(ex_to<numeric>(itin2->coeff));
+ if (itin1->rest.compare(itin2->rest)==0) {
+ itin1->coeff = ex_to<numeric>(itin1->coeff).
+ add_dyn(ex_to<numeric>(itin2->coeff));
if (expair_needs_further_processing(itin1))
needs_further_processing = true;
must_copy = true;
}
sub generate_from_to {
- my ($template,$seq_template1,$seq_template2,$from,$to)=@_;
+ my ($template,$seq_template1,$seq_template2,$seq_template3,$from,$to)=@_;
my ($res,$N,$SEQ);
$res='';
for ($N=$from; $N<=$to; $N++) {
$SEQ1=generate_seq($seq_template1,$N);
$SEQ2=generate_seq($seq_template2,$N);
+ $SEQ3=generate_seq($seq_template3,$N);
$res .= eval('"' . $template . '"');
$SEQ1=''; # to avoid main::SEQ1 used only once warning
$SEQ2=''; # same as above
+ $SEQ3=''; # same as above
}
return $res;
}
sub generate {
- my ($template,$seq_template1,$seq_template2)=@_;
- return generate_from_to($template,$seq_template1,$seq_template2,1,$maxargs);
+ my ($template,$seq_template1,$seq_template2,$seq_template3)=@_;
+ return generate_from_to($template,$seq_template1,$seq_template2,$seq_template3,1,$maxargs);
}
-$declare_function_macro = generate_from_to(
- <<'END_OF_DECLARE_FUNCTION_MACRO','const GiNaC::ex & p${N}','p${N}',1,$maxargs);
+$declare_function_macro = generate(
+ <<'END_OF_DECLARE_FUNCTION_MACRO','typename T${N}','const T${N} & p${N}','GiNaC::ex(p${N})');
#define DECLARE_FUNCTION_${N}P(NAME) \\
extern const unsigned function_index_##NAME; \\
-inline const GiNaC::function NAME(${SEQ1}) { \\
- return GiNaC::function(function_index_##NAME, ${SEQ2}); \\
+template<${SEQ1}> \\
+inline const GiNaC::function NAME(${SEQ2}) { \\
+ return GiNaC::function(function_index_##NAME, ${SEQ3}); \\
}
END_OF_DECLARE_FUNCTION_MACRO
$typedef_eval_funcp=generate(
'typedef ex (* eval_funcp_${N})(${SEQ1});'."\n",
-'const ex &','');
+'const ex &','','');
$typedef_evalf_funcp=generate(
'typedef ex (* evalf_funcp_${N})(${SEQ1});'."\n",
-'const ex &','');
+'const ex &','','');
$typedef_derivative_funcp=generate(
'typedef ex (* derivative_funcp_${N})(${SEQ1}, unsigned);'."\n",
-'const ex &','');
+'const ex &','','');
$typedef_series_funcp=generate(
'typedef ex (* series_funcp_${N})(${SEQ1}, const relational &, int, unsigned);'."\n",
-'const ex &','');
+'const ex &','','');
-$eval_func_interface=generate(' function_options & eval_func(eval_funcp_${N} e);'."\n",'','');
+$eval_func_interface=generate(' function_options & eval_func(eval_funcp_${N} e);'."\n",'','','');
-$evalf_func_interface=generate(' function_options & evalf_func(evalf_funcp_${N} ef);'."\n",'','');
+$evalf_func_interface=generate(' function_options & evalf_func(evalf_funcp_${N} ef);'."\n",'','','');
-$derivative_func_interface=generate(' function_options & derivative_func(derivative_funcp_${N} d);'."\n",'','');
+$derivative_func_interface=generate(' function_options & derivative_func(derivative_funcp_${N} d);'."\n",'','','');
-$series_func_interface=generate(' function_options & series_func(series_funcp_${N} s);'."\n",'','');
+$series_func_interface=generate(' function_options & series_func(series_funcp_${N} s);'."\n",'','','');
$constructors_interface=generate(
' function(unsigned ser, ${SEQ1});'."\n",
-'const ex & param${N}','');
+'const ex & param${N}','','');
$constructors_implementation=generate(
- <<'END_OF_CONSTRUCTORS_IMPLEMENTATION','const ex & param${N}','param${N}');
+ <<'END_OF_CONSTRUCTORS_IMPLEMENTATION','const ex & param${N}','param${N}','');
function::function(unsigned ser, ${SEQ1})
: exprseq(${SEQ2}), serial(ser)
{
END_OF_CONSTRUCTORS_IMPLEMENTATION
$eval_switch_statement=generate(
- <<'END_OF_EVAL_SWITCH_STATEMENT','seq[${N}-1]','');
+ <<'END_OF_EVAL_SWITCH_STATEMENT','seq[${N}-1]','','');
case ${N}:
eval_result = ((eval_funcp_${N})(registered_functions()[serial].eval_f))(${SEQ1});
break;
END_OF_EVAL_SWITCH_STATEMENT
$evalf_switch_statement=generate(
- <<'END_OF_EVALF_SWITCH_STATEMENT','eseq[${N}-1]','');
+ <<'END_OF_EVALF_SWITCH_STATEMENT','eseq[${N}-1]','','');
case ${N}:
return ((evalf_funcp_${N})(registered_functions()[serial].evalf_f))(${SEQ1});
END_OF_EVALF_SWITCH_STATEMENT
$diff_switch_statement=generate(
- <<'END_OF_DIFF_SWITCH_STATEMENT','seq[${N}-1]','');
+ <<'END_OF_DIFF_SWITCH_STATEMENT','seq[${N}-1]','','');
case ${N}:
return ((derivative_funcp_${N})(registered_functions()[serial].derivative_f))(${SEQ1},diff_param);
END_OF_DIFF_SWITCH_STATEMENT
$series_switch_statement=generate(
- <<'END_OF_SERIES_SWITCH_STATEMENT','seq[${N}-1]','');
+ <<'END_OF_SERIES_SWITCH_STATEMENT','seq[${N}-1]','','');
case ${N}:
try {
res = ((series_funcp_${N})(registered_functions()[serial].series_f))(${SEQ1},r,order,options);
END_OF_SERIES_SWITCH_STATEMENT
$eval_func_implementation=generate(
- <<'END_OF_EVAL_FUNC_IMPLEMENTATION','','');
+ <<'END_OF_EVAL_FUNC_IMPLEMENTATION','','','');
function_options & function_options::eval_func(eval_funcp_${N} e)
{
test_and_set_nparams(${N});
END_OF_EVAL_FUNC_IMPLEMENTATION
$evalf_func_implementation=generate(
- <<'END_OF_EVALF_FUNC_IMPLEMENTATION','','');
+ <<'END_OF_EVALF_FUNC_IMPLEMENTATION','','','');
function_options & function_options::evalf_func(evalf_funcp_${N} ef)
{
test_and_set_nparams(${N});
END_OF_EVALF_FUNC_IMPLEMENTATION
$derivative_func_implementation=generate(
- <<'END_OF_DERIVATIVE_FUNC_IMPLEMENTATION','','');
+ <<'END_OF_DERIVATIVE_FUNC_IMPLEMENTATION','','','');
function_options & function_options::derivative_func(derivative_funcp_${N} d)
{
test_and_set_nparams(${N});
END_OF_DERIVATIVE_FUNC_IMPLEMENTATION
$series_func_implementation=generate(
- <<'END_OF_SERIES_FUNC_IMPLEMENTATION','','');
+ <<'END_OF_SERIES_FUNC_IMPLEMENTATION','','','');
function_options & function_options::series_func(series_funcp_${N} s)
{
test_and_set_nparams(${N});
// usually executed before main(), so the exception could not
// caught anyhow
std::cerr << "WARNING: function name " << opt.name
- << " already in use!" << std::endl;
+ << " already in use!" << std::endl;
}
registered_functions().push_back(opt);
if (opt.use_remember) {
remember_table::remember_tables().
push_back(remember_table(opt.remember_size,
- opt.remember_assoc_size,
- opt.remember_strategy));
+ opt.remember_assoc_size,
+ opt.remember_strategy));
} else {
remember_table::remember_tables().push_back(remember_table());
}
// method:
// construct series manually in a dummy symbol s
const symbol s;
- ex ser = zeta(2);
+ ex ser = zeta(_ex2());
// manually construct the primitive expansion
for (int i=1; i<order; ++i)
ser += pow(1-s,i) * (numeric(1,i)*(I*Pi+log(s-1)) - numeric(1,i*i));
// overloading at work: we cannot use the macros here
/** Riemann's Zeta-function. */
extern const unsigned function_index_zeta1;
-inline function zeta(const ex & p1) {
- return function(function_index_zeta1, p1);
+template<typename T1>
+inline function zeta(const T1 & p1) {
+ return function(function_index_zeta1, ex(p1));
}
/** Derivatives of Riemann's Zeta-function. */
extern const unsigned function_index_zeta2;
-inline function zeta(const ex & p1, const ex & p2) {
- return function(function_index_zeta2, p1, p2);
+template<typename T1, typename T2>
+inline function zeta(const T1 & p1, const T2 & p2) {
+ return function(function_index_zeta2, ex(p1), ex(p2));
}
/** Gamma-function. */
// overloading at work: we cannot use the macros here
/** Psi-function (aka digamma-function). */
extern const unsigned function_index_psi1;
-inline function psi(const ex & p1) {
- return function(function_index_psi1, p1);
+template<typename T1>
+inline function psi(const T1 & p1) {
+ return function(function_index_psi1, ex(p1));
}
/** Derivatives of Psi-function (aka polygamma-functions). */
extern const unsigned function_index_psi2;
-inline function psi(const ex & p1, const ex & p2) {
- return function(function_index_psi2, p1, p2);
+template<typename T1, typename T2>
+inline function psi(const T1 & p1, const T2 & p2) {
+ return function(function_index_psi2, ex(p1), ex(p2));
}
/** Factorial function. */
// if we got here we have to care for a simple pole of tgamma(-m):
numeric m = -ex_to<numeric>(arg_pt);
ex recur;
- for (numeric p; p<=m; ++p)
+ for (numeric p = 0; p<=m; ++p)
recur += log(arg+p);
return (lgamma(arg+m+_ex1())-recur).series(rel, order, options);
}
{
if (x.info(info_flags::numeric)) {
// trap integer arguments:
- if (x.info(info_flags::integer)) {
+ const numeric two_x = _num2()*ex_to<numeric>(x);
+ if (two_x.is_even()) {
// tgamma(n) -> (n-1)! for postitive n
- if (x.info(info_flags::posint)) {
+ if (two_x.is_positive()) {
return factorial(ex_to<numeric>(x).sub(_num1()));
} else {
throw (pole_error("tgamma_eval(): simple pole",1));
}
}
// trap half integer arguments:
- if ((x*2).info(info_flags::integer)) {
+ if (two_x.is_integer()) {
// trap positive x==(n+1/2)
// tgamma(n+1/2) -> Pi^(1/2)*(1*3*..*(2*n-1))/(2^n)
- if ((x*_ex2()).info(info_flags::posint)) {
- numeric n = ex_to<numeric>(x).sub(_num1_2());
- numeric coefficient = doublefactorial(n.mul(_num2()).sub(_num1()));
- coefficient = coefficient.div(pow(_num2(),n));
- return coefficient * pow(Pi,_ex1_2());
+ if (two_x.is_positive()) {
+ const numeric n = ex_to<numeric>(x).sub(_num1_2());
+ return (doublefactorial(n.mul(_num2()).sub(_num1())).div(pow(_num2(),n))) * pow(Pi,_ex1_2());
} else {
// trap negative x==(-n+1/2)
// tgamma(-n+1/2) -> Pi^(1/2)*(-2)^n/(1*3*..*(2*n-1))
- numeric n = abs(ex_to<numeric>(x).sub(_num1_2()));
- numeric coefficient = pow(_num_2(), n);
- coefficient = coefficient.div(doublefactorial(n.mul(_num2()).sub(_num1())));;
- return coefficient*power(Pi,_ex1_2());
+ const numeric n = abs(ex_to<numeric>(x).sub(_num1_2()));
+ return (pow(_num_2(), n).div(doublefactorial(n.mul(_num2()).sub(_num1()))))*power(Pi,_ex1_2());
}
}
// tgamma_evalf should be called here once it becomes available
if (!arg_pt.info(info_flags::integer) || arg_pt.info(info_flags::positive))
throw do_taylor(); // caught by function::series()
// if we got here we have to care for a simple pole at -m:
- numeric m = -ex_to<numeric>(arg_pt);
+ const numeric m = -ex_to<numeric>(arg_pt);
ex ser_denom = _ex1();
for (numeric p; p<=m; ++p)
ser_denom *= arg+p;
// treat all problematic x and y that may not be passed into tgamma,
// because they would throw there although beta(x,y) is well-defined
// using the formula beta(x,y) == (-1)^y * beta(1-x-y, y)
- numeric nx(ex_to<numeric>(x));
- numeric ny(ex_to<numeric>(y));
+ const numeric nx = ex_to<numeric>(x);
+ const numeric ny = ex_to<numeric>(y);
if (nx.is_real() && nx.is_integer() &&
ny.is_real() && ny.is_integer()) {
if (nx.is_negative()) {
static ex psi1_eval(const ex & x)
{
if (x.info(info_flags::numeric)) {
- numeric nx = ex_to<numeric>(x);
+ const numeric nx = ex_to<numeric>(x);
if (nx.is_integer()) {
// integer case
if (nx.is_positive()) {
// psi(n) -> 1 + 1/2 +...+ 1/(n-1) - Euler
- numeric rat(0);
- for (numeric i(nx+_num_1()); i.is_positive(); --i)
+ numeric rat = 0;
+ for (numeric i(nx+_num_1()); i>0; --i)
rat += i.inverse();
return rat-Euler;
} else {
// half integer case
if (nx.is_positive()) {
// psi((2m+1)/2) -> 2/(2m+1) + 2/2m +...+ 2/1 - Euler - 2log(2)
- numeric rat(0);
- for (numeric i((nx+_num_1())*_num2()); i.is_positive(); i-=_num2())
+ numeric rat = 0;
+ for (numeric i = (nx+_num_1())*_num2(); i>0; i-=_num2())
rat += _num2()*i.inverse();
return rat-Euler-_ex2()*log(_ex2());
} else {
// to relate psi(-m-1/2) to psi(1/2):
// psi(-m-1/2) == psi(1/2) + r
// where r == ((-1/2)^(-1) + ... + (-m-1/2)^(-1))
- numeric recur(0);
- for (numeric p(nx); p<0; ++p)
+ numeric recur = 0;
+ for (numeric p = nx; p<0; ++p)
recur -= pow(p, _num_1());
return recur+psi(_ex1_2());
}
if (!arg_pt.info(info_flags::integer) || arg_pt.info(info_flags::positive))
throw do_taylor(); // caught by function::series()
// if we got here we have to care for a simple pole at -m:
- numeric m = -ex_to<numeric>(arg_pt);
+ const numeric m = -ex_to<numeric>(arg_pt);
ex recur;
for (numeric p; p<=m; ++p)
recur += power(arg+p,_ex_1());
return log(tgamma(x));
if (n.info(info_flags::numeric) && n.info(info_flags::posint) &&
x.info(info_flags::numeric)) {
- numeric nn = ex_to<numeric>(n);
- numeric nx = ex_to<numeric>(x);
+ const numeric nn = ex_to<numeric>(n);
+ const numeric nx = ex_to<numeric>(x);
if (nx.is_integer()) {
// integer case
if (nx.is_equal(_num1()))
// to relate psi(n,m) to psi(n,1):
// psi(n,m) == psi(n,1) + r
// where r == (-)^n * n! * (1^(-n-1) + ... + (m-1)^(-n-1))
- numeric recur(0);
- for (numeric p(1); p<nx; ++p)
+ numeric recur = 0;
+ for (numeric p = 1; p<nx; ++p)
recur += pow(p, -nn+_num_1());
recur *= factorial(nn)*pow(_num_1(), nn);
return recur+psi(n,_ex1());
// use psi(n,1/2) == (-)^(n+1) * n! * (2^(n+1)-1) * zeta(n+1)
return pow(_num_1(),nn+_num1())*factorial(nn)*(pow(_num2(),nn+_num1()) + _num_1())*zeta(ex(nn+_num1()));
if (nx.is_positive()) {
- numeric m = nx - _num1_2();
+ const numeric m = nx - _num1_2();
// use the multiplication formula
// psi(n,2*m) == (psi(n,m) + psi(n,m+1/2)) / 2^(n+1)
// to revert to positive integer case
// to relate psi(n,-m-1/2) to psi(n,1/2):
// psi(n,-m-1/2) == psi(n,1/2) + r
// where r == (-)^(n+1) * n! * ((-1/2)^(-n-1) + ... + (-m-1/2)^(-n-1))
- numeric recur(0);
- for (numeric p(nx); p<0; ++p)
+ numeric recur = 0;
+ for (numeric p = nx; p<0; ++p)
recur += pow(p, -nn+_num_1());
recur *= factorial(nn)*pow(_num_1(), nn+_num_1());
return recur+psi(n,_ex1_2());
if (!arg_pt.info(info_flags::integer) || arg_pt.info(info_flags::positive))
throw do_taylor(); // caught by function::series()
// if we got here we have to care for a pole of order n+1 at -m:
- numeric m = -ex_to<numeric>(arg_pt);
+ const numeric m = -ex_to<numeric>(arg_pt);
ex recur;
for (numeric p; p<=m; ++p)
recur += power(arg+p,-n+_ex_1());
return _ex1();
}
// exp(n*Pi*I/2) -> {+1|+I|-1|-I}
- ex TwoExOverPiI=(_ex2()*x)/(Pi*I);
+ const ex TwoExOverPiI=(_ex2()*x)/(Pi*I);
if (TwoExOverPiI.info(info_flags::integer)) {
numeric z=mod(ex_to<numeric>(TwoExOverPiI),_num4());
if (z.is_equal(_num0()))
const symbol *s = static_cast<symbol *>(rel.lhs().bp);
const ex point = rel.rhs();
const symbol foo;
- ex replarg = series(log(arg), *s==foo, order).subs(foo==point);
+ const ex replarg = series(log(arg), *s==foo, order).subs(foo==point);
epvector seq;
seq.push_back(expair(-I*csgn(arg*I)*Pi, _ex0()));
seq.push_back(expair(Order(_ex1()), order));
static ex sin_eval(const ex & x)
{
// sin(n/d*Pi) -> { all known non-nested radicals }
- ex SixtyExOverPi = _ex60()*x/Pi;
+ const ex SixtyExOverPi = _ex60()*x/Pi;
ex sign = _ex1();
if (SixtyExOverPi.info(info_flags::integer)) {
numeric z = mod(ex_to<numeric>(SixtyExOverPi),_num120());
static ex cos_eval(const ex & x)
{
// cos(n/d*Pi) -> { all known non-nested radicals }
- ex SixtyExOverPi = _ex60()*x/Pi;
+ const ex SixtyExOverPi = _ex60()*x/Pi;
ex sign = _ex1();
if (SixtyExOverPi.info(info_flags::integer)) {
numeric z = mod(ex_to<numeric>(SixtyExOverPi),_num120());
static ex tan_eval(const ex & x)
{
// tan(n/d*Pi) -> { all known non-nested radicals }
- ex SixtyExOverPi = _ex60()*x/Pi;
+ const ex SixtyExOverPi = _ex60()*x/Pi;
ex sign = _ex1();
if (SixtyExOverPi.info(info_flags::integer)) {
numeric z = mod(ex_to<numeric>(SixtyExOverPi),_num60());
const symbol *s = static_cast<symbol *>(rel.lhs().bp);
const ex point = rel.rhs();
const symbol foo;
- ex replarg = series(atan(arg), *s==foo, order).subs(foo==point);
+ const ex replarg = series(atan(arg), *s==foo, order).subs(foo==point);
ex Order0correction = replarg.op(0)+csgn(arg)*Pi*_ex_1_2();
if ((I*arg_pt)<_ex0())
Order0correction += log((I*arg_pt+_ex_1())/(I*arg_pt+_ex1()))*I*_ex_1_2();
const symbol *s = static_cast<symbol *>(rel.lhs().bp);
const ex point = rel.rhs();
const symbol foo;
- ex replarg = series(atanh(arg), *s==foo, order).subs(foo==point);
+ const ex replarg = series(atanh(arg), *s==foo, order).subs(foo==point);
ex Order0correction = replarg.op(0)+csgn(I*arg)*Pi*I*_ex1_2();
if (arg_pt<_ex0())
Order0correction += log((arg_pt+_ex_1())/(arg_pt+_ex1()))*_ex1_2();
}
if (expanded_seqp)
delete expanded_seqp;
-
+
// Now the only remaining thing to do is to multiply the factors which
// were not sums into the "last_expanded" sum
if (is_ex_exactly_of_type(last_expanded, add)) {
-
add const & finaladd = ex_to<add>(last_expanded);
+
const add & finaladd = ex_to<add>(last_expanded);
exvector distrseq;
int n = finaladd.nops();
distrseq.reserve(n);
{
// First, expand the children
exvector expanded_seq = expandchildren(options);
-
+
// Now, look for all the factors that are sums and remember their
// position and number of terms. One remark is in order here: we do not
// take into account the overall_coeff of the add objects. This is
const add & expanded_addref = ex_to<add>(*cit);
number_of_add_operands[number_of_adds] = expanded_addref.seq.size();
number_of_expanded_terms *= expanded_addref.seq.size();
- number_of_adds++;
+ ++number_of_adds;
}
- current_position++;
+ ++current_position;
}
// If there are no sums, we are done
// other ctors
public:
- explicit numeric(int i);
- explicit numeric(unsigned int i);
- explicit numeric(long i);
- explicit numeric(unsigned long i);
- explicit numeric(long numer, long denom);
- explicit numeric(double d);
- explicit numeric(const char *);
+ numeric(int i);
+ numeric(unsigned int i);
+ numeric(long i);
+ numeric(unsigned long i);
+ numeric(long numer, long denom);
+ numeric(double d);
+ numeric(const char *);
// functions overriding virtual functions from bases classes
public:
// other ctors
//////////
-power::power(const ex & lh, const ex & rh) : inherited(TINFO_power), basis(lh), exponent(rh)
-{
- debugmsg("power ctor from ex,ex",LOGLEVEL_CONSTRUCT);
-}
-
-/** Ctor from an ex and a bare numeric. This is somewhat more efficient than
- * the normal ctor from two ex whenever it can be used. */
-power::power(const ex & lh, const numeric & rh) : inherited(TINFO_power), basis(lh), exponent(rh)
-{
- debugmsg("power ctor from ex,numeric",LOGLEVEL_CONSTRUCT);
-}
+// all inlined
//////////
// archiving
ex power::evalm(void) const
{
- ex ebasis = basis.evalm();
- ex eexponent = exponent.evalm();
+ const ex ebasis = basis.evalm();
+ const ex eexponent = exponent.evalm();
if (is_ex_of_type(ebasis,matrix)) {
if (is_ex_of_type(eexponent,numeric)) {
return (new matrix(ex_to<matrix>(ebasis).pow(eexponent)))->setflag(status_flags::dynallocated);
upper_limit[l] = n;
}
- while (1) {
+ while (true) {
exvector term;
term.reserve(m+1);
for (l=0; l<m-1; l++) {
term.push_back(f);
- /*
- cout << "begin term" << endl;
- for (int i=0; i<m-1; i++) {
- cout << "k[" << i << "]=" << k[i] << endl;
- cout << "k_cum[" << i << "]=" << k_cum[i] << endl;
- cout << "upper_limit[" << i << "]=" << upper_limit[i] << endl;
- }
- for_each(term.begin(), term.end(), ostream_iterator<ex>(cout, "\n"));
- cout << "end term" << endl;
- */
-
- // TODO: optimize this
+ // TODO: Can we optimize this? Alex seemed to think so...
sum.push_back((new mul(term))->setflag(status_flags::dynallocated));
// increment k[]
l = m-2;
- while ((l>=0)&&((++k[l])>upper_limit[l])) {
+ while ((l>=0) && ((++k[l])>upper_limit[l])) {
k[l] = 0;
--l;
}
upper_limit[i] = n-k_cum[i-1];
}
return (new add(sum))->setflag(status_flags::dynallocated |
- status_flags::expanded );
+ status_flags::expanded );
}
// power(+(x,...,z;c),2)=power(+(x,...,z;0),2)+2*c*+(x,...,z;0)+c*c
// first part: ignore overall_coeff and expand other terms
for (epvector::const_iterator cit0=a.seq.begin(); cit0!=last; ++cit0) {
- const ex & r = (*cit0).rest;
- const ex & c = (*cit0).coeff;
+ const ex & r = cit0->rest;
+ const ex & c = cit0->coeff;
GINAC_ASSERT(!is_ex_exactly_of_type(r,add));
GINAC_ASSERT(!is_ex_exactly_of_type(r,power) ||
}
for (epvector::const_iterator cit1=cit0+1; cit1!=last; ++cit1) {
- const ex & r1 = (*cit1).rest;
- const ex & c1 = (*cit1).coeff;
+ const ex & r1 = cit1->rest;
+ const ex & c1 = cit1->coeff;
sum.push_back(a.combine_ex_with_coeff_to_pair((new mul(r,r1))->setflag(status_flags::dynallocated),
_num2().mul(ex_to<numeric>(c)).mul_dyn(ex_to<numeric>(c1))));
}
return (new mul(distrseq,ex_to<numeric>(m.overall_coeff).power_dyn(n)))->setflag(status_flags::dynallocated);
}
-/*
-ex power::expand_noncommutative(const ex & basis, const numeric & exponent,
- unsigned options) const
-{
- ex rest_power = ex(power(basis,exponent.add(_num_1()))).
- expand(options | expand_options::internal_do_not_expand_power_operands);
-
- return ex(mul(rest_power,basis),0).
- expand(options | expand_options::internal_do_not_expand_mul_operands);
-}
-*/
-
// helper function
ex sqrt(const ex & a)
// other ctors
public:
- power(const ex & lh, const ex & rh);
- power(const ex & lh, const numeric & rh);
+ power(const ex & lh, const ex & rh) : inherited(TINFO_power), basis(lh), exponent(rh) {}
+ template<typename T> power(const ex & lh, const T & rh) : inherited(TINFO_power), basis(lh), exponent(rh) {}
// functions overriding virtual functions from bases classes
public:
ex expand_add(const add & a, int n) const;
ex expand_add_2(const add & a) const;
ex expand_mul(const mul & m, const numeric & n) const;
- //ex expand_commutative_3(const ex & basis, const numeric & exponent,
- // unsigned options) const;
- //ex expand_noncommutative(const ex & basis, const numeric & exponent, unsigned options) const;
// member variables
{
return power(b, e);
}
+template<typename T1, typename T2>
+inline ex pow(const T1 & b, const T2 & e)
+{
+ return power(ex(b), ex(e));
+}
/** Square root expression. Returns a power-object with exponent 1/2. */
ex sqrt(const ex & a);
git://www.ginac.de / ginac.git / commitdiff