Polymeric drug delivery platforms offer promising capabilities for controlled drug release thanks to their ability to be custom-designed with specific properties. In this paper, we present a model to simulate the complex interplay between solvent absorption, polymer swelling, drug release and stress development within these types of platforms. A system of nonlinear partial differential equations coupled with an ordinary differential equation is introduced to avoid drawbacks from other models found in the literature. These incorporated a memory effect to account for polymer relaxation but from a numerical point of view, required storing information from all previous time steps, making them computationally expensive. This paper proposes a new numerical method to simulate such drug delivery devices based on nonuniform grids and an implicit midpoint time discretization. Our main results are the proof of second-order convergence of the method for nonsmooth solutions and the scheme's stability under the assumption of quasiuniform grids and a sufficiently small timestep. We also illustrate numerically the second-order convergence result proven in the main result using solutions based on biological information.