In this paper, we deal with the finite difference method for the initial boundary value problem of the nonlinear pseudo-parabolic system $(-1)^Mu_t+A(x,t,u,u_x,\cdots,u_x 2M-1)u_x2M_t=F(x,t,u,u_x,\cdots,u_x 2M)$,$u_xk(o,t)=\psi_{0k}(t), u_xk(L,t)=\psi_{1k}(t),k=0,1,\cdots,M-1,u(x,0)=\phi (x)$ in the rectangular domain $D=[0\leq X\leq L,0\leq t\leq T]$, where $u(x,t)=(u_1(x,t),u_2(x,t),\cdots,u_m(x,t)),\phi (x),\psi_{0k}(t),\psi_{1k}(t),F(x,t,u,u_x,\cdots,u_x 2M)$ are $m$-dimensional vector functions, and $A(x,t,u,u_x,\cdots,u_x2M-1)$ is an $m\times m$ positive definite matrix. The existence and uniqueness of solution for the finite difference system are proved by fixed-point theory. Stability, convergence and error estimates are derived.
