  Complex many-body systems, such as triaxial and reflection-asymmetric nuclei,
weakly-bound halo states, cluster configurations, nuclear fragments produced in
heavy-ion fusion reactions, cold Fermi gases, and pasta phases in neutron star
crust, they are all characterized by large sizes and complex topologies, in
which many geometrical symmetries characteristic of ground-state configurations
are broken. A tool of choice to study such complex forms of matter is an
adaptive multi-resolution wavelet analysis. This method has generated much
excitement since it provides a common framework linking many diversified
methodologies across different fields, including signal processing, data
compression, harmonic analysis and operator theory, fractals, and quantum field
theory.
  To describe complex superfluid many-fermion systems, we introduce an adaptive
pseudo-spectral method for solving self-consistent equations of nuclear density
functional theory in three dimensions, without symmetry restrictions.
  The new adaptive multi-resolution Hartree-Fock-Bogoliubov (HFB) solver
{\madnesshfb} is benchmarked against a two-dimensional coordinate-space solver
{\hfbax} based on B-spline technique and three-dimensional solver {\hfodd}
based on the harmonic oscillator basis expansion. Several examples are
considered, including self-consistent HFB problem for spin-polarized trapped
cold fermions and Skyrme-Hartree-Fock (+BCS) problem for triaxial deformed
nuclei.
  The new {\madnesshfb} framework has many attractive features when applied to
nuclear and atomic problems involving many-particle superfluid systems. Of
particular interest are weakly-bound nuclear configurations close to particle
drip lines, strongly elongated and dinuclear configurations such as those
present in fission and heavy ion fusion, and exotic pasta phases that appear in
the neutron star crust.
