Table of Contents Author Guidelines Submit a Manuscript
Mathematical Problems in Engineering
Volume 2015, Article ID 640305, 13 pages
http://dx.doi.org/10.1155/2015/640305
Research Article

Generalized Finite Difference Time Domain Method and Its Application to Acoustics

1School of Computer Software, Tianjin University, Tianjin 300072, China
2Tianjin Key Laboratory of Cognitive Computing and Application, Tianjin University, Tianjin 300072, China
3Japan Advanced Institute of Science and Technology, Ishikawa 923-1292, Japan

Received 28 August 2014; Revised 11 January 2015; Accepted 11 January 2015

Academic Editor: Shaofan Li

Copyright © 2015 Jianguo Wei et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

A meshless generalized finite difference time domain (GFDTD) method is proposed and applied to transient acoustics to overcome difficulties due to use of grids or mesh. Inspired by the derivation of meshless particle methods, the generalized finite difference method (GFDM) is reformulated utilizing Taylor series expansion. It is in a way different from the conventional derivation of GFDM in which a weighted energy norm was minimized. The similarity and difference between GFDM and particle methods are hence conveniently examined. It is shown that GFDM has better performance than the modified smoothed particle method in approximating the first- and second-order derivatives of 1D and 2D functions. To solve acoustic wave propagation problems, GFDM is used to approximate the spatial derivatives and the leap-frog scheme is used for time integration. By analog with FDTD, the whole algorithm is referred to as GFDTD. Examples in one- and two-dimensional domain with reflection and absorbing boundary conditions are solved and good agreements with the FDTD reference solutions are observed, even with irregular point distribution. The developed GFDTD method has advantages in solving wave propagation in domain with irregular and moving boundaries.

1. Introduction

Partial differential equations (PDEs) modeling problems in science and engineering, such as electromagnetics, acoustics, and hydrodynamics, are usually solved by numerical methods that discretize the computational domain with mesh or grids. Grid-based methods such as finite difference method (FDM), finite element method (FEM), and boundary element method (BEM) [1, 2] have had much achievements and still dominate the field of scientific computing. However, numerical difficulties originating from usage of grids often emerge. For complicated and irregular geometry, implementation of boundary conditions could be a big challenge for FDM. Generation of grids with high quality is not an easy task in FEM and BEM. Moreover, when free surface and moving boundary/interface have to be treated, the transformation of grids will turn the conventional grid-based methods into a difficult, time-consuming process. Numerical accuracy often degenerates and divergence problem occurs.

In recent 20 years, to overcome numerical difficulties due to use of grids or mesh, meshless methods (MMs) based on different techniques have been proposed and widely used in many fields such as hydrodynamics [3], astrophysics [4], and solid mechanics [3, 5]. Among the MMs, generalized finite difference method (GFDM) is the one that evolved from traditional FDM [6, 7] and many different forms have been developed [8]. Benito and his coauthors made great contribution to its recent development [911]. For heat conduction problem, it has been compared with the element-free Garlerkin (EFG) method (one of the most used MMs in solid mechanics) and better performance has been observed [10]. Recently, GFDM was used to solve the wave equations [11] and Burgers’ equations [12] and simulate seismic wave propagation problems in heterogeneous media [13]. An application to the detonation shock dynamics [14] was also carried out. Nevertheless, few work on computational acoustics has been reported.

For acoustic wave propagation problems, the concentration is on the ones in confined domain, for which grid-based methods like FDTD and TDFEM (time-domain finite-element methods) [15], are mostly used. However, moving boundary exists in many acoustic problems like sound wave propagation inside a deforming vocal tract. This problem is hardly solved by conventional grid-based methods and MMs provide a possibility. As one of the MMs, GFDM is extended to transient acoustics in this paper, which is helpful to solve wave propagation problems with moving boundary in the future.

Inspired by the derivation of meshless particle methods, we firstly formulated the GFDM in a way different from the original one that minimizes an energy norm. Such that the relationship between GFDM and meshless particle methods like smoothed particle hydrodynamics (SPH) and its improvements can be conveniently examined. Comparison with the modified dmoothed particle hydrodynamics (MSPH) method, which has better performances than SPH and its corrections [16], shows higher approximation accuracy of the GFDM, especially at the boundary region. By analog with FDTD, a method referred to as generalized finite difference time domain (GFDTD) is proposed, in which GFDM is used to discretize the spatial operators and the leap-frog algorithm is used for time integration. To show its good performance and efficiency, the GFDTD method is applied to transient acoustics. Comparison with conventional FDTD solutions is presented and discussed.

