
Notations
Numerically the integral becomes a sum over all the polygons in the
scene.
This is a system of equations where the unknowns are bi.
The Form Factor
If j corresponds to a light source, fi-j depends on the
angle between the ray of light from the source to the surface and the normal
to the surface.
If the objects i and j are both reflecting surfaces, then
fi-j depends on the angles made by the normals to the surface
in i and j with the line segment from i to j.
fi-i = 0, for i=1 to n.

Gauss-Seidel Method
Iteratively compute bi :
Initialize
Compute bi for the light sources based on the emitted energy.Repeat for a given number of iterations or until the change is not significant anymore.
For all the other polygons, assign bi = 0 (or sum the components from the light sources).
Recalculate each bi based on the previously computed b1, b2, ... bn.For examples and more information, see
If the scene is made of complex surfaces, discretize it as a mesh of polygons.