Radial Basis Functions Methods for Solving Radionuclide Migration, Phase Change and Wood Charring Problems

Modelling the flow through porous media has a great importance for solving the problems of disposal of radioactive waste. When modelling the flow of contaminated material through the geosphere, it is important to consider all internal processes (e.g. advection, dispersion, retardation) within the geosphere, and external processes associated with the near-field and the biosphere. The general reliability and accuracy of transport modelling depend predominantly on input data such as hydraulic conductivity, water velocity on the boundary, radioactive inventory, hydrodynamic dispersion. The output data are concentration, pressure, etc. The most important input data are obtained from field measurement, which are not available for all regions of interest. In such cases, geostatistical science offers a variety of spatial estimation procedures. A vast variety of important physical processes involving heat conduction and materials undergoing a change of phase may be approached as Stefan problems. One of these processes is the heat transfer involving phase changes caused by solidification or melting, which are important in many industrial applications such as the drilling of high ice-content soil, the storage of thermal energy and the safety studies of nuclear reactors. Due to their wide range of applications, the phase change problems have drawn considerable attention of specialists in different fields of science and engineering. Another problem which seems to be completely different but is in mathematical terms very similar to solidification or melting is charring of wood. After wood is exposed to fire it undergoes thermal degradation. The pyrolysis gases undergo flaming combustion as they leave the charred wood surface. The pyrolysis, charring, and combustion of wood have been presented by (Fredlund, 1993) who performed experiments and numerical analyses. For all physical processesmentioned above, themotion of fluids, phase changes, and pyrolysis processes are governed by a set of partial differential equations (PDES). These governing equations are based upon the fundamental conservation laws. The mass, momentum and energy are conserved in any fluid motions. In most cases, the governing equations are too Radial Basis Functions Methods for Solving Radionuclide Migration, Phase Change and Wood Charring Problems


Introduction
Modelling the flow through porous media has a great importance for solving the problems of disposal of radioactive waste.When modelling the flow of contaminated material through the geosphere, it is important to consider all internal processes (e.g.advection, dispersion, retardation) within the geosphere, and external processes associated with the near-field and the biosphere.The general reliability and accuracy of transport modelling depend predominantly on input data such as hydraulic conductivity, water velocity on the boundary, radioactive inventory, hydrodynamic dispersion.The output data are concentration, pressure, etc.The most important input data are obtained from field measurement, which are not available for all regions of interest.In such cases, geostatistical science offers a variety of spatial estimation procedures.A vast variety of important physical processes involving heat conduction and materials undergoing a change of phase may be approached as Stefan problems.One of these processes is the heat transfer involving phase changes caused by solidification or melting, which are important in many industrial applications such as the drilling of high ice-content soil, the storage of thermal energy and the safety studies of nuclear reactors.Due to their wide range of applications, the phase change problems have drawn considerable attention of specialists in different fields of science and engineering.Another problem which seems to be completely different but is in mathematical terms very similar to solidification or melting is charring of wood.After wood is exposed to fire it undergoes thermal degradation.The pyrolysis gases undergo flaming combustion as they leave the charred wood surface.The pyrolysis, charring, and combustion of wood have been presented by (Fredlund, 1993) who performed experiments and numerical analyses.For all physical processes mentioned above, the motion of fluids, phase changes, and pyrolysis processes are governed by a set of partial differential equations (PDES).These governing equations are based upon the fundamental conservation laws.The mass, momentum and energy are conserved in any fluid motions.In most cases, the governing equations are too

