The left-invariant sub-Riemannian problem on the Engel group is considered. The problem gives the nilpotent (non-linear) approximation to generic nonholonomic systems in four-dimensional space with two dimensional control, for instance to a system which describes motion of a mobile robot with a trailer. The global optimality of extremal trajectories is studied via geometric control theory. The global diffeomorphic structure of the exponential mapping is described. As a consequence, the cut time is proved to be equal to the first Maxwell time corresponding to discrete symmetries of the exponential mapping. The Maxwell strata corresponding to the symmetries are described by certain equations in elliptic functions. We study solvability of these equations, give sharp estimates of their roots, and describe their mutual disposition via the analysis of the elliptic functions involved. The problem of finding optimal synthesis in the general case is reduced to a system of three algebraic equations in elliptic functions and elliptic integrals. It seems impossible to analytically solve such equations, therefore, a software for computing optimal trajectories for the nilpotent sub-Riemannian problem on the Engel group is being developed in Wolfram Mathematica. This software has already been used to devise several algorithms for computing approximate paths close to optimal for a mobile robot with a trailer. Those algorithms will be applied for controlling a real mobile robot with a trailer.