In this paper we present details on a formulation and numerical solution of the equations of time-dependent magnetized convection. The ADISM (alternating direction implicit on staggered mesh) method of Chan & Wolff is used to solve these equations in two spatial dimensions. Our motivation stems from the need to understand the interaction between highly stratified compressible convection and magnetic fields under solar and stellar conditions. A strong challenge for theoretical models of both small- and large-scale magnetic fields in stars is provided by high-resolution observations of the solar surface and long-term photometric and spectroscopic observations of stars. We present the results of our simulations in the context of previous work and highlight some qualitative properties of the convection-magnetic field interaction. For example, we find that the thermal relaxation time of a convective layer with magnetic fields is primarily influenced by the mean magnetic flux that permeates the layer (i.e., a measurable quantity). We also find in our simulations the presence of acoustic, magnetoacoustic, and Alfven waves, in agreement with observations. We discuss the importance of an adequate treatment and understanding of these waves, especially with the restrictions imposed by artificial boundary conditions. We present results on the transition between weak and strong magnetic fields, where the highly nonlinear behavior commences, and indicate the sensitivity of our results to the value of the magnetic resistivity, which can significantly influence the details of the interactions. Finally we discuss the consequences of our findings for modeling solar and stellar magnetic fields and indicate what future work will be undertaken.