Radial Basis Functions Methods for Solving
Radionuclide Migration, Phase Change and Wood Charring Problems Leopold Vrankar 1 , Franc Runovc 1 and Goran Turk 2 complex to be solved analytically.Therefore the numerical methods have been extensively used to find an approximate solution of the equations.
Traditionally, the most popular methods have been the finite element methods (FEM), the finite difference methods (FDM), and boundary element method (BEM).In recent years, the radial basis functions (RBF) methods have emerged as a novel computing method in scientific computing community.All the conventional methods mentioned above can be considered as mesh-based methods characterized by their reliance on a computation mesh with certain relationship between the nodes.In spite of their great success in solving scientific and engineering problems over the past five decades, these conventional numerical methods still have some drawbacks that impair their computational efficiency and even limit their applicability to more practical problems, particularly in three-dimensional space.The RBF methods could overcome some of these drawbacks by constructing the approximations entirely in terms of a set of points.The methods can be extended to multi-dimensional problems without significant effort.Over the last 30 years, many researchers have shown a great deal of interest in RBFs.It was used for modelling radionuclide migration (Vrankar et al., 2004a), (Vrankar et al., 2005), structural topology optimization (WangS&W angMY,2006) and many other applications (Ling, 2003).
The chapter is organized as follows.In the first section a short description of geostatistics is given.The sections on radial basis functions follow.The efficiency of the presented method is explained on three very different cases.The first example presents the usage of RBFs in geostatistical analysis for modelling of the radionuclide migration.The advection-dispersion equation is used in either Eulerian or Lagrangian form.
The second case presents a phase change problem.The moving interface is captured by the level set method at all time with the zero contour of a smooth function known as the level set function.A new approach is used to solve a convective transport equation for advancing the level set function in time.This new approach is based on the asymmetric meshless collocation method and the adaptive combination of RBFs method and Greedy algorithm for trial subspaces selection.
The last case explores the char formation in the wood as a function of surrounding temperature.The problem is solved numerically by the radial base function (RBF) methods.
The results are tested on the one-dimensional case in the standard fire conditions, and the same model is used to analyse a two-dimensional behaviour of timber beam exposed to fire from three sides.At the end of the chapter some conclusions are given.

Geostatistics
The term geostatistics is employed here as a generic term, meaning the application of the theory of random fields in the earth sciences (Kitanidis & VoMvoris, 1983).The parameters are distributed in space and can thus be called regionalized variables.The parameters of a given geologic formation can conveniently be represented as realisations of random variables which form random fields.Stochastic simulation is a widely accepted tool in various areas of geostatistics.Simulations are termed globally accurate through the reproduction of one-, two-, or multiple-point statistics representative of the area under study.
The most convenient method for simulation of random fields is sequential Gaussian simulation (Deutsch & Journel, 1998) because all successive conditional distributions from which simulated values are drawn are Gaussian with parameters determined by the solution of a simple kriging system.Sequential Gaussian simulation procedure: 1. First, use a sequential Gaussian simulation to transform the data into a normal distribution.
2. Perform variogram modelling on the data.Select one grid node at random, then krige the value at that location.This also gives the kriged variance.
3. Draw a random number from a normal distribution that has a variance equivalent to the kriged variance and a mean equivalent to the kriged value.This number is the simulated number for that grid node.
4. Select another grid node at random and repeat.For the kriging, include all the previously simulated nodes to preserve the spatial variability as modelled in the variogram.
5. When all nodes have been simulated, back transform to the original distribution.This gives us first realization using a different random number sequence to generate multiple realizations of the map.
Kriging (D.G. Krige, a South African mining engineer) is a collection of generalized linear regression techniques for minimizing an estimation variance defined from a prior model for a covariance (semivariogram) (Olea, 1991).One of the basic statistical measures of geostatistics is the semivariance, which is used to express the rate of chance of a regionalized variable along a specific orientation (or a measure of degree of spatial dependence between samples along a specific support) (Davis, 1986).If we calculate the semivariances for different values of h,we can plot the results in the form of a semivariogram: where Z(x) and h are the realization of random process in the space and norm of vector h = ( h ).Since the semivariogram is a function of distance, the weights change according to the geographic arrangement of the samples.Kriging can be used to make contour maps, but unlike conventional contouring algorithms, it has certain statistically optimal properties.

Meshless methods
So far, research on the numerical method has focused on the idea of using a meshless methodology for the numerical solution of PDEs.One of the common characteristics of all mesh-free methods is their ability to construct functional approximation or interpolation entirely from information at a set of scattered nodes, among which there is no relationship or connectivity needed: Three approaches to meshless methods have been successfully proposed.The first one is based on the finite element method and employs Petrov-Galerkin weak formulation (Atluri & Shen, 2002).The second one is of boundary element type (Li & Aluru, 2002).The third approach employs the RBFs.The base of this approach is its employment of different interpolating functions to approximate solutions of differential equations.Kansa (1990a) introduced multiquadric functions to solve hyperbolic, parabolic and elliptic differential equations with collocation methods.One of the most powerful RBF methods is based on multiquadric basis functions (MQ), first used by (Hardy, 1971).

