Attitude determination using double-differenced GNSS carrier phase measurements is studied. A realistic stochastic model is employed to take the correlations among the double-differenced measurements into full consideration. Two important issues concerning iteratively solving the nonlinear least-squares attitude determination problem are treated, namely, the initial guess and the iteration scheme. An analytical and sub-optimal solution is employed to provide the initial guess. In this solution, the orthogonal and determinant constraints among the elements of the direction cosine matrix (DCM) of the attitude are firstly ignored, and hence a relaxed 3×3 matrix is estimated using the linear weighted least-squares method. Then a mathematically feasible DCM, i.e., orthogonal and with +1 determinant, is extracted from the relaxed matrix estimate, optimally in the sense of minimum Frobenius norm. This analytical initial guess estimation method can be used for all feasible cases, including some generated ones, e.g., the case with only 3 antennas and only 3 satellites, subject possibly to some necessary, yet minor modifications. In each iteration, an error attitude, whose DCM is parameterized using the Gibbs vector, is introduced to relate the previously estimated and the true DCM. By linearizing the measurement model at the zero Gibbs vector, the least-squares estimate of the Gibbs vector is obtained and then used to correct the previously estimated DCM. By repeating this process, the truly least-squares estimate of the attitude can be achieved progressively. These are in fact Gauss-Newton iterations. For the final estimate, the variance covariance matrix (VCM) of the attitude estimation error can be retained to evaluate or predict the estimation accuracy. The extraction of the widely used rollpitch-yaw angles and the VCM of their additive estimation errors from the final solution is also presented. Numerical experiments are conducted to check the performance of the developed theory. For the case with 3 2-meter long and orthogonally mounted baselines, 5 visible satellites, and 5-millimeter standard deviations of the carrier phase measurements, the root mean squared errors (RMSE) of the roll-pitch-yaw angles in the analytical solution are well below 0.5 degrees, and the estimates converge after only one iteration, with all three RMSEs below 0.2 degrees.