The article is focused on a two‐dimensional geometrically nonlinear formulation of a Bernoulli beam element that can accommodate arbitrarily large rotations of cross sections. The formulation is based on the integrated form of equilibrium equations, which are combined with the kinematic equations and generalized material equations, leading to a set of three first‐order differential equations. These equations are then discretized by finite differences and the boundary value problem is converted into an initial value problem using a technique inspired by the shooting method. Accuracy of the numerical approximation is conveniently increased by refining the integration scheme on the element level while the number of global degrees of freedom is kept constant, which leads to high computational efficiency. The element has been implemented into an open‐source finite element code. Numerical examples show a favorable comparison with standard beam elements formulated in the finite‐strain framework and with analytical solutions.