139
Radial Basis Functions Methods for Solving Radionuclide Migration, Phase Change and Wood Charring Problems  (Wang & Liu, 2002).It is also important to mention the local version of the RBF collocation which does not produce the ill conditioning and is less sensitive to the shape parameter selection.The local version of the RBF collocation is based on compactly supported RBF spline: (ϕ(r)=(1 − r) m + p(r)) where p(r) is a polynomial of the Wendland (Wendland, 1995) compactly supported (CS-RBF) spline.The parameter c > 0 is a shape parameter controlling the fitting of a smoothing surface to the data.It has a significant influence on the accuracy of the solution.For this reason, in almost all previous researches the shape parameter must be adjusted with the number of centers in order to produce equation systems that are sufficiently well conditioned to be solved with the standard finite precision arithmetic.The optimal choice of the constant shape parameter is still an open question, and it is most often selected by the trial and error approach.
To overcome the problems of ill-conditioned matrix many efforts have been made to find new computational methods being capable of avoiding the ill-conditioning problems using linear solvers.In the literature, the following methods are reported: • (a) Using variable shape parameters (Kansa, 1990b), • (b) preconditioning the coefficient matrix, see Ling and Kansa (Ling & Kansa, 2004), 140 Advanced Methods for Practical Applications in Fluid Mechanics www.intechopen.com• (c) using domain decomposition methods in overlapping or non-overlapping schemes that decompose a very large ill-conditioned problem into many subproblems with better conditioning (Kansa, 1990a), • (d) optimizing the center locations by the Greedy algorithm.
Combination of the above can be used as optimized solution procedure and accuracy.

Asymmetric meshless collocation methods for stationary problems
We briefly review the RBF asymmetric collocation scheme.We consider a PDE in the general form of where u is the unknown solution, d denotes the dimension, ∂Ω is the boundary of the domain Ω, L is the differential operator on the interior, and B is the operator that specifies the boundary conditions of the Dirichlet, Neumann or mixed type.Both, f and g, are given functions with sufficient smoothness mapping R d → R. Using the Kansa's asymmetric multiquadric collocation method, the unknown PDE solution u is approximated by RBFs in the form: where ϕ is any type of radial basis function, and ν l , ..., ν M ∈ ∏ d q , are polynomials of degree m or less, Q := q − 1 + d d , see (Iske, 1994).Let (x j ) N j=1 be the N collocation points in Ω ∂Ω.We assume the collocation points are arranged in such a way that the first N I points are in Ω, whereas the last N B points are on ∂Ω.To evaluate or approximate the N + Q unknown coefficients, at least, N + Q linearly independent equations are needed.Ensuring that U(x) satisfies equations ( 2) and (3) at the collocation points results in a good approximation of the solution u.The first N equations are given by and the others are given as 141 Radial Basis Functions Methods for Solving Radionuclide Migration, Phase Change and Wood Charring Problems The last Q equations could be obtained by imposing additional conditions on p(•): The above procedure leads to the system of equations where The RBF method for solving PDEs is straightforward and easy to implement.Since the RBF is differentiable, spatial derivatives are computed by simply differentiating the MQ-RBFs; time derivatives are computed by differentiating the time dependent expansions coefficients.The spatial and temporal partial derivatives result in either a linear or nonlinear set of ordinary differential equations.

Modelling of the radionuclide migration
The central issue in modelling is on the one hand consistency between conceptual and mathematical models and, on the other hand between conceptual models and scenarios.
A conceptual model is a qualitative description of the system functioning in a form which corresponds to mathematical representation.Each scenario is a set of features, processes and events which has to be considered together to assess the impact of the disposal in the future.Groundwater models are presented by motion and continuity equations.The majority of the codes currently used or under development are based on the advective-dispersive equation (Bear, 1972)

