Let $OmegasubsetR^n$ be a smooth bounded domain and let $a_1,a_2,dots,a_{i_0}inOmega$, $widehat{Omega}=Omegasetminus{a_1,a_2,dots,a_{i_0}}$ and $widehat{R^n}=R^nsetminus{a_1,a_2,dots,a_{i_0}}$. We prove the existence of solution $u$ of the fast diffusion equation $u_t=Delta u^m$, $u>0$, in $widehat{Omega}times (0,infty)$ ($widehat{R^n}times (0,infty)$ respectively) which satisfies $u(x,t)toinfty$ as $xto a_i$ for any $t>0$ and $i=1,cdots,i_0$, when $0<m<frac{n-2}{n}$, $ngeq 3$, and the initial value satisfies $0le u_0in L^p_{loc}(2{Omega}setminus{a_1,cdots,a_{i_0}})$ ($u_0in L^p_{loc}(widehat{R^n})$ respectively) for some constant $p>frac{n(1-m)}{2}$ and $u_0(x)ge lambda_i|x-a_i|^{-gamma_i}$ for $xapprox a_i$ and some constants $gamma_i>frac{2}{1-m},lambda_i>0$, for all $i=1,2,dots,i_0$. We also find the blow-up rate of such solutions near the blow-up points $a_1,a_2,dots,a_{i_0}$, and obtain the asymptotic large time behaviour of such singular solutions. More precisely we prove that if $u_0gemu_0$ on $widehat{Omega}$ ($widehat{R^n}$, respectively) for some constant $mu_0>0$ and $gamma_1>frac{n-2}{m}$, then the singular solution $u$ converges locally uniformly on every compact subset of $widehat{Omega}$ (or $widehat{R^n}$ respectively) to infinity as $ttoinfty$. If $u_0gemu_0$ on $widehat{Omega}$ ($widehat{R^n}$, respectively) for some constant $mu_0>0$ and satisfies $lambda_i|x-a_i|^{-gamma_i}le u_0(x)le lambda_i|x-a_i|^{-gamma_i}$ for $xapprox a_i$ and some constants $frac{2}{1-m}<gamma_ilegamma_i<frac{n-2}{m}$, $lambda_i>0$, $lambda_i>0$, $i=1,2,dots,i_0$, we prove that $u$ converges in $C^2(K)$ for any compact subset $K$ of $2{Omega}setminus{a_1,a_2,dots,a_{i_0}}$ (or $widehat{R^n}$ respectively) to a harmonic function as $ttoinfty$.