FiniteField

 

f = Function[x, x^2 + 1]  
f = #^2+1 &  
IrreduciblePolynomialQ[f[x], Modulus -> 7] True
F = FiniteField[7, f]  
Information[F, "FieldSize"]

49