Eulerian form of the advection-dispersion equation
The first step of radionuclide transport modelling is to obtain the Darcy velocity, here denoted by V.The continuity equation (e.g. for fluid phase) can be written 142 Advanced Methods for Practical Applications in Fluid Mechanics www.intechopen.com where ̺ f is the density of fluid.If the porosity ǫ is constant, equation ( 14) is simplified to With assumed constant density, Equation ( 15) simplifies to In saturated porous media the fluid flow is described by the Darcy equation where ∇p is the pressure gradient, K is hydraulic conductivity tensor, and μ is the dynamic viscosity.If we assume that the hydraulic conductivity is constant in each point of the considered area, then the combination of equations ( 16) and ( 17) gives a type of Laplace equation.In cases of homogeneous or heterogeneous and anisotropic porous media and incompressible fluid, the calculation of the velocities in principal directions which were determined from the pressure of the fluid obtained from the Laplace differential equation is presented in (Vrankar et al., 2004a), (Vrankar et al., 2005).The Lagrangian form of the advection-dispersion equation is presented in (Vrankar et al., 2004b) The velocities obtained from Laplace equation are used in the advection-dispersion equation.
The advection-dispersion equation for transport through the saturated porous media zone with retardation and decay is: where x is the groundwater flow axis, y is the transverse axis, u is the concentration of contaminant in the groundwater [Bqm −3 ], D x and D y are the components of dispersion tensor (two processes are incorporated: molecular diffusion and dispersion.Molecular diffusion is a thermo-chemical process, where mass is transported due to thermal or solutes gradients.On the other hand, dispersion is mechanical process, where spreading of the substance is caused due to the motion of the fluid.Details you can found in (Bear, 1972)) [m 2 y −1 ] in the saturated zone, ω is porosity of the saturated zone[−], v x is Darcy velocity [my −1 ] at interior points, R is the retardation factor in the saturated zone [−] and λ is the radioactive decay constant [y −1 ].
In these cases [y] means years.
For the parabolic problem, we consider the implicit scheme: 143 Radial Basis Functions Methods for Solving Radionuclide Migration, Phase Change and Wood Charring Problems where δt is the time step and u n and u n+1 are the contaminant concentrations at the time t n and t n+1 .
The approximate solution is expressed as : where α n+1 j , j = 1, ..., N are the unknown coefficients to be determined.ϕ j (x, y) is Hardy's multiquadrics function Hardy (1971): where c is shape parameter.
By substituting (20) into (18), we have: from which we can solve the N × N linear system of ( 22)-( 23) for the unknown coefficients α n+1 j , j = 1, ..., N. Let N = N I + N B be the number of collocation points, N I is the number of interior points and N B is the number of boundary points.Then (20) gives the approximate solution at any point in the domain Ω.The traditional finite difference method (FDM) was also used for solving the Laplace and advection-dispersion equation.For the approximation of the first derivative second-order central difference or one-sided difference were used.But for the approximation of the second derivatives we used the second-order central difference.

2D Solution
The simulation was performed on a rectangular area of 450 by 300 m.The source (initial condition) was Thorium (Th − 230) with activity 1 • 10 6 Bq and half life of 77000 years.The distribution of hydraulic conductivity (values between 66.00 and 83.41 [ m y ]) for one specific simulation is shown in Fig. 2. The groundwater flow field is presented for steady-state conditions.Except for the inflow (left side) and outflow (right side), all boundaries have a no-flow condition.The flow velocity was 1 m/y.At the outflow side, time-constant pressures at the boundaries were set.

Moving-boundary problems and Level Set Method(LSM)
The location of the solid-liquid interface for many phase change problems is not known a priori and must be determined during the course of the analysis.Mathematically, the motion of the interface is expressed implicitly in an equation for the conservation of thermal energy at the interface, the so-called Stefan conditions.These conditions causes the system to become non-linear, and therefore each problem is somewhat unique.Moving boundaries are also associated with time-dependent problems and the position of the boundary has to be determined as a function of time and space.
A level set method has become an efficient tool for tracking, modelling and simulating the motion of free boundaries in fluid mechanics, combustion, computer animation and image processing (Osher & Fedkiw, 2003), (Sethian, 1999).The objective of this part of the chapter is to present the combination of RBF approach and the level set method for solving two-dimensional moving-boundary problems.The LSM is a numerical technique for tracking shapes and interfaces.The advantage of the LSM is that one can perform numerical computations involving curves and surfaces on fixed Cartesian grid.The LSM makes it very easy to follow shapes that change topology (e.g., when a shape splits in two or develops holes).