2. Generalized Finite Difference Method (GFDM)

Other than conventional derivation of GFDM by minimizing an energy norm [10], a different derivation of GFDM is presented in this section. Taylor series expansion of around point remaining up to second-order terms yieldswhere , , , and  .

By multiplying both sides of (1) with and ( is a weighting function with compact support) and integrating the resulted equations over the support domain , we get two equations, and the following, as an example, is the result for :where is a volume measure.

Repeating the same procedure with , , and instead of and , we get other three equations and the following is the result for : To approximate the integrations by Riemann sum, the volume of the support domain is divided into points with associated volumes , . Equations (2), (3), and the other three constitute a system of five equations written in matrix form aswithwhere with , , and , and is a measure of the support size.

The conventional derivation of GFDM is presented in the appendix. It is clear that the difference between the conventional and the current derivation is not only the procedure but also the final form. The conventional derivation loses term (see (A.5)). If all the points in the domain have the same volume, at both sides of (4) will be cancelled, and the two final forms will be the same. However, can hardly be the same when points are irregularly spaced. From this point of view, our derived final form is more general and takes point irregularity into account.

3. Modified Smoothed Particle Hydrodynamics (MSPH)

As a modification to SPH, the MSPH method improves the accuracy of the approximations especially at points near the boundary of the domain [16]. It uses Taylor series expansion of function as in (1). Similar to the derivations of (2) and (3), but with different weight functions , , , , and , the following equations, as examples, for and , are obtained:Again the Riemann sum over the support domain is used to approximate the integrations and a system of five equations is obtained asCompared with formula (4), the only difference is the terms multiplied to both sides of (1). In GFDM, , , , , and are used instead of , , , , and in MSPH. As a result, GFDM avoids computing the derivatives of the weight function and hence saves computational efforts and leads to more choice of the weight function.

4. Numerical Tests for Approximation of Derivatives

In previous sections the deviation of GFDM and MSPH is presented. In this section, to compare the performance of the two methods, they are used to approximate the derivatives of certain 1D and 2D functions. For the convenience of evaluation, a global error measure is defined as follows: where can be , , , and and the superscripts and refer to the exact and numerical solutions, respectively.

The quartic spline function is used as the weight function : where is the kernel radius taken as ( is the space interval) which is usually used in meshless methods.

4.1. One-Dimensional Case

Consider the following function:Figure 1 shows the first- and second-order derivatives estimated by GFDM and MSPH and the exact results when the domain is discretized into 21 equally spaced points. It is seen that GFDM has better performance in both derivatives especially for the points near boundaries. When the number of points increases to 51, the results are similar as exhibited in Figure 2. Error analysis shown in Table 1 indicates that GFDM has higher accuracy. With increasing number of points, the global error decreases.

Table 1: Approximation errors in the derivatives of the function .
Figure 1: Estimates of the first derivative (a) and second derivative (b) of function with 21 equally spaced points on .
Figure 2: Estimates of the first derivative (a) and second derivative (b) of function with 51 equally spaced points on .
4.2. Two-Dimensional Case

For the functionits first- and second-order derivatives together with estimations by GFDM and MSPH are shown in Figure 3. In each direction 21 points are employed. As expected, GFDM has higher approximation accuracy than MSPH for both first- and second-order derivatives as shown in Table 2.

Table 2: Approximation errors in the derivatives of the function .
Figure 3: Exact (top) and estimated (a) and (b) of function by GFDM (middle) and MSPH (bottom).

5. Generalized Finite Difference Time Domain Method for Computational Acoustics

For computational acoustics, the mostly used approach is the FDTD method, which was originally designed for the simulation of electromagnetics [1, 2]. As a finite difference scheme, its applicability to complex problems suffers from aforementioned difficulties, for which the generalized finite difference can be a good alternative. In this section, together with the basics of computational acoustics, a meshless method is proposed, in which GFDM is used to discretize the spatial derivatives and the leap-frog algorithm is used to discretize the temporal derivatives. By analog with FDTD, it is referred to as generalized finite difference time domain (GFDTD) method and is expected to have advantages due to its meshless property.

