Higher order radiosity approximation

 

The original radiosity method is based on finite element techniques. In other words, the radiosity distribution is searched in a piecewise constant function form, reducing the original problem to the calculation of the values of the steps.

The idea of piecewise constant approximation is theoretically simple and easy to accomplish, but an accurate solution would require a large number of steps, making the solution of the linear equation difficult. Besides, the constant approximation can introduce unexpected artifacts in the picture even if it is softened by Gouraud shading.  

This section addresses this problem by applying a variational method for the solution of the integral equation [SK93].

The variational solution consists of the following steps [Mih70]:  

1. It establishes a functional which is extreme for a function (radiosity distribution) if and only if the function satisfies the original integral equation (the basic radiosity equation).
2. It generates the extreme solution of the functional by  Ritz's method, that is, it approximates the function to be found by a function series, where the coefficients are unknown parameters, and the extremum is calculated by making the partial derivatives of the functional (which is a function of the unknown coefficients) equal to zero. This results in a linear equation which is solved for the coefficients defining the radiosity distribution function.

Note the similarities between the second step and the original radiosity method. The proposed variational method can, in fact, be regarded as a generalization of the finite element method, and, as we shall see, it contains that method if the basis functions of the function series are selected as piecewise constant functions being equal to zero except for a small portion of the surfaces. Nevertheless, we are not restricted to these basis functions, and can select other function bases, which can approximate the radiosity distribution more accurately and by fewer basis functions, resulting in a better solution and requiring the calculation of a significantly smaller linear equation.

  
Figure: Geometry of the radiosity calculation

Let the diffuse coefficient be at point p and the visibility indicator between points p and p' be H(p,p'). Using the notations of figure 1.9, and denoting the radiosity and emission at point p by B(p) and E(p) respectively, the basic radiosity equation is:

where f(p,p') is the point-to-point form factor:

 

Dividing both sides by dA, the radiosity equation is then:

 

Let us define a linear operator :

 

Then the radiosity equation can also be written as follows:

 

The solution of the radiosity problem means to find a function B satisfying this equation. The domain of possible functions can obviously be restricted to functions whose square has finite integration over surface A. This function space is usually called space where the scalar product is defined as:

If were a symmetric and positive operator, that is, for any u,v in ,

were an identity and

then according to the minimal theorem of quadratic functionals [Ode76] the solution of equation 1.62 could also be found as the stationary point of the following functional:

 

Note that makes no difference in the stationary point, since it does not depend on B, but it simplifies the resulting formula.

To prove that if and only if some satisfies

for a symmetric and positive operator , then is extreme for the assumption that is positive and symmetric can be used:

Since only the term depends on B and this term is minimal if and only if is zero due to the assumption that is positive, therefore the functional is really extreme for that which satisfies equation 1.62.

Unfortunately is not symmetric in its original form (equation 1.61) due to the asymmetry of the radiosity equation which depends on but not on . One possible approach to this problem is the subdivision of surfaces into finite patches having constant diffuse coefficients, and working with multi-variate functionals, but this results in a significant computational overhead.

Now another solution is proposed that eliminates the asymmetry by calculating B(p) indirectly through the generation of . In order to do this, both sides of the radiosity equation are divided by :

Let us define , and g(p,p') by the following formulae

 

Using these definitions, we get the following form of the original radiosity equation:

 

Since g(p,p') = g(p',p), this integral equation is defined by a symmetric linear operator :

As can easily be proven, operator is not only symmetric but also positive taking into account that for physically correct models:

 

This means that the solution of the modified radiosity equation is equivalent to finding the stationary point of the following functional:

This extreme property of functional I can also be proven by generating the functional's first variation and making it equal to zero:

Using elementary derivation rules and taking into account the following symmetry relation:

the formula of the first variation is transformed to:

The term closed in brackets should be zero to make the expression zero for any variation. That is exactly the original radiosity equation, hence finding the stationary point of functional I is really equivalent to solving integral equation 1.71.

In order to find the extremum of functional , Ritz's method is used. Assume that the unknown function is approximated by a function series:

where form a complete function system (that is, any piecewise continuous function can be approximated by their linear combination), and are unknown coefficients. This assumption makes functional an n-variate function , which is extreme if all the partial derivatives are zero. Having made every equal to zero, a linear equation system can be derived for the unknown -s ( ):

 

This general formula provides a linear equation for any kind of complete function system , thus it can be regarded as a basis of many different radiosity approximation techniques, because the different selection of basis functions, , results in different methods of determining the radiosity distribution.

Three types of function bases are discussed:

• piecewise constant functions which lead to the traditional method, proving that the original approach is a special case of this general framework,
• piecewise linear functions which, as we shall see, are not more difficult than the piecewise constant approximations, but they can provide more accurate solutions. It is, in fact, a refined version of the method of `` vertex-surface form factors'',  
• harmonic (cosine) functions where the basis functions are not of finite element type because they can approximate the radiosity element type because they can approximate the domain, and thus fall into the category of global element methods.  

 
Figure: One-dimensional analogy of proposed basis functions

• Piecewise constant radiosity approximation
• Linear finite element techniques
• Global element approach - harmonic functions