SciELO - Scientific Electronic Library Online

Home Pagelista alfabética de revistas  

Servicios Personalizados




Links relacionados


Maderas. Ciencia y tecnología

versión On-line ISSN 0718-221X

Maderas, Cienc. tecnol. v.8 n.3 Concepción  2006 

Maderas. Ciencia y tecnología. 8(3):149-158, 2006



Jarl-Gunnar Salin
SP Swedish National Testing and Research Institute Building Technology and Mechanics Wood Technology P.O.Box 5609, SE-114 86 Stockholm, Sweden

Corresponding author:

Received: 19.12.2005. Accepted: 14.04.2006.

The sapwood fibers in green softwood are normally filled – or nearly filled – with free water. This water forms a continuous phase through the narrow bordered pits between fiber lumens. In a drying process water evaporates from the menisci formed at these narrow throats and initially located close to the surface of the wood. As the amount of free water decreases in this way, a retraction of the meniscus will happen in the widest throat, which can be located anywhere, not necessarily close to the surface. As the individual size of these narrow throats is a stochastic variable, this process will produce complicated water phase structures that gradually fragmentize. If this drying behavior is viewed as such an invasion percolation process it is obvious that the result differs from normal diffusion, which however most often is used as the basis for wood drying modeling.

A computer based simulation model that utilizes the above principle has been developed. This model has been used for the simulation of a few situations of practical and theoretical importance. It has been found that this invasion percolation approach can reproduce some previously unexplained features seen in softwood drying processes. Some detailed examples are presented in this paper.

Keywords: wood drying, percolation, simulation, drying rate, dry shell


The fibers in the sapwood – i.e. the living part – of softwood species are, in the green state, normally filled with free water. Water in adjacent fibers is interconnected through the bordered pit openings in the cell wall. This free water will thus form a continuous (interconnected) phase covering a multitude of fibers. Such a group of fibers filled with water will be referred to as a water “cluster” in the wood. In a piece of sapwood, there will in the beginning – before drying starts – be only one such water cluster, or very few separate clusters at the most. When this piece of wood is placed in a dryer, then evaporation of water will occur from the outer surface of the cluster. The evaporation can take place from the free water surface (meniscus) in a bordered pit or in the open end of a cut fiber, but also as bound water diffusion through the cell wall into an empty fiber or into surrounding air. If there are empty fibers between the cluster and the wood piece surface, then there will be a vapor and bound water diffusion through this “dry” shell.

Anyway, as a result of this evaporation the volume of the free water in the cluster will decrease. The menisci in the openings at the outer border of the water cluster will move inwards and the radius of curvature of the free water surface will decrease. For a reasonably slow drying process, the suction pressure in the water phase, created by the curved menisci, will spread undamped throughout the cluster so that the same radius of curvature is established at all openings. Finally the widest opening cannot withstand the capillary suction pressure any longer and this water surface collapses into the corresponding fiber which is thus gradually emptied. This will uncover new bordered pit openings where new menisci are formed and then the process is repeated. It is important to notice that the location of the widest opening can, at each stage, be anywhere at the outer border of the water cluster and not necessarily close to the surface of the piece of wood, where presumably the highest evaporation is located. A stochastic element is thus introduced into the process and the water cluster will gradually shrink and also split into smaller clusters in a rather unpredictable way. It is obvious that this process differs substantially from a diffusion process which has been the normal way to describe wood drying also above the fiber saturation point (FSP).

The physical behavior described above can be generalized as a stochastic process in a regular network where changes occur at each grid point according to a rather simple mathematical algorithm (“widest opening” criterion in this case). Such processes are referred to as percolation, in this case especially invasion percolation, and these have been studied by physicists and mathematicians during the last two decades (Stauffer, Aharony, 1992). The invasion percolation approach has been used to analyze drying in previous works. Some that should be mentioned and that have developed pore network models of drying using the percolation concept are Prat (1993, 1995, 2000) Laurindo and Prat (1996, 1998) and Segura and Toledo (2005a, 2005b, 2005c) and references therein. Application of the concept for wood drying has however not been done, and it even seems that only non-hygroscopic media have been considered. This invasion percolation approach will now be used to describe the sapwood drying process above FSP and the results are analyzed in order to find out if this approach can reflect the experimentally observed behavior.


