[proofplan]
We prove the identity law directly from the flex condition at $O$. For associativity, we use the classical nine-point principle for plane cubics: if two plane cubic curves meet a smooth plane cubic $E$ in nine intersection points counted with multiplicity and agree on eight of them, then they agree on the ninth. Applying this to two carefully chosen reducible cubic curves built from three lines shows that the two third points controlling $(P \oplus Q)\oplus R$ and $P\oplus(Q\oplus R)$ are equal. Finally, applying the involution $\iota$ to that equality gives associativity.
[/proofplan]
[step:Establish the inverse-through-$O$ map and the identity law]
For $S \in E(k)$, write $-S := \iota(S)$. By definition, $-O = O$. If $S \ne O$, the line through $S$ and $O$ is defined over $k$, and its third intersection point with the $k$-curve $E$ is again $k$-rational, because the degree-three intersection polynomial has two $k$-rational roots already accounted for by $S$ and $O$.
Since $O$ is a flex point, the tangent line to $E$ at $O$ meets $E$ at $O$ with multiplicity $3$. Therefore
\begin{align*}
O \circ O = O.
\end{align*}
For any $P \in E(k)$, the line through $P$ and $O$ meets $E$ in the three points $P$, $O$, and $-P$, counted with multiplicity. Hence
\begin{align*}
P \circ O = -P.
\end{align*}
Applying $\iota$ gives
\begin{align*}
P \oplus O = \iota(P \circ O) = \iota(-P) = P.
\end{align*}
The same line is used for $O \circ P$, so
\begin{align*}
O \oplus P = P.
\end{align*}
Thus $O$ is the identity element for $\oplus$.
[guided]
The point $O$ is chosen to be a flex because the identity law depends exactly on how the tangent line at $O$ intersects the cubic. Since $O$ is a flex, the tangent line at $O$ has intersection multiplicity $3$ with $E$ at $O$, so the third intersection point of this tangent line with $E$ is again $O$. Thus
\begin{align*}
O \circ O = O.
\end{align*}
Now fix $P \in E(k)$. The projective line through $P$ and $O$ is defined over $k$. Its intersection with the plane cubic $E$ has total multiplicity $3$. By the definition of the inverse-through-$O$ map, the three intersection points are
\begin{align*}
P,\qquad O,\qquad -P,
\end{align*}
where $-P := \iota(P)$. Therefore the third point on the line through $P$ and $O$ is $-P$, meaning
\begin{align*}
P \circ O = -P.
\end{align*}
Applying $\iota$ to this equality gives
\begin{align*}
P \oplus O = \iota(P \circ O) = \iota(-P).
\end{align*}
The same line through $P$ and $O$ also shows that the third point associated to $-P$ and $O$ is $P$, so $\iota(-P)=P$. Hence
\begin{align*}
P \oplus O = P.
\end{align*}
Because the line through $O$ and $P$ is the same line as the line through $P$ and $O$, the identical argument gives
\begin{align*}
O \oplus P = P.
\end{align*}
Thus $O$ is a two-sided identity for $\oplus$.
[/guided]
[/step]
[step:Name the auxiliary third intersection points for associativity]
Fix $P,Q,R \in E(k)$. Define $A \in E(k)$ by
\begin{align*}
A := P \circ Q.
\end{align*}
Define $B \in E(k)$ by
\begin{align*}
B := Q \circ R.
\end{align*}
Then
\begin{align*}
P \oplus Q = -A,\qquad Q \oplus R = -B.
\end{align*}
Define $C \in E(k)$ by
\begin{align*}
C := (-A) \circ R.
\end{align*}
Define $D \in E(k)$ by
\begin{align*}
D := P \circ (-B).
\end{align*}
Then
\begin{align*}
(P \oplus Q)\oplus R = -C,\qquad P \oplus (Q \oplus R) = -D.
\end{align*}
Therefore it is enough to prove $C=D$.
[/step]
[step:Compare the two intersection divisors on $E$]
Let $L_{X,Y} \subset \mathbb{P}^2_k$ denote the projective line through two points $X,Y \in E(k)$, with $L_{X,X}$ interpreted as the tangent line to $E$ at $X$. For a plane curve $G \subset \mathbb{P}^2_k$ not containing $E$, let $\operatorname{div}_E(G)$ denote the effective intersection divisor on the smooth projective curve $E$, with coefficients equal to the local intersection multiplicities. Smoothness of $E$ makes these local intersection multiplicities defined, projectivity gives total degree by [Bezout's Theorem](/page/Bezout%27s%20Theorem), and the [genus formula for smooth plane curves](/page/Genus%20of%20a%20Smooth%20Plane%20Curve) implies that a nonsingular projective plane cubic has genus $1$.
Define the reducible plane cubic curve $\Gamma_1 \subset \mathbb{P}^2_k$ by
\begin{align*}
\Gamma_1 := L_{P,Q} \cup L_{-A,R} \cup L_{B,-B}.
\end{align*}
Define the reducible plane cubic curve $\Gamma_2 \subset \mathbb{P}^2_k$ by
\begin{align*}
\Gamma_2 := L_{Q,R} \cup L_{P,-B} \cup L_{A,-A}.
\end{align*}
The definitions of $\circ$ and $-$ are definitions of intersection divisors with multiplicity, so
\begin{align*}
\operatorname{div}_E(L_{P,Q}) = P + Q + A.
\end{align*}
Also
\begin{align*}
\operatorname{div}_E(L_{-A,R}) = -A + R + C.
\end{align*}
Finally
\begin{align*}
\operatorname{div}_E(L_{B,-B}) = B + (-B) + O.
\end{align*}
Adding these three line-section divisors gives
\begin{align*}
\operatorname{div}_E(\Gamma_1) = P + Q + A + (-A) + R + C + B + (-B) + O.
\end{align*}
Similarly,
\begin{align*}
\operatorname{div}_E(\Gamma_2) = Q + R + B + P + (-B) + D + A + (-A) + O.
\end{align*}
These are equalities of effective divisors on $E$, not merely lists of set-theoretic points, so they remain valid when some of the displayed points coincide.
We now use the precise nine-point principle needed here: if $G_1$ and $G_2$ are plane cubic curves not containing the smooth plane cubic $E$, and if
\begin{align*}
\operatorname{div}_E(G_1) = Z + X
\end{align*}
and
\begin{align*}
\operatorname{div}_E(G_2) = Z + Y
\end{align*}
for the same effective degree-eight divisor $Z$ on $E$ and points $X,Y \in E$, then $X=Y$. Indeed, let $g_1$ and $g_2$ be homogeneous cubic equations for $G_1$ and $G_2$. Since neither $G_1$ nor $G_2$ contains $E$, the restrictions of $g_1$ and $g_2$ to $E$ are nonzero rational sections of the line bundle $\mathcal{O}_E(3)$, and their quotient restricts to a rational function $f:E\dashrightarrow \mathbb{P}^1_k$. Its [principal divisor](/page/Principal%20Divisor) is
\begin{align*}
\operatorname{div}(f) = X - Y.
\end{align*}
If $X \ne Y$, then this divisor is nonzero, so $f$ is nonconstant. Its pole divisor is exactly $Y$, hence has degree $1$, and the degree of the induced morphism $E \to \mathbb{P}^1_k$ is therefore $1$. A [degree-one morphism](/page/Degree%20of%20a%20Morphism) of smooth projective curves is [birational](/page/Birational%20Map), but birational smooth projective curves have the same [genus](/page/Genus), while $E$ has genus $1$ and $\mathbb{P}^1_k$ has genus $0$. This is impossible. Hence $X=Y$.
Apply this principle with
\begin{align*}
Z := P + Q + R + A + (-A) + B + (-B) + O.
\end{align*}
The displayed divisor identities give $\operatorname{div}_E(\Gamma_1)=Z+C$ and $\operatorname{div}_E(\Gamma_2)=Z+D$. Therefore the residual points coincide:
\begin{align*}
C = D.
\end{align*}
[guided]
The delicate point is that we must compare divisors, not just visible points in a picture. Let $L_{X,Y} \subset \mathbb{P}^2_k$ denote the projective line through $X,Y \in E(k)$, with $L_{X,X}$ interpreted as the tangent line to $E$ at $X$. For a plane curve $G \subset \mathbb{P}^2_k$ not containing $E$, define $\operatorname{div}_E(G)$ to be the effective intersection divisor on $E$ whose coefficient at a point is the local intersection multiplicity of $G$ and $E$ there. This is the right object because coincident points are then counted with the correct multiplicity. The smoothness hypothesis on $E$ supplies the local intersection theory on the curve, projectivity allows the degree count from [Bezout's Theorem](/page/Bezout%27s%20Theorem), and the [genus formula for smooth plane curves](/page/Genus%20of%20a%20Smooth%20Plane%20Curve) gives that a nonsingular projective plane cubic has genus $1$.
Define $\Gamma_1 \subset \mathbb{P}^2_k$ by
\begin{align*}
\Gamma_1 := L_{P,Q} \cup L_{-A,R} \cup L_{B,-B}.
\end{align*}
Define $\Gamma_2 \subset \mathbb{P}^2_k$ by
\begin{align*}
\Gamma_2 := L_{Q,R} \cup L_{P,-B} \cup L_{A,-A}.
\end{align*}
Because $A=P\circ Q$, the line $L_{P,Q}$ has intersection divisor
\begin{align*}
\operatorname{div}_E(L_{P,Q}) = P + Q + A.
\end{align*}
This statement includes the tangent case and the case where two of the symbols denote the same point, because the coefficients are intersection multiplicities. Likewise, because $C=(-A)\circ R$,
\begin{align*}
\operatorname{div}_E(L_{-A,R}) = -A + R + C.
\end{align*}
Since $-B$ is defined by the line through $B$ and $O$,
\begin{align*}
\operatorname{div}_E(L_{B,-B}) = B + (-B) + O.
\end{align*}
Adding the three line-section divisors gives
\begin{align*}
\operatorname{div}_E(\Gamma_1) = P + Q + A + (-A) + R + C + B + (-B) + O.
\end{align*}
The same reasoning for the three components of $\Gamma_2$ gives
\begin{align*}
\operatorname{div}_E(\Gamma_2) = Q + R + B + P + (-B) + D + A + (-A) + O.
\end{align*}
Now we state the nine-point principle in the form actually being used. Suppose $G_1$ and $G_2$ are plane cubic curves not containing $E$, and suppose their intersection divisors with $E$ have the form
\begin{align*}
\operatorname{div}_E(G_1) = Z + X
\end{align*}
and
\begin{align*}
\operatorname{div}_E(G_2) = Z + Y
\end{align*}
where $Z$ is the same effective divisor of degree $8$ and $X,Y \in E$. Choose homogeneous cubic equations $g_1$ and $g_2$ defining $G_1$ and $G_2$. Since neither cubic curve contains $E$, the restrictions of $g_1$ and $g_2$ to $E$ are nonzero rational sections of $\mathcal{O}_E(3)$. Quotienting two nonzero rational sections of the same line bundle gives a rational function
\begin{align*}
f:E\dashrightarrow \mathbb{P}^1_k.
\end{align*}
The zero and pole orders of this quotient are exactly the differences of the two intersection divisors, so its [principal divisor](/page/Principal%20Divisor) is
\begin{align*}
\operatorname{div}(f) = X - Y.
\end{align*}
If $X \ne Y$, then $\operatorname{div}(f)\ne 0$, so $f$ is nonconstant. Its pole divisor is exactly the single point $Y$ with multiplicity $1$, so the pole divisor has degree $1$. For a nonconstant rational function on a smooth projective curve, the degree of the induced morphism to $\mathbb{P}^1_k$ equals the degree of its pole divisor; hence $f$ gives a degree-one morphism $E \to \mathbb{P}^1_k$. A [degree-one morphism](/page/Degree%20of%20a%20Morphism) of smooth projective curves is [birational](/page/Birational%20Map). That cannot occur here because birational smooth projective curves have the same [genus](/page/Genus), while $E$ has genus $1$ and $\mathbb{P}^1_k$ has genus $0$. Hence $X=Y$.
We apply this principle with
\begin{align*}
Z := P + Q + R + A + (-A) + B + (-B) + O.
\end{align*}
The divisor identities above say exactly that $\operatorname{div}_E(\Gamma_1)=Z+C$ and $\operatorname{div}_E(\Gamma_2)=Z+D$. Therefore the residual degree-one divisors are equal, and hence
\begin{align*}
C = D.
\end{align*}
This proves the needed geometric equality in all cases, including tangent cases and cases where several of the named points coincide.
[/guided]
[/step]
[step:Apply the inverse-through-$O$ map to obtain associativity]
Since $C=D$, applying $\iota$ gives
\begin{align*}
-C = -D.
\end{align*}
Using the definitions of $C$ and $D$, we compute first
\begin{align*}
(P \oplus Q)\oplus R = (-A)\oplus R = \iota((-A)\circ R) = \iota(C) = -C.
\end{align*}
We also compute
\begin{align*}
P\oplus(Q\oplus R) = P\oplus(-B) = \iota(P\circ(-B)) = \iota(D) = -D.
\end{align*}
Therefore
\begin{align*}
(P \oplus Q)\oplus R = P\oplus(Q\oplus R).
\end{align*}
Since $P,Q,R \in E(k)$ were arbitrary, $\oplus$ is associative.
[/step]