Practice and Matrix

 

Design Issues

The important thing to do are to encapsulate the data needed for the solution. Also, encapsulating the computation of the transfer (basis coefficients) is very convenient.

class RadObj
Members dealing with the meshing and hierarchy bool Subdivide()
int GetNumChildren()
RadObj *GetChild(int) Members dealing with getting geometric data needed for doing the projection. Often corners are needed, which implies here that U and V are in the range [0,1]. GetSample(UV, Point, Normal)
float GetArea()
Members needed for actually updating the radiance of the surface. RadData holds the solution data and other parts of the future linear system. Spectrum GetRho()
RadData *RetRadData()

The RadData class holds the coefficients of the basis functions for both the computed radiance and the emission. In the future it will hold data necessary to represent the matrix. class RadData
GetRadCoefficient(s)()
GetEmissionCoefficient(s)()

The ComputeProjection function computes the coefficients on the reciever due to the environment. In the future, it will take a specific source as well as the reciever. ComputeProjection(Rec)

For an example of what I am using, see here. I don't recommend you use this, or even follow it closely, as it includes a lot more generality than you need now, and probably different generality than you would want after this class is over. The generality comes at a large cost in complexity, and I'm not even sure it's worth it for my code yet.

Creating the Radiosity Matrix

Read CW [10] pp 48-60.

First, intuitively. We have

$$\widehat{L} \sum_{i} c_i B_i$$ = (10)

$$\left< \mbox{$\widetilde{B_i}$ } \, , \, L_e+ ML \right>$$ c_i = (11)

Unfortunately, we don't know L exactly, so why not substitute $\widehat{L}$ for L? Doing this gives

$$% c_i \approx \left< \mbox{$\widetilde{B_i}$ } \, , \, L_e+ M\widehat{L} \right> \\$$

$$% c_i \approx \left< \mbox{$\widetilde{B_i}$ } \, , \, L_e+ M\left(\sum_{j} c_j B_j\right) \right>$$

$$\left< \mbox{$\widetilde{B_i}$ } \, , \, L_e \right> + \left< \mbox{$\widetilde{B_i}$ } \, , \, M\left(\sum_{j} c_j B_j\right) \right>$$ = (12)

$$\left< \mbox{$\widetilde{B_i}$ } \, , \, L_e \right> + \sum_j c_j \left< \mbox{$\widetilde{B_i}$ } \, , \, MB_j \right>$$ = (13)

$$\widehat{E}_i + \sum_j c_j \left< \mbox{$\widetilde{B_i}$ } \, , \, MB_j \right>$$ = (14)

We now have a system of equations that can be represented by

$$\widehat{L}= (I - \widehat{M}) \widehat{E}.$$

where

$$\widehat{M}_{ij} = \left< \mbox{$\widetilde{B_i}$ } \, , \, MB_j \right>$$

The more correct way to derive this is to look at something called the residual. What we want to do is minimize

$$\left\vert \left\vert L- \widehat{L} \right\vert \right\vert$$

except we have no information about L. We do know that if

$$\widehat{L}= L_e+ M\widehat{L}$$

then $\widehat{L}$ is the solution. In general it won't be equal. We define the residual as the difference between the two sides.

$$r(x) = \widehat{L}- L_e+ M\widehat{L}.$$

We can't make $r \equiv 0$, but we can try to make the projection of r onto our basis function be 0.

$$\left< r \, , \, \mbox{$\widetilde{B_i}$ } \right> = 0$$

Working this out gives the same matrix expression as above.

The choice of inner products (duals) determines the type of solution method you are doing. Delta function duals create a solution where the residual function is zero at the sample points. This is called point collocation. When the duals span the same space as the basis functions (e.g. piecewise constant) then the approach is called the Galerkin method.

More Reading

Arvo works this out in all the gory details[3].