For a Euler-Bernoulli beam with quite large flexible deformation， due to its quite complicated nonlinear geometric equation， in common only the displacements of discrete nodes can be numerically predicted， by such as the multiple shooting method and the differential quadrature method. This study proposes the variational method for solving the large flexible deformation problem of Euler-Bernoulli beams. Instead of directly assuming the displacement functions， this study derives the exact forms of coupled displacement functions， based on the polynomial function of curvature and constant strain of central line. Further， based on the derived displacement functions， a group of nonlinear equations are derived following the standard procedure of the variational method. The proposed variational method for the quite large flexible deformation problem of Euler-Bernoulli beam is verified by finite element analysis. The predictions on an Euler-Bernoulli cantilever beam demonstrate that， with only six unknown parameters， the proposed variational method is able to accurately predict the quite large flexible deformations of an Euler-Bernoulli beam.