We consider the so-called GI/GI/N queueing network in which a stream of jobs with independent and identically distributed service times arrive according to a renewal process to a common queue served by $N$ identical servers in a first-come-first-serve manner. We introduce a two-component infinite-dimensional Markov process that serves as a diffusion model for this network, in the regime where the number of servers goes to infinity and the load on the network scales as $1-\beta N^{-1/2}+o(N^{-1/2})$ for some $\beta>0$. Under suitable assumptions, we characterize this process as the unique solution to a pair of stochastic evolution equations comprised of a real-valued Ito equation and a stochastic partial differential equation on the positive half line, which are coupled together by a nonlinear boundary condition. We construct an asymptotic (equivalent) coupling to show that this Markov process has a unique invariant distribution. This invariant distribution is shown in a companion paper [Aghajani and Ramanan (2016)] to be the limit of the sequence of suitably scaled and centered stationary distributions of the GI/GI/N network, thus resolving (for a large class service distributions) an open problem raised by Halfin and Whitt in [Oper. Res. 29 (1981) 567–588]. The methods introduced here are more generally applicable for the analysis of a broader class of networks.