The distributions of event-by-event harmonic flow coefficients v (n) for n = 2- 4 are measured in = 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 |eta| < 2.5 in a dataset of approximately 7 mu 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.