One dimensional non-steady, non-Darcy flow of water in a rock stratum was reduced into a system described by six ordinary differential equations involving five controlling parameters. Through response computations and time series analysis, chaotic behavior in the reduced system was discussed in details. Firstly, the dynamical response of the reduced system under a set of parameters was calculated, and the power spectrum of the attractor was obtained through fast Lagrangian transformation; then the phase space was reconstructed by fixing embedding dimension to be 6 and delay time to range from 1 to 20, and the correlation dimension of the attractor was calculated based on the curves under the coordinates of logarithm of correlation integral vs. logarithm of covering radius; and lastly, the Lyapunov indices of the attractor were calculated by using Gram-Schmit’s orthogonalization method. The results show that the power spectrum of the attractor is continuous; the correlation dimension of the attractor is equal to 2.36; among the Lyapunov indices, LE1, LE2, LE3 are positive, LE5, LE6 are negative, and LE4 fluctuates near zero. All the analysis indicates that there may exist chaos in the system of non-steady, non-Darcy flow. One dimensional non-steady, non-Darcy flow of water in a rock stratum was reduced into a system described by six ordinary differential equations involving five controlling parameters. Through response computations and time series analysis, chaotic behavior in the reduced system was discussed in details Firstly, the dynamical response of the reduced system under a set of parameters was calculated, and the power spectrum of the attractor was obtained through fast Lagrangian transformation; then the phase space was reconstructed by fixing embedding dimension to be 6 and delay time to range from 1 to 20, and the correlation dimension of the attractor were calculated based on the curves under the coordinates of logarithm of correlation integral vs. logarithm of covering radius; and lastly, the Lyapunov indices of the attractor were calculated by using Gram-Schmit’s orthogonalization method. The results show that the power spectrum of the attractor is continuous; the correlation dimension of LE1, LE2, LE3 are positive, LE5, LE6 are negative, and LE4 fluctuates near zero. All the analysis indicates that there may exist chaos in the system of non-steady, non -Darcy flow.