Understanding the evolution of solute concentration gradients underpins the prediction of porous media processes limited by mass transfer. Here, we present the development of a mathematical model that describes the dissolution of spherical bubbles in two-dimensional regular pore networks. The model is solved numerically for lattices with up to 169 bubbles by evaluating the role of pore network connectivity, vacant lattice sites and the initial bubble size distribution. In dense lattices, diffusive shielding prolongs the average dissolution time of the lattice, and the strength of the phenomenon depends on the network connectivity. The extension of the final dissolution time relative to the unbounded (bulk) case follows the power-law function, Bk/ℓ, where the constant ℓ is the inter-bubble spacing, B is the number of bubbles, and the exponent k depends on the network connectivity. The solute concentration field is both the consequence and a factor affecting bubble dissolution or growth. The geometry of the pore network perturbs the inward propagation of the dissolution front and can generate vacant sites within the bubble lattice. This effect is enhanced by increasing the lattice size and decreasing the network connectivity, yielding strongly nonuniform solute concentration fields. Sparse bubble lattices experience decreased collective effects, but they feature a more complex evolution, because the solute concentration field is nonuniform from the outset.