In this paper we consider a furtivity problem in the context of time-dependent three-dimensional acoustic obstacle scattering. The scattering problem for a ‘passive’ obstacle is the following: an incoming acoustic wavepacket is scattered by a bounded simply connected obstacle with locally Lipschitz boundary having a known boundary acoustic impedance. The scattered wave is the solution of an exterior problem for the wave equation. To make the obstacle furtive we leave ‘passive’ obstacles and we consider ‘active’ obstacles, that is obstacles that, when hit by the incoming wavepacket, react with a pressure current circulating on their boundary. The furtivity problem consists of making the acoustic field scattered by the obstacle ‘as small as possible’ by choosing a control function, that is a pressure current on the boundary of the obstacle, in the function space of the admissible controls. It consists of finding the control function that minimizes a cost functional that will be made precise later. This furtivity problem is of great relevance in many applications. The mathematical model for this furtivity problem is a control problem for the wave equation. In the boundary condition for the wave equation on the boundary of the obstacle we introduce a control function, the so-called pressure current. The cost functional depends on the control function, and on the scattered acoustic field. Note that the scattered field depends on the control function via the boundary conditions. Using the Pontryagin maximum principle we show that, for a suitable choice of the cost functional, the first-order optimality conditions for the furtivity problem considered can be formulated as an exterior problem defined outside the obstacle for a system of two coupled wave equations. This is the main purpose of the paper. Moreover, to solve this exterior problem numerically we develop a highly parallelizable method based on a ‘perturbative series’ of the type proposed in 1. This method obtains the time-dependent scattered field and the control function as superpositions of time harmonic functions. The space-dependent parts of each time harmonic component of the scattered field and of the control function are obtained by solving an exterior boundary value problem for two coupled Helmholtz equations. The mathematical model and the numerical method proposed are validated by studying some test problems numerically. The results obtained with a parallel implementation of the numerical method proposed on the test problems are shown and discussed from the numerical and the physical point of view. The quantitative character of the results obtained is established. Animations (audio, video) relative to the numerical experiments can be found at stacks.iop.org/WRM/11/549. MThis article features online multimedia enhancements