Print["MBnum v.0.1, last modified: 18.06.09"]; MBnum::usage = "MBnum[repr_, epsilon_, points_, powers_, loop_]" Options[MBnum]={ RulesConstrain->60, ShowMBrep->False, ShowPositions->True, ShowRules->False, ShowComments->True, ShowMBcontinue->False, ShowMBexpand->False, ShowMBintegrate->False, ShowETArules->False, ShowETAcontinue1->False, ShowETAexpand1->False, ShowETAcontinue2->False, ShowETAexpand2->False, ShowMBintegrateETA->False, Analytical->False, Numerical->True }; OptChange[fun_,new___Rule]:= Module[{i}, Map[If[Length[i=Position[{new},First[#]->_,1,1]]===0, #,{new}[[Sequence @@ First[i]]]]&, Options[fun]] ] MBnum[repr_, epsilon_, points_, powers_, loop_,options___Rule] := Module[{option,AllRules, positions, cond, NoEtaRepr, NoEtaRules, NoEtaInts, NoEtaSer, num1list, PowersEta, ReprEta, RulesEta, num2list, AnalyticalResult,NumericalResult}, option=OptChange[MBnum,options]; (* after= BarnesLemma[rep,1]; *) If[ShowMBrep /.option, WriteString[$Output,"repr=",InputForm[repr], "\n","Length=",Length[repr],"\n\n"]]; Rules=Map[fun[#,eps->0,{},{eps}]&, repr /.powers]; AllRules = TimeConstrained[Rules /.fun->MBoptimizedRules,RulesConstrain /.option,Rules /.fun->MBrules]; positions = Position[AllRules, {}, 1]; If[ShowPositions /.option, Print["ETA's will be aplied on positions: ",Flatten[positions]]]; If[ShowRules /.option, Print["List of all rules:\n",AllRules]]; (*---------------------*) (*----- 1) NO ETA -----*) (*---------------------*) If[ShowComments /.option, Print["1. Calculating 'no eta' parts..."]]; NoEtaRepr = Delete[repr, positions] /.powers; If[Flatten[NoEtaRepr]==={}, {If[ShowComments /.option, Print[" No 'no eta' parts found!!!"]]; NoEtaSer={}; num1list={}}, {NoEtaRules= Delete[AllRules, positions]; If[ShowComments /.option, Print[" Running MBcontinue..."]]; NoEtaInts = MapThread[MBcontinue[#1,eps->0,#2,Verbose->False]&, {NoEtaRepr,NoEtaRules}]; If[ShowMBcontinue /.option, Print[NoEtaInts]]; If[ShowComments /.option, Print[" Running MBexpand..."]]; NoEtaSer = Map[MBmerge[MBexpand[{#},Exp[loop*eps*EulerGamma],{eps,0,epsilon}]]&, NoEtaInts]; NoEtaSer = MBmerge[NoEtaSer]; If[ShowMBexpand /.option, Print[NoEtaSer]]; If[Numerical/.option, If[ShowComments /.option, Print[" Running MBintegrate..."]]; num1list = Map[MBintegrate[{#},points,Verbose->False]&, DeleteCases[NoEtaSer,{}]]; num1list = First[#]& /@ num1list; If[ShowMBintegrate /.option, Print[num1list]]] }]; (*-----------------------*) (*----- 2) WITH ETA -----*) (*-----------------------*) Label["eta"]; If[ShowComments /.option, Print["2. Calculating 'eta' parts..."]]; If[positions==={}, {If[ShowComments /.option, Print[" No 'eta' parts found!!!"]]; EtaExp={}; num2list={}}, {Do[ PowersEta = ReplacePart[powers, powers[[i]] /.Rule[n_,x_]->Rule[n,x+eta], i]; If[ShowComments /.option, Print[" Using eta nr: ",i, "\n ",PowersEta]]; ReprEta = Map[Part[repr,#] /.PowersEta &, Flatten[positions]]; RulesEta = Map[MBoptimizedRules[#,eps->0,{},{eps,eta}]&, ReprEta]; If[ShowComments /.option, Print[" Calculating eta rules on nr: ",i]]; If[ShowETArules /.option, Print["List of eta rules:\n",RulesEta]]; If[ShowComments /.option, Print[" Running MBcontinue1 for eta nr: ",i]]; EtaInts = MapThread[MBcontinue[#1,eta->0,#2,Verbose->False]&, {ReprEta,RulesEta}]; cond=Union[Map[MemberQ[#, MBint[___]]&,EtaInts]]//First; If[cond, {Break[]}, {If[ShowComments /.option, Print[" Couldn't find rules, checking next eta..."]]}]; If[i==Length[powers], {If[ShowComments /.option, Print[" Couldn't find rules for etas!!!"]]; Goto["end"]}], {i,1,Length[powers]}]; If[ShowETAcontinue1 /.option, Print[EtaInts]]; If[ShowComments /.option, Print[" Running MBexpand1 for eta..."]]; EtaSer = MBexpand[#,1,{eta,0,0}]& /@ EtaInts; If[ShowETAexpand1 /.option, Print[EtaSer]]; If[ShowComments /.option, Print[" Running MBcontinue2 for eta..."]]; EtaCont = EtaSer /.MBint[integrand_,rules_]:> MBcontinue[integrand,eps->0,rules,Verbose->False]; If[ShowETAcontinue2 /.option, Print[EtaCont]]; If[ShowComments /.option, Print[" Running MBexpand2 for eta..."]]; EtaExp = MBmerge[MBexpand[#,Exp[loop*eps*EulerGamma],{eps,0,epsilon}]]& /@ EtaCont; EtaExp = MBmerge[EtaExp]; If[ShowETAexpand2 /.option, Print[EtaExp]]; If[Numerical/.option, If[ShowComments /.option, Print[" Running MBintegrate eta..."]]; num2list = Map[MBintegrate[{#},points,Verbose->False]&, EtaExp]; num2list = First[#]& /@ num2list; If[ShowMBintegrateETA /.option, Print[num2list]]] }]; Label["end"]; AnalyticalResult = Join[NoEtaSer,EtaExp]; If[Numerical/.option,NumericalResult = Plus @@ Join[num1list,num2list]]; If[((Analytical /.option) && (Numerical /.option)), {result = Join[{NumericalResult},{AnalyticalResult}]}, {If[(Analytical /.option), result = AnalyticalResult, result = NumericalResult]; }]; Return[result]]; (* 3.06, ShowMBrep option added 18.06, Numerical->False improved, adding If *)