Because 3-dimensional networks are difficult to visualize, a 2-dimensional network is first used as an example of the invasion percolation model adapted to a drying process. Consider a (rectangular) area made up of a regular square network, Figure 1, where each square represents the cross section of a wood fiber. Each square has openings to each of the four neighboring squares and the width of these openings is a stochastic variable with a given distribution. There are in reality many openings between two adjacent fibers but we only need to consider the widest opening.

Figure 1. Two- dimensional square network

Assume that all squares in Figure 1 are filled with water and that the amount of water in each square (lumen volume) is a stochastic variable with a given distribution. Assume further that evaporation takes place from the upper horizontal border of this area into the surrounding air. Transfer of vapor or water through the vertical or bottom borders is not taken into account in this case. As a result of the evaporation, the square in the uppermost horizontal row which has the widest opening – for instance “A” in Figure 1 – will be emptied first. The time required for complete removal of its free water depends on the rate of evaporation (from the whole water cluster) and the lumen volume of “A”. As square A is emptied three new openings are uncovered as well as three new surfaces from which evaporation can take place. These three surfaces are however at a certain distance from the outer border, so there will be a diffusion of moisture through A towards the outer border. This transient moisture field is governed by the diffusion equation. However, solution of this partial differential equation is a heavy task and a simplified approach is preferred. Instead the steady-state field is calculated according to the Laplace equation. The moment a square is completely emptied a new solution which incorporates this square is recalculated. We are generally interested in relatively large networks and thus the Laplace equation is preferably solved numerically by finite differences using a grid directly based on the square network itself.

When the moisture field for the “dry” region is solved in this way, the individual fluxes from all border surfaces of the water cluster can be calculated and the sum of these fluxes equals the rate with which the next square will be emptied. The next “active” square is selected according to the “widest opening” criteria (in Figure 1 it would be B, C or D) a new steady state field is calculated for the “dry” region, and so on. At a certain point the water cluster will most likely split into two or more separate clusters. The same procedure continues of course for each cluster but it is now very important to calculate the time to the point of no free water in the “active” square for each cluster separately in order to establish in which order squares of different clusters are emptied. This requires an efficient method to determine which water filled squares that belong to the same water cluster. Fortunately, an algorithm exists – the so called Hoshen-Kopelman algorithm (Stauffer, Aharony, 1992) – which can handle this task rather efficiently.

Figure 2 shows the development of such an invasion percolation process on a 100 x 100 square network with evaporation from the upper horizontal surface. Blue color represents squares (fibers) filled with water and red color empty squares (no free water). Dark red color corresponds to a high relative humidity (high bound water content) and light red to a low RH. The very complex, almost fractal, structure seen in Figure 2 is typical. There is no well defined receding front which clearly indicates a deviation from traditional diffusion based models. A very typical feature is the “finger" that in the first picture penetrates deep into the area. This behavior can be explained in the following way. As squares with a wide opening are emptied, the distribution of the openings of remaining water filled squares is shifted towards smaller openings. But when new openings are uncovered, these obey the original distribution. There is thus a slightly increased probability that the next square emptied is a neighbor to the previous one. It can also be seen from Figure 2 that the number of separate clusters is very high at certain stages of the process.

Figure 2. Drying in a 100 x 100 two-dimensional square network

The nature of the statistical distributions for the width of openings between squares (fibers) and the volume of free water in each square (lumen volume) should be discussed. In Figure 2 both distributions where assumed to be uniform, i.e. the width varied uniformly between 0 and 1, and the volume between 0,5 and 1,5. As far as the width is concerned, the type of distribution does not matter in a homogeneous, isotropic network. The governing factor is the order of the values, not whether they are concentrated close to the average value or not. This order determines for each separate cluster how squares are emptied. The type of water volume distribution is on the contrary important. The volumes of the squares being emptied, determine which one is finished first and thus in which order separate clusters are selected in the process. In the case of wood drying it is thus important to determine the lumen volume distribution more precisely. Based on data found in (Atmer, Thörnqvist, 1982) a normal (Gauss) distribution with an average value of 0,0046 mm3 and a standard deviation of 50% seems adequate for both Scots pine (Pinus sylvestris) and Norway spruce (Picea sp.). When volume values are selected with a random number generator in the model, extreme values (for instance negative volumes) have to be rejected. This distribution has been used in the following simulations.


