  We propose a new, efficient multi-scale method to decompose a map (or signal
in general) into components maps that contain structures of different sizes. In
the widely-used wave transform, artifacts containing negative values arise
around regions with sharp transitions due to the application of band-limited
filters. In our approach, the decomposition is achieved by solving a modified,
non-linear version of the diffusion equation. This is inspired by the
anisotropic diffusion methods, which establish the link between image filtering
and partial differential equations. In our case, the artifact issue is
addressed where the positivity of the decomposed images is guaranteed. Our new
method is particularly suitable for signals which contain localized, non-linear
features, as typical of astronomical observations. It can be used to study the
multi-scale structures of astronomical maps quantitatively and should be useful
in observation-related tasks such as background removal. We thus propose a new
measure called the "scale spectrum", which describes how the image values
distribute among different components in the scale space, to describe maps. The
method allows for input arrays of an arbitrary number of dimensions, and a
python3 implementation of the algorithms is included in the Appendix and
available at https://gxli.github.io/Constrained-Diffusion-Decomposition/.
