What causes objects to become tidally locked? I'm trying to write a gravity simulation (suns planets etc), and was hoping tidal locking could be one feature demonstrated.

Using a simple equation for gravity has produced some interesting results, but (unless its emergent behaviour) I see nothing that would encourage tidal locking. But, after some reading it appears tidal locking is quite common, planets and their satellites, planets and suns, suns and other suns (binary stars).

Is it a result of the formation stage of these objects, or is it somehow a function of the equation of gravity?
Tidal locking occurs because the planet deforms the satellite into an oval, with long axis pointing towards the planet. If the satellite is rotating the long axis will move away from being pointing towards the planet, and the gravity of the planet will tend to pull it back, slowing the rotation until one face is permanently facing the planet. Tidal locking isn't a result of the formation processes, but a consequence of satellites not being perfectly rigid.

In order to model the effects of tides on the orbits and rotation periods of satellites you need to know several important pieces of information.

First you obviously need to know the size of the planet and the satellite (both in terms of mass and radius) the shape of the orbit and the rotation rate of both planet and satellite. For many objects, these values are well known.

Next, and this is the tricky bit, you need to know how the satellite and planet will be deformed by the other's gravity, and how much tidal heating will occur. These are the so-called "love number" (after Augustus Love) and the dissipation function, Q. 

It is hard to estimate these. For the Earth Moon system the ratio k/Q is known to be 0.0011. (but the Earth is a poor model for other planets, which don't have a substantial ocean, or a liquid core) 

For other planets the value of Q varies between 10 and 10000, with larger values for the gas giants, and k can be estimated from the rigidity of the bodies. 

A simple gravitation model is not able to capture the subtleties of the gravitational interaction between two mutually deforming bodies, indeed for most simulations, the planets are modelled as points, or at most as spheres, and this is good enough for all but the highest precision calculations.

Tidal locking takes a long time (by human standards) but a relatively short time compared with the age of the solar system. The time taken is very strongly dependent (order 6) on the radius of the orbit.

Direct simulation would be more or less impossible: the deformations are too small, and the time scale of locking is too large. It would be possible (though difficult) to model tidal locking in a simulation with unrealistic values for the rigidity of the satellite, and the size of the planet (think jelly world, orbiting a (Newtonian) black hole) so the deformation is greater and the locking time shorter. However modelling the  elastic deformation of a body under gravity is far from trivial.