The governing equations for acoustic wave propagation problems arewhere is pressure, is particle velocity, is the density of the medium, and is the speed of sound.

5.1. Spatial Derivative Approximations by GFDM

The spatial derivatives on the right-hand side of (12) are approximated by GFDM. By solving (4) we get the approximations of and . That is, the derivatives of variable at point can be approximated by function values at points inside the support domain centered at aswhere and are the difference coefficients for center point and and are coefficients for other points in the support domain. As GFDM can reproduce constant functions [4], we haveBy taking both and in (12) as , the approximated spatial operators in (12) are accordingly obtained.

5.2. Explicit Leap-Frog Scheme in GFDTD

Generally, the temporal derivatives on the left-hand side of (12) can be integrated by any time marching algorithms. Inspired by the conventional FDTD method, the second-order accurate explicit leap-frog scheme is used herein, in which two variables and are alternatively calculated. The velocity is computed at the half time step and the pressure is calculated at the integer time step [2]. After temporal approximations the semi-discretization of (12) becomeswhere superscript represents the time step. The leap-frog scheme is conditionally stable and the time step should satisfy the Courant-Friedrichs-Lewy (CFL) condition; that is, , where is the dimension of the problem.

Substituting (13) into the right-hand side of (15) as the approximation of the first-order spatial derivatives, the full discretized system of equations becomeswhere the particle velocity is a 2D vector .

By analog with FDTD, the full discretization scheme is referred to as generalized finite difference time domain (GFDTD) method.

6. Numerical Results

To validate the proposed GFDTD method, it is applied to one- and two-dimensional wave propagation problems. Three cases are presented. The first two examine the acoustic wave propagation in one- and two-dimensional domain, respectively, and the third one is a real case with different types of boundary conditions. All the cases use (9) as the weight function. The other parameters are and  m/s and the time interval is set to be μs to satisfy the CFL condition. FDTD solutions are chosen as the reference.

6.1. One-Dimensional Case

In this case 501 points are equally spaced in the domain . In the middle of it, there is a wave source in Gaussian pulse form:

As shown in Figure 4, the simulated results at two time levels μs and 500 μs have good agreement with the FDTD solutions. Compared with FDTD, the relative errors concerning the pressure and are

Figure 4: Gaussian wave propagation at two time points: (a) μs and (b) 500 μs.
6.2. Two-Dimensional Case

The Gaussian wave propagation in two-dimensional domain is simulated in this section. The length of the square domain is 0.1 m and 101 points are uniformly distributed in each direction. The wave source starts from the middle of the domain. Figure 5 compares the solutions of FDTD and GFDTD at μs. The results along  m are shown in Figure 5(c). Again, good agreement is observed and the relative error is less than 2%.

Figure 5: Gaussian wave propagation in a square domain at time μs. (a) The FDTD result, (b) the GFDTD result, and (c) the comparison along  m.

To show the advantage of the proposed GFDTD over the conventional FDTD, irregular point distribution is examined. All the points used above are allowed to have perturbation around their original locations to make the distribution irregular, part of which is shown in Figure 6(a) and the result is shown in Figure 6(b). The comparison of the result along  m shown in Figure 6(c) indicates that, with irregular distribution of computational points, the Gaussian wave propagates as well as before.

Figure 6: Gaussian wave propagation in a square domain with irregular point distribution. (a) Part of the irregular distribution, (b) simulated results, and (c) comparison with regular point distribution after μs.

In the GFDTD simulation of wave propagation with irregular point distribution, the volume associated to each point had better to be considered as analyzed at the end of Section 2. In 2D case, the volume associated with a given point is the area that the point dominates. Here we use Delaunay triangulation and Voronoi diagram [17] to calculate the area and the results are shown in Figure 7. The volume of each point is shown in Figure 7(a). Due to the designed perturbation, the volume fluctuates around  m2. Based on Delaunay triangulation and Voronoi diagram, the area associated with a red point is indicated by the area surrounded by the blue lines as shown in Figure 7(b). To compare the results with regular and irregular point distributions, the values on regularly distributed points have to be firstly interpolated by the calculated results with irregularly spaced points. The third-order accurate cubic interpolation method [17] is used herein. The relative error is The relatively small error is due to the gentle irregular point distribution. That is, if the point irregularity is small, the effect of is negligible. This is consistent with Benito’s results [911].

