Galactic-scale winds are a generic feature of massive galaxies with high star formation rates across a broad range of redshifts . Despite their importance , a detailed physical understanding of what drives these mass-loaded global flows has remained elusive . In this paper , we explore the dynamical impact of cosmic rays by performing the first three-dimensional , adaptive mesh refinement simulations of an isolated starbursting galaxy that includes a basic model for the production , dynamics and diffusion of galactic cosmic rays . We find that including cosmic rays naturally leads to robust , massive , bipolar outflows from our 10 ^ { 12 } M _ { \odot } halo , with a mass-loading factor \dot { M } / { SFR } = 0.3 for our fiducial run . Other reasonable parameter choices led to mass-loading factors above unity . The wind is multiphase and is accelerated to velocities well in excess of the escape velocity . We employ a two-fluid model for the thermal gas and relativistic CR plasma and model a range of physics relevant to galaxy formation , including radiative cooling , shocks , self-gravity , star formation , supernovae feedback into both the thermal and CR gas , and isotropic CR diffusion . Injecting cosmic rays into star-forming regions can provide significant pressure support for the interstellar medium , suppressing star formation and thickening the disk . We find that CR diffusion plays a central role in driving superwinds , rapidly transferring long-lived CRs from the highest density regions of the disk to the ISM at large , where their pressure gradient can smoothly accelerate the gas out of the disk .