Hi!

Series is a possibility, but remember that Series expansions are LOCAL and focus on the neighborhood of some point of interest. Such series do not converge everywhere! In your case, you can compute two series, focusing either on the center of the domain, or on the right side, for instance :

S[\[Phi]_] := (4/
     0.014) (-Log[1 - 0.8])^(1/2) (1 - \[Phi]) (-Log[1 - \[Phi]])^(1/
      2);
DLSc = Normal[
  Series[S[\[Phi]], {\[Phi], 0.4, 10}]]; (* focus on center *)
DLSr = Normal[
  Series[S[\[Phi]], {\[Phi], 0.8, 10}]];(* focus on right side *)
TayC[x_] := DLSc /. \[Phi] -> x;
TayR[x_] := DLSr /. \[Phi] -> x;

Focusing on the right side, one can see that the three approximations are correct, but that TayR works much better than TayC :

Plot[{S[x], TayC[x], TayR[x]}, {x, 0.5, 1}, PlotRange -> All, 
 Frame -> True, PlotLegends -> Automatic]

Now, focusing on the center (0.4, say), you can see that TayC works quite well, excepted near the boundary!

Plot[{S[x], TayC[x](*,TayR[x]*)}, {x, 0, 1}, PlotRange -> All, 
 Frame -> True, PlotLegends -> Automatic]

Bur you could as well use polynomial approximation, either with orthogonal polynomials or Bernstein polynomial,... It depends on what you need!

Hope this helps.

Krishan,

maybe Pade-approximation might be an option, e.g.:

func[\[Phi]_] := (2000/7) (-Log[1/5])^(1/2) (1 - \[Phi]) (-Log[1 - \[Phi]])^(1/2)
Simplify[PadeApproximant[func[\[Phi]], {\[Phi], 0, 3}]]