The paper considers shape optimization of Euler-Bernoulli beams with circular, square and rectangular cross-sections made of axially functionally graded materials at a prescribed fundamental frequency. Optimization is carried out by the beam mass minimization. Considerations involve the case of coupled bending and axial vibrations, where complex boundary conditions are the cause of coupling. Pontryagin's maximum principle is used to solve shape optimization, where a limited diameter or a beam cross-sectional width is used for control. Diameter limit is considered so that the optimized shape of a beam is within the limits of the validity of Euler-Bernoulli theory, and its strength does not decrease for smaller cross-sectional dimensions. The resulting system of differential equations is a two-point boundary value problem, and the shooting method is applied to solve it. The property of self-coupled systems is utilized, where all adjoint variables, except for one variable, are expressed through state variables, which facilitates solving the appropriate differential equations. Theoretical considerations are illustrated by an example. Also, the savings of beam mass in percent are determined, using the cantilever beam with optimal variable cross-section against the cantilever beam of a constant cross-section, where both beams have the same prescribed fundamental frequency.