We perform hybrid Monte Carlo simulation of (2+1+1)-flavors lattice QCD with the optimal domain-wall fermion (which has the effective 4D Dirac operator exactly equal to the Zolotarev optimal rational approximation of the overlap Dirac operator). The gauge ensemble is generated on the 323×64 lattice with the extent Ns=16 in the fifth dimension, and with the plaquette gauge action at β=6/g2=6.20. The lattice spacing (a≃0.063 fm) is determined by the Wilson flow, using the value t0=0.1416(8) fm obtained by the MILC Collaboration for the (2+1+1)-flavors QCD. The masses of s and c quarks are fixed by the masses of the vector mesons ϕ(1020) and J/ψ(3097) respectively; while the mass of the u/d quarks is heavier than their physical values, with the unitary pion mass Mπ≃280 MeV (and MπL≃3). We compute the point-to-point quark propagators, and measure the time-correlation functions of meson and baryon interpolators. Our results of the mass spectra of the lowest-lying hadrons containing s and c quarks are in good agreement with the high energy experimental values, together with the predictions of the charmed baryons which have not been observed in experiments.