The 2-dimensional model described above and illustrated in Figure 2 is easily expanded to a 3-dimensional model using the same principles. If the squares in Figure 2 represent cross sections of fibers, then this 100 x 100 network corresponds to an area about 3 x 3 mm2 in size in reality. If the ultimate target is to simulate drying of timber, it is quite obvious that the number of fibers involved is too high to handle by this method. This means that only representative parts of the real solid can be handled. Then the question arises how the borders of this part should be described in the model. From Figure 2 it is quite clear that the vertical and the lower horizontal borders introduce “disturbances” into the simulation. In almost all cases of timber drying, evaporation takes place from at least two surfaces. If we consider the case of Figure 2 then evaporation from the lower horizontal surface removes the corresponding disturbance and gives a possibility to simulate drying of at least thin pieces. This approach has been used in the following 3-dimensional cases. The disturbances from the vertical surfaces can be reduced in two different ways. First, by ignoring squares close to the surface – for instance using only 60 central squares out of 100 in Figure 2 – then this problem is minimized. A second alternative method is to couple the left border to the right border (Figure 2). In this way an infinite – but periodic - network is obtained. The effect of this periodicity is difficult to determine.

Figure 3 illustrates two stages of a percolation process on a 3-dimensional 30 x 30 x 98 cubic network where evaporation takes place from the upper and lower horizontal surfaces into air (with the same climate at both ends). The red part in the middle shows cubes (fibers) that are emptied. The cubes that still are filled with free water are seen in the blue part on the right. To the left a horizontal cross section, at the height indicated, is shown. Generally the same features as found in the 2-dimensional case of Figure 2 are clearly seen, i.e. a very complicated structure and no well defined receding front. There seems in this case to be more emptied fibers in the upper part than in the lower, despite the symmetry.

The distance between the upper and lower surfaces would in reality be about 3 mm if this direction is seen as the radial direction. The length of the fiber (tracheid) is normally about 100 times longer than the diameter. In the model we present the fiber as a cube which means that the length scale in the longitudinal direction is considerably compressed. The piece of wood in Figure 3 can thus be interpreted as a 0,9 x 90 mm2 sheet of 3 mm thick veneer.


Figure 3: Drying in a 30x30x98 three-dimensional cubic network

It is now especially interesting to calculate the average moisture content at different heights as a function of time. For Scots pine (Pinus sylvestris)fully water filled lumens correspond to an MC of roughly 160% and the MC at the FSP is about 25%. Figure 4 is based on these values and shows the moisture profiles at different stages. The curves are equidistant in the sense that equally many fibers were emptied between curves. Each curve represents the average of ten different simulations and as the result should approach symmetry, the curves are actually made symmetric by taking the average of the left and right part.

Figure 4. Moisture profiles in the “radial” direction during drying of an initially completely water filled cubic network.

Two surprising properties are seen in Figure 4. In the beginning when “dry fingers” are penetrating into the body, the profiles are slightly rounded, but very soon profiles that are completely flat in the middle part are produced. This agrees qualitatively exactly with experimental data obtained by computer tomography scanning of sapwood during drying (Wiberg, Morén 1999a). This so-called “gradient-free” drying is seen in Figure 5. The second interesting feature is that these flat profiles seem to stop between 60 and 70% MC and receding fronts start to develop instead. In mathematical percolation theory the threshold probability is a well known concept. If sites in a very large (infinite) network are randomly occupied with a given probability then there will always be a continuous cluster from one side of the system to the other for probabilities above the threshold, and never for probabilities below the threshold. For cubic networks this threshold probability is 0.3116 (Stauffer, Aharony 1992) which converted to MC corresponds to 67,1% which is in very good agreement with Figure 4. This MC level is also in rather good agreement with experimental data reported by (Wiberg, Morén, 1999a). This threshold corresponds also to the concept of “irreducible saturation” introduced to describe the loss of continuity in the liquid phase.

Figure 5. Consecutive total density profiles in a pine sample (cross section 40 x 40 mm2) during drying, obtained by CT-scanning. (From Wiberg, Morén, 1999b).

