### Trouble verifying modularity properties

Hey all,

I'm having a bit of trouble doing some numerical things with modular forms, and I simply can't figure out where I'm going wrong.

The $j$ function should satisfy $j(\gamma \tau) = j(\tau)$ for every $\tau$ in the upper half plane and every $\gamma\in SL_2(\mathbb{Z})$. I wrote some code to numerically compute some values of the $j$ function (there may be a better way to do it for the $j$ function, but I hope to eventually migrate this code to work for other modular forms which I define). Here it is below.

```
## Starts Laurent series in q
R.<q> = LaurentSeriesRing(QQ)
I = CC.0 # imaginary unit
precision = 75
##evaluates a function using its q-expansion
def evaluate(f,z):
result = 0
coeffs = f.coefficients()
exps = f.exponents()
for i in range(0,len(coeffs)):
result = result + coeffs[i]*z^(exps[i])
return result
## computes the action of a member of the modular group on tau in the upper half plane
def action(gamma,tau):
return ((gamma[0]*tau + gamma[1])/(gamma[2]*tau + gamma[3]))
## Produce Eisenstein series with specified weight and precision of q-expansion
def eis(weight,precision):
t = EisensteinForms(1,weight)
t.set_precision(precision)
t = t.eisenstein_series()
e = t[0].q_expansion()
return e*(1/e.list()[0])
## gives you q which corresponds to tau
def qt(tau):
return exp(2*pi*I*tau)
## Defining delta cusp form
delta = CuspForms(1,12).0
delta = delta.q_expansion(precision)
# Computes j function
g2 = 60*eis(4,precision)/240
j = 1728*g2^3/(27*delta)
```

Now when I run the following code:

```
tau = 1+I
gamma = [3,-1,4,-1]
```~~print(evaluate(j,qt(tau))) ~~print(evaluate(j,qt(tau)).n()) #j(tau)
~~print(evaluate(j,qt(action(gamma,tau)))) ~~print(evaluate(j,qt(action(gamma,tau))).n()) #j(gamma tau)

the values $j(\tau)$ and $j(\gamma\tau)$ are not equal! I would appreciate any help.