X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=blobdiff_plain;f=ginac%2Fmul.cpp;h=dd96a70fafdb15114985a1dce12d6d4aeba3569b;hp=ab46ae0e91c8a589a8599ca3300de76352636eed;hb=04aace56d3d2ca82be59460bf1922600ae43449a;hpb=1be6e8f6ed13bcfb278d13e01f1b4e4851a91114 diff --git a/ginac/mul.cpp b/ginac/mul.cpp index ab46ae0e..dd96a70f 100644 --- a/ginac/mul.cpp +++ b/ginac/mul.cpp @@ -3,7 +3,7 @@ * Implementation of GiNaC's products of expressions. */ /* - * GiNaC Copyright (C) 1999-2005 Johannes Gutenberg University Mainz, Germany + * GiNaC Copyright (C) 1999-2007 Johannes Gutenberg University Mainz, Germany * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by @@ -51,7 +51,7 @@ GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(mul, expairseq, mul::mul() { - tinfo_key = TINFO_mul; + tinfo_key = &mul::tinfo_static; } ////////// @@ -62,7 +62,7 @@ mul::mul() mul::mul(const ex & lh, const ex & rh) { - tinfo_key = TINFO_mul; + tinfo_key = &mul::tinfo_static; overall_coeff = _ex1; construct_from_2_ex(lh,rh); GINAC_ASSERT(is_canonical()); @@ -70,7 +70,7 @@ mul::mul(const ex & lh, const ex & rh) mul::mul(const exvector & v) { - tinfo_key = TINFO_mul; + tinfo_key = &mul::tinfo_static; overall_coeff = _ex1; construct_from_exvector(v); GINAC_ASSERT(is_canonical()); @@ -78,32 +78,32 @@ mul::mul(const exvector & v) mul::mul(const epvector & v) { - tinfo_key = TINFO_mul; + tinfo_key = &mul::tinfo_static; overall_coeff = _ex1; construct_from_epvector(v); GINAC_ASSERT(is_canonical()); } -mul::mul(const epvector & v, const ex & oc) +mul::mul(const epvector & v, const ex & oc, bool do_index_renaming) { - tinfo_key = TINFO_mul; + tinfo_key = &mul::tinfo_static; overall_coeff = oc; - construct_from_epvector(v); + construct_from_epvector(v, do_index_renaming); GINAC_ASSERT(is_canonical()); } -mul::mul(std::auto_ptr vp, const ex & oc) +mul::mul(std::auto_ptr vp, const ex & oc, bool do_index_renaming) { - tinfo_key = TINFO_mul; + tinfo_key = &mul::tinfo_static; GINAC_ASSERT(vp.get()!=0); overall_coeff = oc; - construct_from_epvector(*vp); + construct_from_epvector(*vp, do_index_renaming); GINAC_ASSERT(is_canonical()); } mul::mul(const ex & lh, const ex & mh, const ex & rh) { - tinfo_key = TINFO_mul; + tinfo_key = &mul::tinfo_static; exvector factors; factors.reserve(3); factors.push_back(lh); @@ -218,8 +218,12 @@ void mul::do_print_csrc(const print_csrc & c, unsigned level) const c.s << "("; if (!overall_coeff.is_equal(_ex1)) { - overall_coeff.print(c, precedence()); - c.s << "*"; + if (overall_coeff.is_equal(_ex_1)) + c.s << "-"; + else { + overall_coeff.print(c, precedence()); + c.s << "*"; + } } // Print arguments, separated by "*" or "/" @@ -461,6 +465,41 @@ ex mul::evalf(int level) const return mul(s, overall_coeff.evalf(level)); } +void mul::find_real_imag(ex & rp, ex & ip) const +{ + rp = overall_coeff.real_part(); + ip = overall_coeff.imag_part(); + for (epvector::const_iterator i=seq.begin(); i!=seq.end(); ++i) { + ex factor = recombine_pair_to_ex(*i); + ex new_rp = factor.real_part(); + ex new_ip = factor.imag_part(); + if(new_ip.is_zero()) { + rp *= new_rp; + ip *= new_rp; + } else { + ex temp = rp*new_rp - ip*new_ip; + ip = ip*new_rp + rp*new_ip; + rp = temp; + } + } + rp = rp.expand(); + ip = ip.expand(); +} + +ex mul::real_part() const +{ + ex rp, ip; + find_real_imag(rp, ip); + return rp; +} + +ex mul::imag_part() const +{ + ex rp, ip; + find_real_imag(rp, ip); + return ip; +} + ex mul::evalm() const { // numeric*matrix @@ -559,7 +598,7 @@ bool tryfactsubs(const ex & origfactor, const ex & patternfactor, int & nummatch return true; } -/** Checks wheter e matches to the pattern pat and the (possibly to be updated +/** Checks wheter e matches to the pattern pat and the (possibly to be updated) * list of replacements repls. This matching is in the sense of algebraic * substitutions. Matching starts with pat.op(factor) of the pattern because * the factors before this one have already been matched. The (possibly @@ -595,10 +634,28 @@ bool algebraic_match_mul_with_mul(const mul &e, const ex &pat, lst &repls, return false; } +bool mul::has(const ex & pattern, unsigned options) const +{ + if(!(options&has_options::algebraic)) + return basic::has(pattern,options); + if(is_a(pattern)) { + lst repls; + int nummatches = std::numeric_limits::max(); + std::vector subsed(seq.size(), false); + std::vector matched(seq.size(), false); + if(algebraic_match_mul_with_mul(*this, pattern, repls, 0, nummatches, + subsed, matched)) + return true; + } + return basic::has(pattern, options); +} + ex mul::algebraic_subs_mul(const exmap & m, unsigned options) const { std::vector subsed(seq.size(), false); exvector subsresult(seq.size()); + ex divide_by = 1; + ex multiply_by = 1; for (exmap::const_iterator it = m.begin(); it != m.end(); ++it) { @@ -606,36 +663,35 @@ ex mul::algebraic_subs_mul(const exmap & m, unsigned options) const retry1: int nummatches = std::numeric_limits::max(); std::vector currsubsed(seq.size(), false); - bool succeed = true; lst repls; if(!algebraic_match_mul_with_mul(*this, it->first, repls, 0, nummatches, subsed, currsubsed)) continue; - bool foundfirstsubsedfactor = false; - for (size_t j=0; jsecond.subs(ex(repls), subs_options::no_pattern) / it->first.subs(ex(repls), subs_options::no_pattern), nummatches); - } + for (size_t j=0; jfirst.subs(ex(repls), subs_options::no_pattern); + divide_by *= power(subsed_pattern, nummatches); + ex subsed_result + = it->second.subs(ex(repls), subs_options::no_pattern); + multiply_by *= power(subsed_result, nummatches); goto retry1; } else { -retry2: - int nummatches = std::numeric_limits::max(); - lst repls; for (size_t j=0; jnops(); j++) { - if (!subsed[j] && tryfactsubs(op(j), it->first, nummatches, repls)) { + int nummatches = std::numeric_limits::max(); + lst repls; + if (!subsed[j] && tryfactsubs(op(j), it->first, nummatches, repls)){ subsed[j] = true; - subsresult[j] = op(j) * power(it->second.subs(ex(repls), subs_options::no_pattern) / it->first.subs(ex(repls), subs_options::no_pattern), nummatches); - goto retry2; + ex subsed_pattern + = it->first.subs(ex(repls), subs_options::no_pattern); + divide_by *= power(subsed_pattern, nummatches); + ex subsed_result + = it->second.subs(ex(repls), subs_options::no_pattern); + multiply_by *= power(subsed_result, nummatches); } } } @@ -651,15 +707,7 @@ retry2: if (!subsfound) return subs_one_level(m, options | subs_options::algebraic); - exvector ev; ev.reserve(nops()); - for (size_t i=0; isetflag(status_flags::dynallocated); + return ((*this)/divide_by)*multiply_by; } // protected @@ -715,8 +763,8 @@ unsigned mul::return_type() const if ((rt == return_types::noncommutative) && (!all_commutative)) { // another nc element found, compare type_infos if (noncommutative_element->rest.return_type_tinfo() != i->rest.return_type_tinfo()) { - // diffent types -> mul is ncc - return return_types::noncommutative_composite; + // different types -> mul is ncc + return return_types::noncommutative_composite; } } ++i; @@ -725,10 +773,10 @@ unsigned mul::return_type() const return all_commutative ? return_types::commutative : return_types::noncommutative; } -unsigned mul::return_type_tinfo() const +tinfo_t mul::return_type_tinfo() const { if (seq.empty()) - return tinfo_key; // mul without factors: should not happen + return this; // mul without factors: should not happen // return type_info of first noncommutative element epvector::const_iterator i = seq.begin(), end = seq.end(); @@ -738,17 +786,17 @@ unsigned mul::return_type_tinfo() const ++i; } // no noncommutative element found, should not happen - return tinfo_key; + return this; } -ex mul::thisexpairseq(const epvector & v, const ex & oc) const +ex mul::thisexpairseq(const epvector & v, const ex & oc, bool do_index_renaming) const { - return (new mul(v, oc))->setflag(status_flags::dynallocated); + return (new mul(v, oc, do_index_renaming))->setflag(status_flags::dynallocated); } -ex mul::thisexpairseq(std::auto_ptr vp, const ex & oc) const +ex mul::thisexpairseq(std::auto_ptr vp, const ex & oc, bool do_index_renaming) const { - return (new mul(vp, oc))->setflag(status_flags::dynallocated); + return (new mul(vp, oc, do_index_renaming))->setflag(status_flags::dynallocated); } expair mul::split_ex_to_pair(const ex & e) const @@ -920,11 +968,11 @@ ex mul::expand(unsigned options) const exvector add1_dummy_indices, add2_dummy_indices, add_indices; for (epvector::const_iterator i=add1begin; i!=add1end; ++i) { - add_indices = get_all_dummy_indices(i->rest); + add_indices = get_all_dummy_indices_safely(i->rest); add1_dummy_indices.insert(add1_dummy_indices.end(), add_indices.begin(), add_indices.end()); } for (epvector::const_iterator i=add2begin; i!=add2end; ++i) { - add_indices = get_all_dummy_indices(i->rest); + add_indices = get_all_dummy_indices_safely(i->rest); add2_dummy_indices.insert(add2_dummy_indices.end(), add_indices.begin(), add_indices.end()); } @@ -971,7 +1019,7 @@ ex mul::expand(unsigned options) const size_t n = last_expanded.nops(); exvector distrseq; distrseq.reserve(n); - exvector va = get_all_dummy_indices(mul(non_adds)); + exvector va = get_all_dummy_indices_safely(mul(non_adds)); sort(va.begin(), va.end(), ex_is_less()); for (size_t i=0; i