This is the python script used to do the symbolic calculations that were mentioned in a footnote in article 79649
#!/usr/bin/env python3
# ================================================================
# This script does symbolic calculations for article 79649
#
# Written by Randy S (2023)
# ================================================================
import copy
from math import factorial
# https://docs.sympy.org/latest/index.html
from sympy import symbols, Eq, solveset, expand, simplify
dims = 5
maxn = 14
# ===================================================================
def show_term_in_LaTeX( term ):
output = '&'
if term['coef'] > 0:
output += '+'
output += str( term['coef'] ) + '\\frac{'
Vexp = 1
for k in range(len(term['Vfactors'])):
if Vexp > 1:
Vexp -= 1
else:
output += 'g_{' + str( term['Vfactors'][k] ) + '}'
Vexp = 1
while k+Vexp < len(term['Vfactors']) and term['Vfactors'][k+Vexp] == term['Vfactors'][k]:
Vexp += 1
if Vexp > 1:
output += '^{' + str(Vexp) + '}'
output += '}{(1+g_2)^{' + str( term['Dexp'] ) + '}}'
print( output )
print( '\cr')
# ===================================================================
def show_bterms_in_LaTeX( n, terms ):
print( '\\begin{align*}' )
print( '\\beta_{' + str(n) + '} &= ' )
print( '(' + str(n) + ' - ' + str(int(n/2-1)) + 'd)g_{' + str(n) + '} ' )
print( '\cr')
for k in range(len(terms)):
show_term_in_LaTeX( terms[k] )
print( '\end{align*}' )
# ===================================================================
def has_odds( term ):
any_odd = False
for k in range(len(term['Vfactors'])):
if term['Vfactors'][k] % 2 == 1:
any_odd = True
return any_odd
# ===================================================================
def solve_for_highest_g( n, terms ):
result = []
for k in range(len(terms)):
newterm = copy.deepcopy( terms[k] )
if newterm['Vfactors'][0] > n:
# replace with the tree-level term
newterm['Dexp'] = -1
newterm['coef'] = -(n - int(n/2-1)*dims)
newterm['Vfactors'] = []
newterm['Vfactors'].append(n)
else:
newterm['Dexp'] -= 1
newterm['coef'] *= -1
result.append( newterm )
return result
# ===================================================================
def discard_odds( terms ):
result = []
for k in range(len(terms)):
if not has_odds( terms[k] ):
result.append( copy.deepcopy( terms[k] ) )
return result
# ===================================================================
def take_deriv_of_one_term( term ):
result = []
# append the term that comes from the derivative hitting the denominator
newfactors = copy.deepcopy( term['Vfactors'] )
newfactors.append(3)
result.append( {'Dexp': term['Dexp']+1, \
'Vfactors': newfactors, \
'coef': -term['Dexp'] * term['coef'] } )
# append the terms that come from the derivative hitting each factor in the numerator
for k in range(len(term['Vfactors'])):
newterm = copy.deepcopy( term )
newterm['Vfactors'][k] = newterm['Vfactors'][k] + 1
result.append( newterm )
return result
# ===================================================================
def consolidate_result( result, result_single ):
for k in range(len(result_single)):
result_single[k]['Vfactors'].sort( reverse = True )
for k in range(len(result_single)):
found_match = False
for j in range(len(result)):
if result_single[k]['Vfactors'] == result[j]['Vfactors']:
# these terms have the same factors, so consolidate them
result[j]['coef'] += result_single[k]['coef']
found_match = True
if not found_match:
result.append( result_single[k] )
# ===================================================================
def take_deriv( terms ):
result = []
for k in range(len(terms)):
result_single = take_deriv_of_one_term( terms[k] )
consolidate_result( result, result_single )
return result
# ===================================================================
def construct_beta_for_SymPy( n, bterms, g ):
k = int(n/2-1)
result = (n-k*dims) * g[n]
for k in range(len(bterms)):
thisterm = bterms[k]['coef']
for j in range(bterms[k]['Dexp']):
thisterm = thisterm/(1+g[2])
for j in range(len(bterms[k]['Vfactors'])):
thisterm = thisterm * g[bterms[k]['Vfactors'][j]]
result = result + thisterm
return result
# ===================================================================
def do_calculations_with_SymPy( beta, g ):
# abbreviation
d = dims
print( '' )
for n in range(2,maxn+1,2):
print( 'beta[' + str(n) + '] = ', beta[n] )
# solve the equation beta[n] = 0 for g[n+2]
geqtmp = ['']*(maxn+3)
for n in range(2,maxn+1,2):
geqtmp[n+2] = solveset( beta[n], g[n+2] )
# extract the one-and-only solution in each case
geq = ['']*(maxn+3)
geq[2] = g[2]
for n in range(2,maxn+1,2):
geq[n+2] = list(geqtmp[n+2])[0]
print( '' )
for n in range(2,maxn+1,2):
print( 'beta[' + str(n) + '] = 0 implies g' + str(n+2) + ' = ', geq[n+2] )
# express all of the geq's entirely in terms of g2
for n in range(2,maxn+1,2):
for k in range(n,2,-2):
geq[n+2] = geq[n+2].subs( g[k], geq[k] )
for n in range(2,maxn+1,2):
geq[n+2] = expand( geq[n+2] )
print( '' )
for n in range(2,maxn+1,2):
print( 'beta[n] = 0 for all n implies g' + str(n+2) + ' = ', geq[n+2] )
print( '' )
print( 'h[n]\equiv g[n]/n!' )
for n in range(2,maxn+1,2):
print( 'h[' + str(n+2) + '] = ', simplify(geq[n+2]/factorial(n+2)) )
# calculate ratios to study convergence of the series
ratio = ['']*(maxn+3)
for n in range(2,maxn+1,2):
ratio[n] = simplify( geq[n+2]/(geq[n]*(n+2)*(n+1)) )
print( '' )
print( 'r[n]\equiv (g[n+2]/(n+2)!)/(g[n]/n!)' )
for n in range(2,maxn+1,2):
print( 'r[' + str(n) + '] = ', ratio[n] )
# ===================================================================
# Main program
# ---------------------------------
# Initialize variables for SymPy
#
# manual version: g = symbols( 'g0 g1 g2 g3 g4 g5 g6' )
# ---------------------------------
g = ['']*(maxn+3)
for n in range(2,maxn+3,2):
g[n] = symbols( 'g' + str(n) )
beta = ['']*(maxn+1)
# ---------------------------------
# Initialize my home-grown symbolic calculation
#
# start with the first derivative of log(D)
# ---------------------------------
terms = [{'Dexp': 1, 'coef': 1, 'Vfactors': [3]}]
# ---------------------------------
# Construct the beta-functions by
# calculating higher derivatives iteratively
# using my home-grown symbolic calculations
# ---------------------------------
for n in range(2,maxn+1):
terms = copy.deepcopy( take_deriv( terms ) )
if n%2 == 0:
# bterms = terms in \beta_n(\vec g)
bterms = copy.deepcopy( discard_odds( terms ) )
show_bterms_in_LaTeX( n, bterms )
beta[n] = construct_beta_for_SymPy( n, bterms, g )
print( '' )
print( '\\begin{verbatim}' )
print( 'beta[' + str(n) + '] = ', beta[n] )
print( '\\end{verbatim}' )
print( '' )
# ---------------------------------
# Use SymPy to solve the conditions beta[n] = 0
# ---------------------------------
print( '\\begin{verbatim}' )
do_calculations_with_SymPy( beta, g )
print( '\\end{verbatim}' )