Level set function and equation
In two dimensions, the LSM represents a closed curve Γ ∈ R 2 in the plane as the zero level set of a two-dimensional auxiliary function Φ(x, t) : R 2 × R → R, where x is a position of the interface, t is a moment in time.Therefore, the closed curve is presented as: The function Φ is also called a level set function and is assumed to take positive values inside the region delimited by the curve Γ and negative values elsewhere (Osher & Fedkiw, 2003), (Sethian, 1999).The level set function can be defined as a signed distance function to the interface.The moving interface is then captured at all time by locating the set of Γ(t) for which Φ vanishes.The movement of the level set function can be described as the following Cauchy problem (Tsai & Osher, 2003): where Φ 0 (x) means the initial position of the interface and v T =[ν 1 , ν 2 ] is the continuous field, which is a function of position x.The above partial differential equation is often solved by using the finite difference method (FDM) on Cartesian grids (Osher & Fedkiw, 2003), (Sethian, 1999).

Level set equation construction with RBFs
The level set equation can be constructed in many ways.In our case, the Crank-Nicolson implicit scheme will be presented.The classical time stepping schemes are stabilized with the adaptive greedy algorithm, see (Vrankar et al., 2010).

146
Advanced Methods for Practical Applications in Fluid Mechanics www.intechopen.com We consider the Crank-Nicolson implicit scheme of (25): where t n+1 = t n + △t, Φ n+1 and Φ n are the level set variable at time t n+1 and t n .

Oriented flow
The first example is translation of circular interface/bubble (the term bubble is connected to the study of multiphase flows (ETH, 2003).In the multi-fluid model, the motion of the specific dispersed phase elements (bubbles, drops) is not followed; rather, the dispersed phase elements which interact with the continuous phase flow field are observed.Therefore, computational modelling is assuming a greater role in the study of multiphase flows.) in oriented flow in Fig. 4. The circular interface of radius r = 0.15, initially centered at (0.2, 0.7) and moved by the oriented flow in a cavity of size 1 × 1 with the velocity field (u, v) defined as follows: For capturing the moving circular interface in Fig. 3, LMS and RBFs are used.The calculation is performed up to t = 2, in computational grid 20 × 20 (upper part of the Fig. 4) and 30 × 30 (lower part of the Fig. 4) .The function Φ is evaluated on a 38 × 38 grid in order to find zero contours.We choose the following shape parameters c = 0.033, c = 0.5, and c = 2 to ensure smoothness of the circular interfaces over time and the initial time step is 0.01.
In such a velocity field, the circular interface is passively transported in the form of circulation and stretching.The LMS and RBFs are used to capture the moving interface with time, up to t = 4 with a time step △t = 0.01, in computational grids 20 × 20 and 30 × 30.The function Φ is then evaluated on a 38 × 38 grid in order to find zero contours.We again choose the large shape parameter c = 2 to demonstrate stability over ill-conditioned linear systems.In Fig. 5 , we can see that with or without the adaptive algorithm, the circular interfaces are more or less at the same location.However, the circular interface obtained with the adaptive algorithm is more skewed; without adaptive algorithm, a more circular bubble is obtained.We can see that some artifacts appear in the corner of the computational domain.

Wood charring problems
The model which explains a very complex phenomena of wood charring consists of the differential equation for heat transfer with corresponding boundary conditions, which prescribe the heat flow on the exposed boundaries of the cross-section.In our case, different types of boundary conditions were used, e.g.Dirichlet, Neumann, radiation term, etc.The char formation in the timber beam as a function of its temperature is taken into account by the model.
Since the analytical solution is seldom obtainable, the problem is solved numerically by the RBF methods (e. g. non-symmetric RBF collocation).Picard's or Newton's methods are used for the solution of the second order non-linear partial differential equations.

