The distributions of event-by-event harmonic flow coefficients v n for n = 2- 4 are measured in sNN−−−−√ = 2.76 TeV Pb + Pb collisions using the ATLAS detector at the LHC. The measurements are performed using charged particles with transverse momentum p T > 0.5 GeV and in the pseudorapidity range |η| < 2.5 in a dataset of approximately 7 μb−1 recorded in 2010. The shapes of the v n distributions suggest that the associated flow vectors are described by a two-dimensional Gaussian function in central collisions for v 2 and over most of the measured centrality range for v 3 and v 4. Significant deviations from this function are observed for v 2 in mid-central and peripheral collisions, and a small deviation is observed for v 3 in mid-central collisions. In order to be sensitive to these deviations, it is shown that the commonly used multi-particle cumulants, involving four particles or more, need to be measured with a precision better than a few percent. The v n distributions are also measured independently for charged particles with 0.5 < p T < 1 GeV and p T > 1 GeV. When these distributions are rescaled to the same mean values, the adjusted shapes are found to be nearly the same for these two p T ranges. The v n distributions are compared with the eccentricity distributions from two models for the initial collision geometry: a Glauber model and a model that includes corrections to the initial geometry due to gluon saturation effects. Both models fail to describe the experimental data consistently over most of the measured centrality range.