There is a rather deep polynomial‑time algorithm for counting the $\mathbb F_q$‑rational points of an elliptic curve published by René Schoof in 1985 (with subsequent improvements by Noam Elkies and A. O. L. Atkin). It is based on two core ideas:
- The number of points is closely linked to a functional equation
$$ \varphi^2-t\varphi+q = 0 \qquad\in\operatorname{End}(E)$$
that the Frobenius endomorphism
$$ \varphi\colon\;E\to E,\;\begin{cases}\mathcal O&\mapsto \mathcal O\\(x,y)&\mapsto (x^q,y^q) \end{cases} $$
satisfies in the endomorphism ring of $E$. If $t\in\mathbb Z$ is chosen such that this equation holds, it is called the trace of Frobenius and one can show that
$$ \#E(\mathbb F_q)=q+1-t \text. $$
(The reason $\varphi$ has anything to do with point counting is that it leaves exactly the points with coordinates in $\mathbb F_q$ invariant.)
- For odd $\ell$, there exist division polynomials $\psi_\ell\in\mathbb F_q[x]$ which vanish precisely on the $x$‑coordinates of the finite $\ell$‑torsion points of $E$. Therefore, one can compute $t\bmod\ell$ by checking for which $k\in\{0,\dots,\ell{-}1\}$ the functional equation $\varphi^2-k\varphi+q=0$ holds on a symbolic point $(x,y)$ where $x$ is a hypothetical root of $\psi_\ell$; in other words, this involves evaluating the endomorphism modulo $\psi_\ell$. The modular reduction is what makes this step polynomial‑time: The evaluation now works on symbolic points whose size is polynomially bounded in $\ell$, rather than (without any reduction) in $q$. Having computed $t\bmod\ell$ for a sufficiently large set $L$ of odd primes (according to Hasse's theorem, $\prod_{\ell\in L}\ell>4\sqrt q$ is enough), the Chinese remainder theorem can be used to reconstruct the trace of Frobenius $t$, thereby yielding the number of $\mathbb F_q$‑rational points $\#E(\mathbb F_q)$ on $E$.
Of course, there are implementations of (the improved variants of) this algorithm in popular computer algebra systems. Using SageMath, the number of $\mathbb F_{q^n}$‑rational points of an elliptic curve
$$ y^2+a_1xy+a_3y=x^3+a_2x^2+a_4x+a_6 $$
defined over $\mathbb F_q$ can be computed as:
EllipticCurve(GF(q), [a1,a2,a3,a4,a6]).cardinality(extension_degree=n)