As is well known, the sapwood of Scots pine (Pinus sylvestris) and Norway spruce (Picea sp.) is not fully saturated with free water even in the green state. To investigate this case, a certain number of cubes were randomly selected and emptied before the percolation algorithm was started. The MC distribution is presented in Figure 6. As seen, the flat profile feature starts immediately and the transit from flat profile to receding front is even more exactly defined than in Figure 4. One problem regarding this simulation should however be mentioned. At the starting point there will be several separate empty cubes that form gas “bubbles” in the structure. When a cube adjacent to the bubble is emptied the gas pressure in the bubble will decrease as the bubble expands. This suction will in reality slightly affect the “widest opening” mechanism and this is not accounted for in the calculations for Figure 6.

The drying rate is of course a very important variable that drying models should predict as accurately as possible. If sapwood drying is modeled as a diffusion process, then the model predicts that the surface is wet for a considerable time. This implies an initial phase with constant rate drying – as is found for many other materials. Experimental data for wood drying processes show however serious deviations from this general picture (Salin, 2002). The percolation model described here however predicts that a number of sites at the surface are emptied from the very beginning. Figure 7 presents the relative drying rate for the configuration seen in Figure 3. The solid curve is the average of several simulations. Curves of this type have previously been reported by for instance (Prat, 2002). The solid curve of Figure 7 drops initially very rapidly, not a horizontal line as for a diffusion based situation. This drop is so rapid that it is perhaps not always seen in practice due to uncontrolled pre-drying of the surface before the wood (timber) is loaded into the dryer.

Figure 6. Moisture profiles in the “radial” direction during drying of an initially partly water filled cubic network.

After this phase an almost constant drying rate follows. It is possible that some researchers have interpreted this as a real constant drying rate phase and concluded that the surface has been wet. A constant rate does not, however, prove that the surface is wet, although the opposite connection is correct. This almost constant drying rate ends when the receding front (see above) phase starts and the rate drops then gradually. The curve of Figure 7 is in qualitative agreement with correction factors that have been introduced into diffusion based wood drying models in order to make them fit experimental data better. It is obvious that the percolation concept gives a basis for a natural explanation of these correction factors.

Figure 7. Relative drying rates corresponding to the moisture profiles of Figure 4 (solid line) and Figure 8 (dots).


So far some important properties of the wood structure have not been accounted for in the percolation model above. For instance, isotropic properties have been assumed although it is well known that properties are different in the three principal directions – radial, tangential and longitudinal. Different diffusion coefficients in different directions have thus been introduced in the “dry” region by using the values reported by (Olek et al, 2005). The mass transfer rate at the horizontal surfaces is defined by a surface emission coefficient, with a numerical value corresponding to normal timber drying. Further isothermal conditions within the solid are assumed.. One important feature in this context is the fact that bordered pit openings are concentrated to the radial surfaces of the tracheids. This makes free water flow in the tangential (and longitudinal) direction more easily than in the radial direction. In our specific case, however, the number of openings is not important, only the width of the widest opening. In order to investigate this feature, the configuration of Figure 3 has been repeated but with no openings in the tangential (horizontal in Figure 3) surface between fibers for 90% of these sites. The remaining randomly selected 10% are assumed to obey the same statistical width distribution as in the radial and longitudinal directions. The network is further extended to a 60 x 60 x 98 network. The calculated MC distribution is shown in Figure 8. In general the moisture profiles of Figure 8 resemble those of Figure 4 but the “flat profile” phase is shorter. The receding front phase starts earlier which is natural as the partial absence of openings in the tangential direction will break the continuity in the water phase at an earlier stage.

Figure 8. Moisture profiles in the radial direction during drying when there are no openings between fibers in the radial direction in 90% of the cases (60 x 60 x 98 network).

The drying rate corresponding to the profiles of Figure 8 is presented in Figure 7 as dots. A clear constant drying rate phase is seen but closer to unity than for the isotropic case. This indicates that the “dry” shell at the surface is developing slower in this particular case, but the rate depends on several parameters. Figures 2 and 3 show clearly that there is no well defined dry shell, but it is of course possible to calculate an “equivalent” dry shell thickness from the drying rate results. The drying rate results of Figure 7 are converted into equivalent dry shell thickness in Figure 9. The conclusion from Figure 9 is that the dry shell thickness is very small for a considerable time. It is thus no surprise that it has not been directly observed, except in a few cases. Despite that, the effect on the drying rate is clearly seen almost from the beginning of the drying process. This explains why the measured mass transfer coefficient at the wood surface has appeared to be lower than expected from classical heat and mass transfer theory.

