# File: GradeAntiderivative.mpl # Original version thanks to Albert Rich emailed on 03/21/2017 #Nasser 03/22/2017 Use Maple leaf count instead since buildin #Nasser 03/23/2017 missing 'ln' for ElementaryFunctionQ added #Nasser 03/24/2017 corrected the check for complex result #Nasser 10/27/2017 check for leafsize and do not call ExpnType() # if leaf size is "too large". Set at 500,000 #Nasser 12/22/2019 Added debug flag, added 'dilog' to special functions # see problem 156, file Apostol_Problems GradeAntiderivative := proc(result,optimal) local leaf_count_result, leaf_count_optimal,ExpnType_result,ExpnType_optimal,debug:=false; leaf_count_result:=leafcount(result); #do NOT call ExpnType() if leaf size is too large. Recursion problem if leaf_count_result > 500000 then return "B"; fi; leaf_count_optimal:=leafcount(optimal); ExpnType_result:=ExpnType(result); ExpnType_optimal:=ExpnType(optimal); if debug then print("ExpnType_result",ExpnType_result," ExpnType_optimal=",ExpnType_optimal); fi; # If result and optimal are mathematical expressions, # GradeAntiderivative[result,optimal] returns # "F" if the result fails to integrate an expression that # is integrable # "C" if result involves higher level functions than necessary # "B" if result is more than twice the size of the optimal # antiderivative # "A" if result can be considered optimal #This check below actually is not needed, since I only #call this grading only for passed integrals. i.e. I check #for "F" before calling this. But no harm of keeping it here. #just in case. if not type(result,freeof('int')) then return "F"; end if; if ExpnType_result<=ExpnType_optimal then if debug then print("ExpnType_result<=ExpnType_optimal"); fi; if is_contains_complex(result) then if is_contains_complex(optimal) then if debug then print("both result and optimal complex"); fi; #both result and optimal complex if leaf_count_result<=2*leaf_count_optimal then return "A"; else return "B"; end if else #result contains complex but optimal is not if debug then print("result contains complex but optimal is not"); fi; return "C"; end if else # result do not contain complex # this assumes optimal do not as well if debug then print("result do not contain complex, this assumes optimal do not as well"); fi; if leaf_count_result<=2*leaf_count_optimal then if debug then print("leaf_count_result<=2*leaf_count_optimal"); fi; return "A"; else if debug then print("leaf_count_result>2*leaf_count_optimal"); fi; return "B"; end if end if else #ExpnType(result) > ExpnType(optimal) if debug then print("ExpnType(result) > ExpnType(optimal)"); fi; return "C"; end if end proc: # # is_contains_complex(result) # takes expressions and returns true if it contains "I" else false # #Nasser 032417 is_contains_complex:= proc(expression) return (has(expression,I)); end proc: # The following summarizes the type number assigned an expression # based on the functions it involves # 1 = rational function # 2 = algebraic function # 3 = elementary function # 4 = special function # 5 = hyperpergeometric function # 6 = appell function # 7 = rootsum function # 8 = integrate function # 9 = unknown function ExpnType := proc(expn) if type(expn,'atomic') then 1 elif type(expn,'list') then apply(max,map(ExpnType,expn)) elif type(expn,'sqrt') then if type(op(1,expn),'rational') then 1 else max(2,ExpnType(op(1,expn))) end if elif type(expn,'`^`') then if type(op(2,expn),'integer') then ExpnType(op(1,expn)) elif type(op(2,expn),'rational') then if type(op(1,expn),'rational') then 1 else max(2,ExpnType(op(1,expn))) end if else max(3,ExpnType(op(1,expn)),ExpnType(op(2,expn))) end if elif type(expn,'`+`') or type(expn,'`*`') then max(ExpnType(op(1,expn)),max(ExpnType(rest(expn)))) elif ElementaryFunctionQ(op(0,expn)) then max(3,ExpnType(op(1,expn))) elif SpecialFunctionQ(op(0,expn)) then max(4,apply(max,map(ExpnType,[op(expn)]))) elif HypergeometricFunctionQ(op(0,expn)) then max(5,apply(max,map(ExpnType,[op(expn)]))) elif AppellFunctionQ(op(0,expn)) then max(6,apply(max,map(ExpnType,[op(expn)]))) elif op(0,expn)='int' then max(8,apply(max,map(ExpnType,[op(expn)]))) else 9 end if end proc: ElementaryFunctionQ := proc(func) member(func,[ exp,log,ln, sin,cos,tan,cot,sec,csc, arcsin,arccos,arctan,arccot,arcsec,arccsc, sinh,cosh,tanh,coth,sech,csch, arcsinh,arccosh,arctanh,arccoth,arcsech,arccsch]) end proc: SpecialFunctionQ := proc(func) member(func,[ erf,erfc,erfi, FresnelS,FresnelC, Ei,Ei,Li,Si,Ci,Shi,Chi, GAMMA,lnGAMMA,Psi,Zeta,polylog,dilog,LambertW, EllipticF,EllipticE,EllipticPi]) end proc: HypergeometricFunctionQ := proc(func) member(func,[Hypergeometric1F1,hypergeom,HypergeometricPFQ]) end proc: AppellFunctionQ := proc(func) member(func,[AppellF1]) end proc: # u is a sum or product. rest(u) returns all but the # first term or factor of u. rest := proc(u) local v; if nops(u)=2 then op(2,u) else apply(op(0,u),op(2..nops(u),u)) end if end proc: #leafcount(u) returns the number of nodes in u. #Nasser 3/23/17 Replaced by build-in leafCount from package in Maple leafcount := proc(u) MmaTranslator[Mma][LeafCount](u); end proc: