def_cg_simp_add(e):#TODO: Improve simplification method"""Takes a sum of terms involving Clebsch-Gordan coefficients and simplifies the terms. First, we create two lists, cg_part, which is all the terms involving CG coefficients, and other_part, which is all other terms. The cg_part list is then passed to the simplification methods, which return the new cg_part and any additional terms that are added to other_part """cg_part=[]other_part=[]e=expand(e)forargine.args:ifarg.has(CG):ifisinstance(arg,Sum):other_part.append(_cg_simp_sum(arg))elifisinstance(arg,Mul):terms=1forterminarg.args:ifisinstance(term,Sum):terms*=_cg_simp_sum(term)else:terms*=termifterms.has(CG):cg_part.append(terms)else:other_part.append(terms)else:cg_part.append(arg)else:other_part.append(arg)cg_part,other=_check_varsh_871_1(cg_part)other_part.append(other)cg_part,other=_check_varsh_871_2(cg_part)other_part.append(other)cg_part,other=_check_varsh_872_9(cg_part)other_part.append(other)returnAdd(*cg_part)+Add(*other_part)def_check_varsh_871_1(term_list):# Sum( CG(a,alpha,b,0,a,alpha), (alpha, -a, a)) == KroneckerDelta(b,0)a,alpha,b,lt=map(Wild,('a','alpha','b','lt'))expr=lt*CG(a,alpha,b,0,a,alpha)simp=(2*a+1)*KroneckerDelta(b,0)sign=lt/abs(lt)build_expr=2*a+1index_expr=a+alphareturn_check_cg_simp(expr,simp,sign,lt,term_list,(a,alpha,b,lt),(a,b),build_expr,index_expr)def_check_varsh_871_2(term_list):# Sum((-1)**(a-alpha)*CG(a,alpha,a,-alpha,c,0),(alpha,-a,a))a,alpha,c,lt=map(Wild,('a','alpha','c','lt'))expr=lt*CG(a,alpha,a,-alpha,c,0)simp=sqrt(2*a+1)*KroneckerDelta(c,0)sign=(-1)**(a-alpha)*lt/abs(lt)build_expr=2*a+1index_expr=a+alphareturn_check_cg_simp(expr,simp,sign,lt,term_list,(a,alpha,c,lt),(a,c),build_expr,index_expr)def_check_varsh_872_9(term_list):# Sum( CG(a,alpha,b,beta,c,gamma)*CG(a,alpha',b,beta',c,gamma), (gamma, -c, c), (c, abs(a-b), a+b))a,alpha,alphap,b,beta,betap,c,gamma,lt=map(Wild,('a','alpha','alphap','b','beta','betap','c','gamma','lt'))# Case alpha==alphap, beta==betap# For numerical alpha,betaexpr=lt*CG(a,alpha,b,beta,c,gamma)**2simp=1sign=lt/abs(lt)x=abs(a-b)y=abs(alpha+beta)build_expr=a+b+1-Piecewise((x,x>y),(0,Eq(x,y)),(y,y>x))index_expr=a+b-cterm_list,other1=_check_cg_simp(expr,simp,sign,lt,term_list,(a,alpha,b,beta,c,gamma,lt),(a,alpha,b,beta),build_expr,index_expr)# For symbolic alpha,betax=abs(a-b)y=a+bbuild_expr=(y+1-x)*(x+y+1)index_expr=(c-x)*(x+c)+c+gammaterm_list,other2=_check_cg_simp(expr,simp,sign,lt,term_list,(a,alpha,b,beta,c,gamma,lt),(a,alpha,b,beta),build_expr,index_expr)# Case alpha!=alphap or beta!=betap# Note: this only works with leading term of 1, pattern matching is unable to match when there is a Wild leading term# For numerical alpha,alphap,beta,betapexpr=CG(a,alpha,b,beta,c,gamma)*CG(a,alphap,b,betap,c,gamma)simp=KroneckerDelta(alpha,alphap)*KroneckerDelta(beta,betap)sign=sympify(1)x=abs(a-b)y=abs(alpha+beta)build_expr=a+b+1-Piecewise((x,x>y),(0,Eq(x,y)),(y,y>x))index_expr=a+b-cterm_list,other3=_check_cg_simp(expr,simp,sign,sympify(1),term_list,(a,alpha,alphap,b,beta,betap,c,gamma),(a,alpha,alphap,b,beta,betap),build_expr,index_expr)# For symbolic alpha,alphap,beta,betapx=abs(a-b)y=a+bbuild_expr=(y+1-x)*(x+y+1)index_expr=(c-x)*(x+c)+c+gammaterm_list,other4=_check_cg_simp(expr,simp,sign,sympify(1),term_list,(a,alpha,alphap,b,beta,betap,c,gamma),(a,alpha,alphap,b,beta,betap),build_expr,index_expr)returnterm_list,other1+other2+other4def_check_cg_simp(expr,simp,sign,lt,term_list,variables,dep_variables,build_index_expr,index_expr):""" Checks for simplifications that can be made, returning a tuple of the simplified list of terms and any terms generated by simplification. Parameters ========== expr: expression The expression with Wild terms that will be matched to the terms in the sum simp: expression The expression with Wild terms that is substituted in place of the CG terms in the case of simplification sign: expression The expression with Wild terms denoting the sign that is on expr that must match lt: expression The expression with Wild terms that gives the leading term of the matched expr term_list: list A list of all of the terms is the sum to be simplified variables: list A list of all the variables that appears in expr dep_variables: list A list of the variables that must match for all the terms in the sum, i.e. the dependant variables build_index_expr: expression Expression with Wild terms giving the number of elements in cg_index index_expr: expression Expression with Wild terms giving the index terms have when storing them to cg_index """other_part=0i=0whilei<len(term_list):sub_1=_check_cg(term_list[i],expr,len(variables))ifsub_1isNone:i+=1continueifnotsympify(build_index_expr.subs(sub_1)).is_number:i+=1continuesub_dep=[(x,sub_1[x])forxindep_variables]cg_index=[None]*build_index_expr.subs(sub_1)forjinrange(i,len(term_list)):sub_2=_check_cg(term_list[j],expr.subs(sub_dep),len(variables)-len(dep_variables),sign=(sign.subs(sub_1),sign.subs(sub_dep)))ifsub_2isNone:continueifnotsympify(index_expr.subs(sub_dep).subs(sub_2)).is_number:continuecg_index[index_expr.subs(sub_dep).subs(sub_2)]=j,expr.subs(lt,1).subs(sub_dep).subs(sub_2),lt.subs(sub_2),sign.subs(sub_dep).subs(sub_2)ifall(iisnotNoneforiincg_index):min_lt=min(*[abs(term[2])fortermincg_index])indicies=[term[0]fortermincg_index]indicies.sort()indicies.reverse()[term_list.pop(i)foriinindicies]fortermincg_index:ifabs(term[2])>min_lt:term_list.append((term[2]-min_lt*term[3])*term[1])other_part+=min_lt*(sign*simp).subs(sub_1)else:i+=1returnterm_list,other_partdef_check_cg(cg_term,expr,length,sign=None):"""Checks whether a term matches the given expression"""# TODO: Check for symmetriesmatches=cg_term.match(expr)ifmatchesisNone:returnifsignisnotNone:ifnotisinstance(sign,tuple):raiseTypeError('sign must be a tuple')ifnotsign[0]==(sign[1]).subs(matches):returniflen(matches)==length:returnmatchesdef_cg_simp_sum(e):e=_check_varsh_sum_871_1(e)e=_check_varsh_sum_871_2(e)e=_check_varsh_sum_872_4(e)returnedef_check_varsh_sum_871_1(e):a=Wild('a')alpha=symbols('alpha')b=Wild('b')match=e.match(Sum(CG(a,alpha,b,0,a,alpha),(alpha,-a,a)))ifmatchisnotNoneandlen(match)==2:return((2*a+1)*KroneckerDelta(b,0)).subs(match)returnedef_check_varsh_sum_871_2(e):a=Wild('a')alpha=symbols('alpha')c=Wild('c')match=e.match(Sum((-1)**(a-alpha)*CG(a,alpha,a,-alpha,c,0),(alpha,-a,a)))ifmatchisnotNoneandlen(match)==2:return(sqrt(2*a+1)*KroneckerDelta(c,0)).subs(match)returnedef_check_varsh_sum_872_4(e):a=Wild('a')alpha=Wild('alpha')b=Wild('b')beta=Wild('beta')c=Wild('c')cp=Wild('cp')gamma=Wild('gamma')gammap=Wild('gammap')match1=e.match(Sum(CG(a,alpha,b,beta,c,gamma)*CG(a,alpha,b,beta,cp,gammap),(alpha,-a,a),(beta,-b,b)))ifmatch1isnotNoneandlen(match1)==8:return(KroneckerDelta(c,cp)*KroneckerDelta(gamma,gammap)).subs(match1)match2=e.match(Sum(CG(a,alpha,b,beta,c,gamma)**2,(alpha,-a,a),(beta,-b,b)))ifmatch2isnotNoneandlen(match2)==6:return1returnedef_cg_list(term):ifisinstance(term,CG):return(term,),1,1cg=[]coeff=1ifnot(isinstance(term,Mul)orisinstance(term,Pow)):raiseNotImplementedError('term must be CG, Add, Mul or Pow')ifisinstance(term,Pow)andsympify(term.exp).is_number:ifsympify(term.exp).is_number:[cg.append(term.base)for_inrange(term.exp)]else:return(term,),1,1ifisinstance(term,Mul):forarginterm.args:ifisinstance(arg,CG):cg.append(arg)else:coeff*=argreturncg,coeff,coeff/abs(coeff)