Governing equations
The heat and mass transfer is governed by the two second order non-linear partial differential equations (Luikov, 1966).Only one equation which describes heat conduction governed predominantly by temperature gradients was considered: where k x and k y represent thermal conductivity (W/mK) in directions x and y of the cross-section of the beam, ̺ is density (kg/m 3 ), c p specific heat (J/kgK) and T temperature (K).The second equation describes moisture diffusion governed by moisture potential and is not considered here.The problem is complete when initial and boundary conditions are specified.The initial condition prescribes the temperature in the cross-section of the beam at the initial time t = 0 T(x, y,0)=T 0 (x, y). (34) The boundary conditions prescribe the heat flow on the exposed boundaries of a cross-section.Thus, the Neumann boundary conditions at the exposed surface are given by balancing heat where e nx and e ny are components of the normal to the boundary surface and h c is convective heat transfer coefficient (W/m 2 K).T A is the ambient temperature.T R is the temperature of the radiating surface, ε R is the effective surface emissivity and σ is the Stefan-Boltzmann constant σ = 5.671.10 −8 W/m 2 K 4 .

Implicit discrete scheme
For the solution of eqn.( 33) with the corresponding initial and boundary conditions, the RBF method is used.We consider the implicit time scheme of eqns.( 33) and ( 35) : where t n+1 = t n + △t, T n+1 and T n are the temperature at time t n+1 and t n .
The approximate solution is expressed as: where α n+1 j , j = 1, ..., N, are the unknown coefficients to be determined and ϕ j (x)= (x − x j ) 2 +(y − y j ) 2 + c 2 are Hardy's multiquadrics functions (Hardy, 1971).By substituting eqn.(38) into eqns.( 36) and (37) and using factorization for the radiation term (T 4 − T 4 R ), we obtain: where N B and N present the number at boundary and all discretized points.

151
Radial Basis Functions Methods for Solving Radionuclide Migration, Phase Change and Wood Charring Problems www.intechopen.com The system of nonlinear equations which result from the space discretization of a nonlinear PDEs were solved by Picard's methods.

Numerical examples
The results are tested on the one-dimensional case in the standard fire conditions, for which comparison is made with the results of one-dimensional charring rate models for wood presented in the literature (White & Nordheim, 1992).The same model is used to analyse a two-dimensional behaviour of timber beam exposed to fire from three sides.

One-dimensional charring
The charring rate of wood usually refers to the dimensional rate, e.g.millimetres per minute, at which wood transforms to char.Many factors are involved in wood charring.No completely satisfactorily charring model has yet been developed.Since the material properties at elevated temperatures are difficult to obtain, constant material properties of wood and char are used.The following data have been used: Most known models suggest constant charring rates.In our case, we used a (White & Nordheim, 1992) non-linear empirical model for charring rate of eight different wood species.
The comparison to the present model in the case of spruce is shown in Fig. 6.In all empirical models it is assumed that the charring of woods starts instantaneously after exposure to fire.In reality, the charring is not immediate.In our model, the charring starts when the temperature of wood reaches the temperature of pyrolysis, which is around 300 • C.This temperature is reached nearly 3 minutes after the fire starts.

A two-dimensional charring
In the two-dimensional case, the formation of char in a timber beam exposed to the standard fire conditions (ISO 834, 1999) on three sides is considered.The upper surface is thermally isolated.The original beam cross-section is rectangular with dimensions 10 × 15 cm.The beam 152 Advanced Methods for Practical Applications in Fluid Mechanics www.intechopen.comcross-section is discretized by the mesh of 10 × 10 points.Material properties are assumed to be the same as the one-dimensional case.The results of the simulation at 10 and 30 minutes after the exposure to fire are given in Fig. 7.

