A constitutive relationship of Kelvin model’s type is developed by sequentially combing a Maxwell component and any number of Kelvin components. Through compiling the head and source files of generalized Kelvin model by using C++ on Microsoft Visual Studio platform, the dynamic link library (DLL) file of the generalized Kelvin constitutive model is developed, and a user-defined generalized Kelvin model is implemented by calling the DLL file in two-dimensional Particle Flow Code (PFC2D). Through stress relaxation tests, in which different numbers of Kelvin models are included in the generalized Kelvin model, the relationship between contact force and time is obtained, which turns out to coincide with the analytical solution. Based on these research results, the programmed generalized Kelvin model is used to simulate the failure process of a silty clay slope of an open pit coal mine. The numerical results illustrate that the failure surface is brim-shaped, whereas the elastic model does not predict such a failure mode, and instead it yields a circular sliding surface, implying that the developed generalized Kelvin model can be used to study the failure mechanism of silty clay slope.