3 ¯ hv ∼ e 2 we can rewrite this second condition simply as l$_{E}$ ≪ d . Thus, the Thomas-Fermi equation (9) for the equilibrium charge density and the hydrodynamic equa- tion (5) for its variation are applicable as long as flake is not useful. Instead, the equation for the n -th fre- quency gap can be obtained by setting q = 0 in Eq. (5). We observe that where β$_{n}$ are the eigenvalues of the equation, ω 2 $_{n}$(0) = β$_{n}$ e $^{2}$v ¯ hl$_{E}$d , (15) l$_{E}$ ≪ d, q ≪ 1 /l$_{E}$. (10) 2 √$_{π}$ d dξ √ | ξ | (1 − ξ $^{2}$)1 / 4 ∫ 1 − 1 dξ $_{′}$δρ ( n $^{)}$( ξ $^{′}$) ξ − ξ ′ = β$_{n}$δρ ( n $^{)}$( ξ ) . (16) ρ$_{0}$ ( x ) = E$_{0}$x √$_{d}$ 2 − x $_{2}$. (11) β$_{1}$ = 1 . 41 , β$_{2}$ = 6 . 49 , β$_{3}$ = 6 . 75 , ... (17) However, the ratio of q and 1 /d can be arbitrary. For a moderate external electric field ∼ 10 $^{4}$V/m the value of electric length l$_{E}$ ∼ 0 . 4 µ m and the first of the conditions (10) is satisfied easily for micron-sized samples. Analytic solution of Eq. (9) is possible when the second term is small, in which case the charge density is [15] The zeroth mode β$_{0}$ = 0, see Eq. (14), is found ana- lytically: δρ (0) ∝ 1 / √ 1 − ξ $^{2}$. It describes charge dis- tribution in the strip in response to a (uniform along x direction and smooth along y -direction) change of its chemical potential [16]. Other solutions of Eq. (16) are found numerically: Substituting this expression back into Eq. (9) we ob- serve that the second term is indeed negligible as long as x ≫ l 2 $_{E}$/d . This is assured whenever the condi- tions (10) are satisfied. It is also worth pointing out that Eq. (11) justifies the linear approximation for the charge density used in deriving Eq. (1) for q ≫ 1 /d , with ρ ′ $_{0}$/e = 1 / ( l 2 $_{E}$d ). With increasing n the eigenmodes of integro-differential equation (16) oscillate faster, but in general do not follow the oscillation theorem familiar from quantum mechan- ics. In particular, the solutions with n = 0 and n = 3 are even while n = 1, n = 2 are odd [17]. We now turn to the analysis of plasma oscillations propagating on top of the density distribution, Eq. (11). For small plasmon momenta, q ≪ 1 /d , electric field ex- tends beyond the width of the flake and the equation (5) needs to be supplemented with the boundary condition, which ensures that electric field (and thus the current) vanishes at the edges, x = ± d : Finally, we mention the case of a gate-controlled p - n junction, Fig. 1b. The equilibrium density profile is linear near x = 0 and saturates for large | x | [18]. Eq. (1) is still applicable for q > 1 /d . In the limit q < 1 /d one should take into account the screening of long-range Coulomb interaction by metallic gates. In this case the logarithm in the spectrum of the gapless plasmon disappears, and the lowest mode Eq. (14) becomes sound-like. Magnetoplasmons. If external magnetic field B is ap- plied perpendicularly to the plane of graphene the plas- mon spectra acquire new modes. The equation of motion (2) should now be modified to include the Lorentz force, P ∫ d − d dx δρ ( x ) x ± d = 0 . (12) ˙ J ( r , t ) = e 2 π ¯ h 2 | µ ( x ) | E ( r , t ) − ev 2 cµ ( x ) J × B . (18) The spectrum of the lowest symmetric mode can be most easily found by integrating Eq. (5) across the width of the flake. The first term in the brackets will then van- ish exactly due to the boundary condition (12). The remaining integral can now be calculated to the log- arithmic accuracy with the help of the approximation K$_{0}$ ( q | x − x $^{′}$| ) = − ln q | x − x $^{′}$| : ∫ d − d dx √ | ρ$_{0}$ ( x ) | e ln ( q | x − x $^{′}$| ) ≈ 2 d Γ $^{2}$(3 / 4) $_{l$_{E}$}$√$_{π}$ ln ( qd ) . (13) Eqs. (5) and (13) combine to give the equation, [ ω 2 − ω 2 $_{0}$( q )] ∫ d − $_{d}$dxδρ ( x ) = 0, that yields the dispersion of the gapless symmetric plasmon, reminiscent of the plasmon spectrum in quasi-one- dimensional wires, The remaining modes, n ≥ 1, are gapped. For these modes ∫ d − $_{d}$dxδρ ( x ) = 0 and simple procedure of integrating Eq. (5) over the width of the ω 2 $_{0}$( q ) = Γ $^{2}$(3 / 4) 4 e $^{2}$vd π ¯ hl$_{E}$ q $^{2}$ln (1 /qd ) , (14) δρ ( x ) + 2 e 2 π { q $^{2}$Z − q ω ( ω$_{B}$ Z ) ′ − d dx Z d dx } × ∫ d − d dx $^{′}$δρ ( x $^{′}$) K$_{0}$ ( | q || x − x $^{′}$| ) = 0 , (19) The relative coefficient between electric and magnetic terms in this equation follows from the expression for the Lorentz force acting on a single particle. The last term has opposite sign for electrons and holes. Note that the frequency of cyclotron motion ω$_{B}$ ( x ) = ev $^{2}$B/cµ ( x ) in graphene p - n junctions is position-dependent. The remaining equations (3)-(4) are intact in the presence of magnetic field. The boundary condition requires now the vanishing of the normal component of electric current at the boundary, rather than simply vanishing of the elec- tric field, as in Eq. (12). Eliminating J and E we arrive at the generalization of equation (5),