Force-free electrodynamics (FFE) describes a particular regime of magnetically dominated relativistic plasmas, which arise on several astrophysical scenarios of interest such as pulsars or active galactic nuclei. In those regimes, the electromagnetic field obeys a modified nonlinear version of Maxwell equations, while the plasma only accommodates to locally cancel out the Lorentz force. The aim of the present talk is to discuss the initial/boundary value formulation of FFE at some given astrophysical settings. We start by showing that, when restricted to the correct constraint submanifold, the system is symmetric hyperbolic: we introduce here a particular hyperbolization for the FFE equations [1], following a covariant approach due to R. Geroch [2]. Then, we analyze the characteristic structure of the resulting evolution system and use this information to construct appropriate boundary conditions [3,4]. In particular, we focus on the treatment to mimic the perfectly conducting surface of a neutron star, where incoming and outgoing physical modes needs to be combined on a very precise way. We shall illustrate this procedure with the simpler vacuum (linear) electrodynamics. Also, we discuss the methods employed to deal with the constraints of the theory at the boundaries. And finally, we show some results from our 3D numerical simulations based on this approach [4,5,6]. [1] F. Carrasco, O. Reula. PRD (93), 2016. DOI: 10.1103/PhysRevD.93.085013 [2] R. Geroch. “Partial differential equations of physics”. In General Relativity, pp. 19-60. Routledge, 1996. [3] F. Carrasco, O. Reula. PRD (96), 2017. [4] F. Carrasco, C. Palenzuela, O. Reula. PRD (98), 2018. DOI: 10.1103/PhysRevD.98.023010 [5] F. Carrasco, D. Viganò. C. Palenzuela, J. Pons. MNRAS Letters (484), 2019. [6] R. Cayuso, F. Carrasco, B. Sbarato, O. Reula., 2019.