We present a finite difference algorithm for integrating the reactive flow model of contractile biological polymer networks on a fixed Eulerian mesh. We discuss the accuracy and limits of the algorithm. To illustrate the application of the algorithm, we carry out a family of computations involving an unreactive contractile network contained in a two-dimensional square reaction vessel. By numerical experiments using different values of the physical parameters of the model, we find that for this simple sort of system two major dynamical modes of contraction are predicted to occur. There is the squeezing type contraction in which the network contracts to a single small clump with gradual expulsion of solution material, and the rending type contraction in which the network tears itself into a number of separate pieces. We find that to a good approximation the transition between the squeezing mode and the rending mode is controlled by a single nondimensional number (the rending number). We discuss the relevance of these results for the analysis of various experimental observations.