Advanced Computing Platform for Theoretical Physics

Commit 79b3749b authored by Feng Feng's avatar Feng Feng
Browse files

fix a bug in fermat_numer_denom

parent 247c84d3
......@@ -1873,6 +1873,15 @@ namespace HepLib {
auto ostr = fermat.Execute(ss.str());
ss.clear();
ss.str("");
// note the order, before exfactor (fermat_normal will be called again here)
ss << "&(U=0);" << endl; // disable ugly printing
if(fermat_use_array) ss << "@(res,[m]);" << endl;
else ss << "@(res,item);" << endl;
ss << "&_G;" << endl;
fermat.Execute(ss.str());
ss.clear();
ss.str("");
// make sure last char is 0
if(ostr[ostr.length()-1]!='0') throw Error("fermat_together: last char is NOT 0.");
......@@ -1883,7 +1892,7 @@ namespace HepLib {
cpos = ostr.find(estr);
if(cpos==string::npos) throw Error(estr+" NOT Found.");
ostr = ostr.substr(0,cpos);
string_trim(ostr);
string_trim(ostr);
symtab st;
Parser fp(st);
......@@ -1891,14 +1900,6 @@ namespace HepLib {
num *= ret.op(0);
if(factor) den *= exfactor(ret.op(1));
else den *= ret.op(1);
ss << "&(U=0);" << endl; // disable ugly printing
if(fermat_use_array) ss << "@(res,[m]);" << endl;
else ss << "@(res,item);" << endl;
ss << "&_G;" << endl;
fermat.Execute(ss.str());
ss.clear();
ss.str("");
}
//fermat.Exit();
......
......@@ -100,6 +100,26 @@ namespace HepLib::FC {
})(ret);
return ret;
}
/**
* @brief convert F(ps, ns) to normal ex, ns is like FIRE convention
* @param expr_in expression contains F
* @return F(ps, ns) converted into normal expression wrapped in F
*/
ex F2F(const ex & expr_in) {
ex ret = expr_in;
ret = MapFunction([](const ex & e, MapFunction &self)->ex{
if(e.match(F(w1, w2))) {
auto ps = e.op(0);
auto ns = e.op(1);
ex res = 1;
for(int i=0; i<ps.nops(); i++) res *= pow(ps.op(i), ex(0)-ns.op(i));
return F(res);
} else if(!e.has((F(w1,w2)))) return e;
else return e.map(self);
})(ret);
return ret;
}
/**
* @brief Apart on matrix
......@@ -785,7 +805,7 @@ namespace HepLib::FC {
if(Verbose>0) cout << " \\--Refined Ints/Pros: " << nints << "/" << ibp_vec_re.size() << " @ " << now(false) << endl;
MapFunction F2ex([ibp_vec](const ex &e, MapFunction &self)->ex {
MapFunction _F2ex([ibp_vec](const ex &e, MapFunction &self)->ex {
if(!e.has(F(w1,w2))) return e;
else if(e.match(F(w1,w2))) {
int idx = ex_to<numeric>(e.op(0)).to_int();
......@@ -799,7 +819,7 @@ namespace HepLib::FC {
auto air = air_vec[idx];
air = subs_naive(air,IR2F);
air = subs_naive(air,rules_ints.first);
air = F2ex(air);
air = _F2ex(air);
auto cv_lst = mma_collect_lst(air, F(w1,w2));
air = 0;
for(auto cv : cv_lst) air += cv.op(1) * fermat_normal(cv.op(0));
......@@ -867,7 +887,7 @@ namespace HepLib::FC {
air = subs_naive(air,rules_ints.first);
air = subs_naive(air,mi_rules.first);
air = subs_naive(air,F2F);
air = F2ex(air);
air = _F2ex(air);
auto cv_lst = mma_collect_lst(air, F(w1,w2));
air = 0;
for(auto cv : cv_lst) air += cv.op(1) * fermat_normal(cv.op(0));
......
......@@ -316,8 +316,9 @@ namespace HepLib::FC {
ex ApartIR2ex(const ex & expr_in);
ex ApartIR2F(const ex & expr_in);
ex F2ex(const ex & expr_in);
ex F2F(const ex & expr_in);
ex ApartIRC(const ex & expr_in, const ex & cut_props=lst{});
void ApartIBP(int IBPmethod, exvector &air_vec, const lst & loops_exts=lst{}, const lst & cut_props=lst{}, std::function<lst(const Base &, const ex &)> uf=IBP::LoopUF);
void ApartIBP(int IBPmethod, exvector &io_vec, const lst & loops_exts=lst{}, const lst & cut_props=lst{}, std::function<lst(const Base &, const ex &)> uf=IBP::LoopUF);
/**
* @brief ApartIR function with 1 argument
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment