  We provide comprehensive regularity results and optimal conditions for a
general class of functionals involving Orlicz multi-phase of the type
\begin{align}
  \label{abst:1}
  v\mapsto \int_{\Omega} F(x,v,Dv)\,dx, \end{align} exhibiting non-standard
growth conditions and non-uniformly elliptic properties.
  The model functional under consideration is given by the Orlicz multi-phase
integral \begin{align}
  \label{abst:2}
  v\mapsto \int_{\Omega} f(x,v)\left[ G(|Dv|) +
\sum\limits_{k=1}^{N}a_k(x)H_{k}(|Dv|) \right]\,dx,\quad N\geqslant 1,
\end{align} where $G,H_{k}$ are $N$-functions and $ 0\leqslant a_{k}(\cdot)\in
L^{\infty}(\Omega) $ with $0 < \nu \leqslant f(\cdot) \leqslant L$. Its
ellipticity ratio varies according to the geometry of the level sets
$\{a_{k}(x)=0\}$ of the modulating coefficient functions $a_{k}(\cdot)$ for
every $k\in \{1,\ldots,N\}$.
  We give a unified treatment to show various regularity results for such
multi-phase problems with the coefficient functions
$\{a_{k}(\cdot)\}_{k=1}^{N}$ not necessarily H\"older continuous even for a
lower level of the regularity. Moreover, assuming that minima of the functional
above belong to better spaces such as $C^{0,\gamma}(\Omega)$ or
$L^{\kappa}(\Omega)$ for some $\gamma\in (0,1)$ and $\kappa\in (1,\infty]$, we
address optimal conditions on nonlinearity for each variant under which we
build comprehensive regularity results.
  On the other hand, since there is a lack of homogeneity properties in the
nonlinearity, we consider an appropriate scaling with keeping the structures of
the problems under which we apply Harmonic type approximation in the setting
varying on the a priori assumption on minima. We believe that the methods and
proofs developed in this paper are suitable to build regularity theorems for a
larger class of non-autonomous functionals.
