Using the largest cosmological reionization simulation to-date (~24 billion particles), we use the genus curve to quantify the topology of neutral hydrogen distribution on scales > 1 Mpc as it evolves during cosmological reionization. We find that the reionization process proceeds primarily in an inside-out fashion, where higher density regions become ionized earlier than lower density regions. There are four distinct topological phases: (1) Pre-reionization at z ~ 15, when the genus curve is consistent with a Gaussian density distribution. (2) Pre-overlap at 10 < z < 15, during which the number of HII bubbles increases gradually with time, until percolation of HII bubbles starts to take effect, characterized by a very flat genus curve at high volume fractions. (3) Overlap at 8 < z < 10, when large HII bubbles rapidly merge, manifested by a precipitous drop in the amplitude of the genus curve. (4) Post-overlap at 6 < z < 8, when HII bubbles have mostly overlapped and the genus curve is consistent with a diminishing number of isolated neutral islands. After the end of reionization (z < 6), the genus of neutral hydrogen is consistent with Gaussian random phase, in agreement with observations.