Figure 7: (a) Area associated with irregular distributed points and (b) Voronoi diagram.
6.3. Two-Dimensional Case with Effect of Boundary Conditions

In Sections 6.1 and 6.2, wave propagation inside a domain is simulated. In this section, to show the effect of boundaries, which is of high importance in transient acoustics, sound wave propagation in a rectangular tube is studied. Two different kinds of boundary conditions including reflection and absorbing boundary are considered. At the left edge of the computational domain there is a Gaussian pulse as the source term. The upper and bottom boundaries are reflection layers and the second-order Mur’s absorbing boundary condition [18] exists at the right side boundary. In this case, the same and are used as before and the time interval is μs. Inside the computational domain points with spatial interval  mm are evenly spaced.

6.3.1. Source Term

The left is a wave source with pressure given by the Gaussian pulse:where and  KHz.

6.3.2. Reflection Boundary Condition

To simulate reflections at the upper and bottom wall boundaries, the model proposed by Yokota et al. [19, 20] and widely used in room acoustics is employed herein. In this model, the normal component of particle velocity and the pressure of the points on the boundary are supposed to satisfy the following condition: where is the normal acoustic impedance on the boundary given byHere the normal sound absorption coefficient is taken as 0.2 as in [19].

6.3.3. Absorbing Boundary Condition

At the right boundary, second-order Mur’s absorbing boundary condition [18] is applied:By applying (12) to (23) and performing time integration, (23) degenerates toWhen GFDTD is used, the discrete form of (24) is obtained as

6.3.4. Results

After applying the source term and the two boundary conditions into our case, the wave is considered to be propagating from left to right inside a tube and gets absorbed at the end of it. Figure 8 shows the simulated results after 200 μs with a color map image that clearly depicts the pressure distribution. In Figure 9, the results after 350 μs are depicted and the absorbing boundary at the right edge leads to no reflection.

Figure 8: Sound wave propagation in a two-dimension tube. (a) FDTD results, (b) GFDTD results, and (c) the comparison along m at μs.
Figure 9: Sound wave propagation in a two-dimension tube. (a) FDTD results, (b) GFDTD results, and (c) the comparison along  m at μs.

7. Conclusion

A new derivation of the generalized finite difference method (GFDM) with Taylor series expansion generates the same formulation as its conventional derivation and clearly demonstrates its relationship with meshless particle methods. GFDM has better performance in derivative approximations than the particle methods. The proposed generalized finite difference time domain (GFDTD) method has been successfully applied to one- and two-dimensional acoustic wave propagation problems with reflection and absorbing boundary conditions. The numerical results are in line with the FDTD reference solutions even with irregular point distribution. The GFDTD method has high potentials in solving transient acoustic problems with moving boundaries, which deserves further studies.

Appendix

Conventional Derivation of GFDM

Considering the 2D case, for the same Taylor expansion in (1), we consider an energy norm :where , , , , and is weighing function with compact support.

The solution of the derivatives is obtained by minimizing the norm , that is, withFor example, the first equation is

Equation (A.4) and the other four give the following system:

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgment

This work is supported in part by the National Natural Science Foundation of China (nos. 51478305 and 61175016) and Key Program (no. 61233009). The authors thank the anonymous reviewers for their most useful suggestions.

