Asymptotic properties of three estimators of probability density function of sample maximum $f_{(m)}:=mfF^{m-1}$ are derived, where $m$ is a function of sample size $n$. One of the estimators is the parametrically fitted by the approximating generalized extreme value density function. However, the parametric fitting is misspecified in finite $m$ cases. The misspecification comes from mainly the following two: the difference $m$ and the selected block size $k$, and the poor approximation $f_{(m)}$ to the generalized extreme value density which depends on the magnitude of $m$ and the extreme index $γ$. The convergence rate of the approximation gets slower as $γ$ tends to zero. As alternatives two nonparametric density estimators are proposed which are free from the misspecification. The first is a plug-in type of kernel density estimator and the second is a block-maxima-based kernel density estimator. Theoretical study clarifies the asymptotic convergence rate of the plug-in type estimator is faster than the block-maxima-based estimator when $γ> -1$. A numerical comparative study on the bandwidth selection shows the performances of a plug-in approach and cross-validation approach depend on $γ$ and are totally comparable. Numerical study demonstrates that the plug-in nonparametric estimator with the estimated bandwidth by either approach overtakes the parametrically fitting estimator especially for distributions with $γ$ close to zero as $m$ gets large.