Conclusions
Modelling of radionuclide migration, moving-boundary and wood charring were presented.The first example presents modelling of radionuclide migration through the geosphere using a combination of radial basis function methods in Eulerian coordinates with geostatistics.In the case of radionuclide migration two steps of evaluations were performed.In the first step the velocities were determined from the pressure of the fluid p by solving the Laplace differential equation.In the second step the advection-dispersion equation was solved to find the concentration of the contaminant.Comparison of the results between the average of contaminant concentrations show that the results are also very similar to the results obtained by finite difference methods (Fig. 3).In the case of calculating the advection-dispersion equation we can conclude that the Kansa method could be an appropriate alternative to the FDM due to its simpler implementation.The second example outlines our work done on an alternative approach to the conventional level set methods for solving two-dimensional moving-boundary problems.This approach is set up from MQ RBFs and the adaptive greedy algorithm.The examples suggest that the solution is more stable by employing the adaptive algorithm.Two examples are presented: oriented flow and shear flow (bubble starts circular but because of the velocity fields the 153 Radial Basis Functions Methods for Solving Radionuclide Migration, Phase Change and Wood Charring Problems bubble stretches).The results of the last example do not depend significantly on the usage of the adaptive greedy algorithm.In the conventional level set methods, the level set equation is solved to evolve the interface using a capturing Eulerian approach.The solving procedure requires appropriate choice of the upwind schemes, reinitialization algorithms and extension velocity methods, which may require excessive amount of computational efforts.In this case we do not choose the reinitialization, because we try to control the smoothness of the moving boundary with appropriate choice of the type of MQ RBFs, shape parameters, time stepping schemes and greedy algorithm.The proposed alternative approach offers to use smaller computational grids, with no reinitialization in order to beat the upwind scheme.In the case of radionuclide migration modelling, the FDM FORTRAN program and the MQ FORTRAN program use comparable CPU times.Usually, fine discretization is applied, but this approch increases CPU time and can cause ill-conditioning.Even when solving two-or three-dimensional problems one tries to obtain the best accuracy with the least amount of CPU time.There are several methods used with RBF methods and other approaches that should be examined before attempting to solve two-or three-dimensional problems with many points, see (Scott & Kansa, 2009).The last example shows the alternative approach for solving wood charring problems with the RBFs methods.Fig. 6 shows that presented approach and the model proposed by (White & Nordheim, 1992) produce comparable results.The same model was used to analyse a two-dimensional behaviour of timber beam exposed to fire from three sides.It shows that the results are comparable to the results obtained in literature (Schnabl & Turk, 2006).

Fig. 1 .
Fig. 1.Graphical presentation of Classical vs. Meshless methods A radial basis function is the function ϕ j (x) := ϕ(||x − x j ||), which depends only on the distance between x ∈ R d and a fixed point x j ∈ R d .Here, ϕ j is continuous and bounded on any bounded sub-domain Ω ⊆ R d whereas ϕ : R d → R. Let r ≥ 0 denote the Euclidean distance between any pair of points in the domain Ω.The commonly used radial basis functions are: linear (ϕ(r)=r), cubic (ϕ(r)=r 3 ), thin-plate spline (ϕ(r)= r 2 log r) and Gaussian (ϕ(r)=e −αr 2 ).The most popular globally supported C ∞ RBFs are MQ(ϕ(r)=( 1 +(r/c) 2 ) β ), β = 1/2),(Wang & Liu, 2002).It is also important to mention the local version of the RBF collocation which does not produce the ill conditioning and is less sensitive to the shape parameter selection.The local version of the RBF collocation is based on compactly supported RBF spline: (ϕ(r)=(1 − r) m + p(r)) where p(r) is a polynomial of the Wendland(Wendland, 1995) compactly supported (CS-RBF) spline.The parameter c > 0 is a shape parameter controlling the fitting of a smoothing surface to the data.It has a significant influence on the accuracy of the solution.For this reason, in almost all previous researches the shape parameter must be adjusted with the number of centers in order to produce equation systems that are sufficiently well conditioned to be solved with the standard finite precision arithmetic.The optimal choice of the constant shape parameter is still an open question, and it is most often selected by the trial and error approach.To overcome the problems of ill-conditioned matrix many efforts have been made to find new computational methods being capable of avoiding the ill-conditioning problems using linear solvers.In the literature, the following methods are reported: with various physical phenomena added.According to this equation, mass transport is controlled by two mechanisms: advection and dispersion.Advection governs the movement of the solute, linked to the fluid, with the water velocity.Water velocity can be assessed through Darcy's law.Dispersion represents mixing of diffusion and random variations from the mean stream.The simulation area will be 2D rectangular with the Neumann and Dirichlet boundary conditions.The Neumann boundary conditions represent flow while Dirichlet boundary conditions represent constant pressure and concentration.

Fig. 6 .
Fig. 6.Comparison of White and Norheim charring model to the present.

Fig. 7 .
Fig. 7. Temperature distribution in the cross-section of spruce beam and the transformation of wood into char at 10 and 30 minutes calculated with the Kansa approach, relaxation parameter: 0.14.