References

  1. M. N. O. Sadiku, Numerical Techniques in Electromagnetics, CRC Press, 2000.
  2. D. M. Sullivan, Electromagnetic Simulation Using the FDTD Method, John Wiley & Sons, 2013.
  3. J. W. Swegle and S. W. Attaway, “On the feasibility of using Smoothed Particle Hydrodynamics for underwater explosion calculations,” Computational Mechanics, vol. 17, no. 3, pp. 151–168, 1995. View at Publisher · View at Google Scholar · View at Scopus
  4. T. Belytschko, Y. Krongauz, D. Organ, M. Fleming, and P. Krysl, “Meshless methods: an overview and recent developments,” Computer Methods in Applied Mechanics and Engineering, vol. 139, no. 1–4, pp. 3–47, 1996. View at Publisher · View at Google Scholar · View at Scopus
  5. G. R. Johnson, R. A. Stryk, and S. R. Beissel, “SPH for high velocity impact computations,” Computer Methods in Applied Mechanics and Engineering, vol. 139, no. 1-4, pp. 347–373, 1996. View at Publisher · View at Google Scholar · View at Scopus
  6. N. Perrone and R. Kao, “A general finite difference method for arbitrary meshes,” Computers & Structures, vol. 5, no. 1, pp. 45–57, 1975. View at Publisher · View at Google Scholar · View at Scopus
  7. P. S. Jensen, “Finite difference techniques for variable grids,” Computers & Structures, vol. 2, no. 1-2, pp. 17–29, 1972. View at Publisher · View at Google Scholar · View at Scopus
  8. T. Liszka and J. Orkisz, “The finite difference method at arbitrary irregular grids and its application in applied mechanics,” Computers and Structures, vol. 11, no. 1-2, pp. 83–95, 1980. View at Publisher · View at Google Scholar · View at Scopus
  9. J. J. Benito, F. Urea, and L. Gavete, “Influence of several factors in the generalized finite difference method,” Applied Mathematical Modelling, vol. 25, no. 12, pp. 1039–1053, 2001. View at Publisher · View at Google Scholar · View at Scopus
  10. L. Gavete, M. L. Gavete, and J. J. Benito, “Improvements of generalized finite difference method and comparison with other meshless method,” Applied Mathematical Modelling, vol. 27, no. 10, pp. 831–847, 2003. View at Publisher · View at Google Scholar · View at Scopus
  11. J. J. Benito, F. Urena, and L. Gavete, “Solving parabolic and hyperbolic equations by the generalized finite difference method,” Journal of Computational and Applied Mathematics, vol. 209, no. 2, pp. 208–233, 2007. View at Publisher · View at Google Scholar · View at MathSciNet · View at Scopus
  12. C.-M. Fan and P.-W. Li, “Generalized finite difference method for solving two-dimensional Burgers’ equations,” Procedia Engineering, vol. 79, pp. 55–60, 2014. View at Publisher · View at Google Scholar
  13. J. J. Benito, F. Ureña, L. Gavete, E. Salete, and A. Muelas, “A GFDM with PML for seismic wave equations in heterogeneous media,” Journal of Computational and Applied Mathematics, vol. 252, pp. 40–51, 2013. View at Publisher · View at Google Scholar · View at MathSciNet · View at Scopus
  14. Y. L. Chen, K. B. Huang, and X. Yu, “Numerical study of detonation shock dynamics using generalized finite difference method,” Science China: Physics, Mechanics and Astronomy, vol. 54, no. 10, pp. 1883–1888, 2011. View at Publisher · View at Google Scholar · View at Scopus
  15. J. F. Lee, R. Lee, and A. Cangellaris, “Time-domain finite-element methods,” IEEE Transactions on Antennas and Propagation, vol. 45, no. 3, pp. 430–442, 1997. View at Publisher · View at Google Scholar · View at MathSciNet · View at Scopus
  16. G. M. Zhang and R. C. Batra, “Modified smoothed particle hydrodynamics method and its application to transient problems,” Computational Mechanics, vol. 34, no. 2, pp. 137–146, 2004. View at Publisher · View at Google Scholar · View at Scopus
  17. M U Guide, The mathworks Inc., Natick, Mass, USA, 2013.
  18. G. Mur, “Absorbing boundary conditions for the finite-difference approximation of the time-domain electromagnetic-field equations,” IEEE Transactions on Electromagnetic Compatibility, vol. 23, no. 4, pp. 377–382, 1981. View at Google Scholar · View at Scopus
  19. T. Yokota, S. Sakamoto, and H. Tachibana, “Visualization of sound propagation and scattering in rooms,” Acoustical Science and Technology, vol. 23, no. 1, pp. 40–46, 2002. View at Publisher · View at Google Scholar · View at Scopus
  20. S. Sakamoto, H. Nagatomo, A. Ushiyama, and H. Tachibana, “Calculation of impulse responses and acoustic parameters in a hall by the finite-difference time-domain method,” Acoustical Science and Technology, vol. 29, no. 4, pp. 256–265, 2008. View at Publisher · View at Google Scholar · View at Scopus