The four-dimensional variational method (4DVar) is one of the data assimilation techniques to integrate numerical simulations and observation data, which especially plays an important role in the cases of large-scale simulation models such as weather forecasting. The conventional 4DVar estimates only the optimum initial conditions and/or parameters for the simulation to fit the data, but never evaluates their uncertainties. We developed a new algorithm for 4DVar that enables us to quantify the uncertainties by the adoption of the second-order adjoint method, which computes the product of a Hessian matrix and an arbitrary given vector without knowing the Hessian matrix. We introduce the foundation of 4DVar and our algorithm for the uncertainty quantification, demonstrating the application results in the cases of the phase-field model that simulates the time-evolution of grain growth in media, and a simple model that simulates seismic wave propagation in the underground.