The increasing use of neural operators and scientific machine learning has enabled accurate surrogate models for parameterized partial differential equations (PDEs). However, most existing approaches focus on learning solution operators for individual PDEs in isolation, often requiring retraining when problem formulations, parameters, or physical components change. In this work, we propose a modular framework for constructing a library of canonical PDE solutions learned via neural networks, which can be systematically combined to approximate solutions of more complex PDEs. The proposed framework is based on independently trained surrogate models for simple PDEs such as the Laplace, Poisson, advection, and diffusion equations. Each PDE is approximated using neural-network-based surrogate models, primarily operator-learning architectures such as DeepONets and Fourier neural operators, which may be trained to predict either solution operators or full PDE solutions, while allowing for the use of alternative physics-informed or data-driven models when appropriate. These learned surrogate models form a reusable library of PDE building blocks. Solutions of more complex PDEs are then constructed by composing the outputs of these surrogate models using established analytical principles, such as superposition, numerical composition techniques including operator splitting methods, or even a unified operator-learning backbone trained to learn the composition map directly. This modular design allows complex PDE solutions to be assembled without training a single monolithic neural model that simultaneously accounts for parameters, boundary conditions, and physical effects. The framework is demonstrated on one-dimensional linear PDEs, including non-homogeneous elliptic problems and advection–diffusion equations through combinations of learned Laplace, Poisson, advection, and diffusion operators. The approach highlights how decomposition can guide surrogate design, reduce model complexity, and enhance flexibility in handling varying parameters and source terms. Ongoing work focuses on extending the framework to higher-dimensional problems, alternative composition strategies such as kernel-based operator composition, and the treatment of nonlinear PDEs.