Figure 9. Thickness of the equivalent dry shell corresponding to the moisture profiles of Figure 4 (solid line) and Figure 8 (dots).


The invasion percolation approach seems to be a promising method to investigate pore level behavior related to drying of free water in the sapwood of softwoods. Several poorly understood phenomena can obviously be explained by this approach. It is however also clear that a more accurate adaptation to the specific properties of wooden microstructures, than done in this survey, is needed. A hexagonal-rectangular fiber structure according to the Comstock model – instead of the cuboids in this presentation – has been introduced in another investigation, as well as the effect of ray cells (Salin 2006). Further an analysis of the “kiln brown stain” mechanism has been performed.


The work presented in this paper has been financed by the Swedish Research Council for Environment, Agricultural Sciences and Spatial Planning, which is gratefully acknowledged.

This paper was originally presented at the 3rd Nordic Drying Conference in Karlstad, Sweden, June 15-17, 2005 but it has been peer reviewed for Maderas. Ciencia y tecnologia journal.


Atmer, B and Thörnqvist, T. 1982. The properties of tracheids in spruce and pine. (In Swedish). The Swedish University of Agricultural Sciences, Department of Forest Products, Uppsala. Report No 134.        [ Links ]

Laurindo, J.B and Prat, M. 1996. Numerical and experimental network study of evaporation in capillary porous media. Phase distributions. Chem. Eng. Sci. 51(23):5171-5185.        [ Links ]

Laurindo, J.B and Prat, M. 1998. Numerical and experimental network study of evaporation in capillary porous media. Drying rates. Chem. Eng. Sci. 53(12): 2257-2269.        [ Links ]

Olek, W.; Perré, P.; Weres, J. 2005. Inverse analysis of the transient bound water diffusion in wood. Holzforschung 59:38-45.        [ Links ]

Prat, M. 1993. Percolation model of drying under isothermal conditions in porous media. Int. J. Multiphase Flow 19(4):691-704.        [ Links ]

Prat, M. 1995. Isothermal drying of non-hygroscopic capillary-porous materials as an invasion percolation process. Int. J. Multiphase Flow 21 (5):875-892.        [ Links ]Prat, M. 2002. Recent advances in pore-scale models for drying of porous media. Chem. Eng. J. 86,153.         [ Links ]

Salin, J-G. 2002. Theoretical analysis of mass transfer from wooden surfaces. Proceedings of the 13th International Drying Symposium (IDS2002). Beijing, China, Aug. 27-30. pp.1826-1834.        [ Links ]

Salin, J-G. 2006. Modelling of the behaviour of free water in sapwood during drying. Part I. A new approach. Accepted for publication in Wood Material Science and Engineering vol.1(1)

Segura, L.A and Toledo, P.G. 2005a. Pore-level modeling of isothermal drying of pore networks. Evaporation and viscous flow. Latin American Applied Research 35:43-50.        [ Links ]

Segura, L.A and Toledo, P.G. 2005b. Pore-level modeling of isothermal drying of pore networks accounting for evaporation, viscous flow, and shrinking. Drying Technology 23:2007-2019.        [ Links ]

Segura, L.A and Toledo, P.G. 2005c. Pore-level modeling of isothermal drying of pore networks. Effect of gravity and pore shape and size distribution on saturation and transport parameters. Chem. Eng. J. 111:237-252.        [ Links ]

Stauffer, D and Aharony, A. 1992. Introduction to Percolation Theory, Taylor & Francis Ltd, London.        [ Links ]

Wiberg, P and Morén, T. 1999a. Moisture flux determination in wood during drying above fibre saturation point using CT-scanning and digital image processing. Holz als Roh und Werkstoff 57:137-144.        [ Links ]

Wiberg, P and Morén, T. 1999b. New information on sapwood drying; CT-scanning profiles of moisture content during drying. First COST Action E15 Wood Drying Workshop, Edinburgh, UK, Oct.13-14.        [ Links ]

Creative Commons License Todo el contenido de esta revista, excepto dónde está identificado, está bajo una Licencia Creative Commons