Discrete Dynamics in Nature and Society

Volume 2014 (2014), Article ID 645726, 12 pages

http://dx.doi.org/10.1155/2014/645726

Research Article

Modified Differential Transform Method for Solving the Model of Pollution for a System of Lakes

Brahim Benhammouda,¹ Hector Vazquez-Leal,² and Luis Hernandez-Martinez³

¹Higher Colleges of Technology, Abu Dhabi Men’s College, P.O. Box 25035, Abu Dhabi, UAE

²Electronic Instrumentation and Atmospheric Sciences School, Universidad Veracruzana, Circuito Gonzalo Aguirre Beltrán S/N, 91000 Xalapa, VER, Mexico

³National Institute for Astrophysics, Optics, and Electronics, Luis Enrique Erro No. 1, Santa Maria 72840 Tonantzintla, PUE, Mexico

Received 19 May 2014; Revised 19 August 2014; Accepted 1 September 2014; Published 15 September 2014

Academic Editor: Carlo Piccardi

Copyright © 2014 Brahim Benhammouda 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.

This work presents the application of the differential transform method (DTM) to the model of pollution for a system of three lakes interconnected by channels. Three input models (periodic, exponentially decaying, and linear) are solved to show that DTM can provide analytical solutions of pollution model in convergent series form. In addition, we present the posttreatment of the power series solutions with the Laplace-Padé resummation method as a useful strategy to extend the domain of convergence of the approximate solutions. The Fehlberg fourth-fifth order Runge-Kutta method with degree four interpolant (RKF45) numerical solution of the lakes system problem is used as a reference to compare with the analytical approximations showing the high accuracy of the results. The main advantage of the proposed technique is that it is based on a few straightforward steps and does not generate secular terms or depend of a perturbation parameter.

1. Introduction

Semianalytical methods like differential transform method (DTM) [1–4], reduced differential transform method (RDTM) [5–7], homotopy perturbation method (HPM) [8–16], homotopy analysis method (HAM) [17], variational iteration method (VIM) [18], and generalized homotopy method [19], multivariate Padé series [20], among others, are powerful tools to approximate linear and nonlinear problems in physics and engineering. Analytical solutions aid researchers to study the effect of different variables or parameters on the function under study easily [21]. Among the above-mentioned methods, the DTM is highlighted by its simplicity and versatility to solve nonlinear differential equations. This method does not rely on a perturbation parameter or a trial function as other popular approximative methods. In [1], the DTM was introduced to the engineering field as a tool to find approximate solutions of electrical circuits. DTM produces approximations based on an iterative procedure derived from the Taylor series expansion. This method is very effective and powerful for solving various kinds of differential equations as nonlinear biochemical reaction model [2], two point boundary-value problems [22], differential-algebraic equations [23], the KdV and mKdV equations [24], the Schrodinger equations [25], fractional differential equations [26], and the Riccati differential equation [27], among others.

Therefore, in this paper, we present the application of a hybrid technique combining DTM, Laplace transform, and Padé approximant [28] to find approximate analytical solutions for a pollution model [29–34]. The aim of the model is to describe the pollution of a system of three lakes [35–39] as depicted in Figure 1. Each lake is considered to be as large compartment and the interconnecting channels as pipes between the compartments with given flow directions. Initially, a pollutant is introduced into the first lake at a given rate which may be constant or may vary with time. Therefore, we are interested in knowing the level of pollution in each lake at any time. We assume the pollutant in each lake to be uniformly distributed throughout the lake by some mixing process, and the volume of water in any of the three lakes remains constant. Also we assume that the type of pollution is persistent and not degrading to other forms.

Figure 1: Three lakes system with interconnecting channels.

We will consider three input models for the pollution source: periodic, exponentially decaying, and linear. To the best of the authors’ knowledge, the periodic and exponential decaying input models are not reported in literature as potential behaviours for the pollution source. Solutions to this problem are first obtained in convergent series form using the DTM. To improve the solution obtained from DTM’s truncated series, we apply Laplace transform to it, and then convert the transformed series into a meromorphic function by forming its Padé approximant. Finally, we take the inverse Laplace transform of the Padé approximant to obtain the approximate analytical solution. This hybrid method (LPDTM) which combines DTM with Laplace-Padé posttreatment greatly improves DTM’s truncated series solutions in convergence rate. In fact, the Laplace-Padé resummation method enlarges the domain of convergence of the truncated power series and often leads to an accurate approximation or the exact solution. It is worth mentioning that the constant input model is not treated here, because it was successfully solved by the DTM in [40].

The proposed method does not generate noise terms also known as secular terms in the solution as the homotopy perturbation based techniques [14]. This property of DTM greatly reduces the volume of computation and improves the efficiency of LPDTM in comparison to the perturbation based methods. What is more, LPDTM does not require a perturbation parameter as the perturbation based techniques (including HPM). Finally, LPDTM is straightforward and can be programmed using computer algebra packages like Maple or Mathematica.

The rest of this paper is organized as follows. In the next section we describe how the DTM can be applied to solve systems of ordinary differential equations. The main idea behind the Padé approximant is given in Section 3. In Section 4, we give the basic concept of the Laplace-Padé resummation method. In Section 5, we give a description of the lakes pollution problem. In Section 6, we apply LPDTM to solve three pollution problems. In Section 7, we give a brief discussion. Finally, a conclusion is drawn in the last section.

2. Differential Transform Method (DTM)

The basic definitions and fundamental operations of differential transform are given in [1, 22–26]. For convenience of the reader, we will give a review of the DTM. We will also describe the DTM to solve systems of ordinary differential equations.

Definition 1. If a function is analytical with respect to in the domain of interest , then

is the transformed function of .

Definition 2. The differential inverse transforms of the set is defined by

Substituting (1) into (2), one deduces that

From Definitions 1 and 2, it is easy to see that the concept of the DTM is obtained from the power series expansion. To illustrate the application of the proposed DTM to solve systems of ordinary differential equations, one considers the nonlinear system

where is a nonlinear smooth function.

System (4) is supplied with some initial conditions

DTM establishes that the solution of (4) can be written as

where are unknowns to be determined by DTM.

Applying the DTM to the initial conditions (5) and system (4), respectively, one obtains the transformed initial conditions

and the recursion system

where is the differential transforms of .

Using (7) and (8), we determine the unknowns ,   Then, the differential inverse transformation of the set of values gives the approximate solution

where is the approximation order of the solution. The exact solution of problem (4)-(5) is then given by

If and are the differential transforms of and , respectively, then the main operations of DTM are shown in Table 1.

Table 1: Main operations of DTM.

Function
Differential transform

The process of DTM can be described as follows.

(1)Apply the differential transform to the initial conditions (5).(2)Apply the differential transform to the differential system (4) to obtain a recursion system for the unknowns (3)Use the transformed initial conditions (7) and the recursion system (8) to determine the unknowns (4)Use the differential inverse transform formula (9) to obtain an approximate solution for the initial value problem (4)-(5).

The solutions series obtained from DTM may have limited regions of convergence, even if we take a large number of terms. Therefore, we propose to apply the Laplace-Padé resummation method to DTM truncated series to enlarge the convergence region as depicted in the next sections.

3. Padé Approximant

Given an analytical function with Maclaurin’s expansion

The Padé approximant to of order which we denote by is defined by [28]

where we considered , and the numerator and denominator have no common factors.

The numerator and the denominator in (12) are constructed so that and and their derivatives agree at up to . That is,

From (13), we have

From (14), we get the following algebraic linear systems:

From (15), we calculate first all the coefficients , . Then, we determine the coefficients , from (16).

Note that for a fixed value of , error (13) is smallest when the numerator and denominator of (12) have the same degree or when the numerator has degree one higher than the denominator.

4. Laplace-Padé Resummation Method

Several approximate methods provide power series solutions (polynomial). Nevertheless, sometimes, this type of solutions lacks of large domains of convergence. Therefore, Laplace-Padé [29–34] resummation method is used in literature to enlarge the domain of convergence of solutions or inclusive to find the exact solutions.

The Laplace-Padé method can be summarized as follows.

(1)First, Laplace transformation is applied to power series (9).(2)Next, is substituted by in the resulting equation.(3)After that, we convert the transformed series into a meromorphic function by forming its Padé approximant of order . and are arbitrarily chosen, but they should be smaller than the order of the power series. In this step, the Padé approximant extends the domain of the truncated series solution to obtain a better accuracy and convergence.(4)Then, is substituted by .(5)Finally, by using the inverse Laplace transformation, we obtain the exact or approximate solution.

5. Description of Pollution Problem

A system of lakes is a set of lakes interconnected by channels. These lakes are modelled by large compartments interconnected by pipes [36]. Figure 1 shows a system of three lakes. At , a pollutant is introduced, for example, from a factory into one of the lakes (here lake 1) at rate . Then the polluted water flows into the other lakes through the channels or pipes as indicated by the arrows [41]. We also assume that the volume of water in each lake does not change and that the pollutant is persistent and uniformly distributed in each lake. With these assumptions, we want to predict the level of pollution in each lake for .

To model the dynamic behavior of the system of lakes, let and , denote the volume of water and the amount of pollutant in lake , respectively. Then the concentration of the pollutant in lake at time is given by

If we assume further that the flow rate from lake to lake is constant, then the flux of the pollutant flowing from lake into lake for is given by

Thus, measures the rate at which the concentration of the pollutant in lake flows into lake at time .

Applying the principle

to each lake, we obtain the following system of first order ordinary differential equations:

If we assume that the lakes are initially free from pollutant, then the initial conditions for (20) are

Since the volume of water in each lake is constant for , then the rate of incoming flow is equal to the rate of outgoing flow for each lake. This leads to the following conditions on flow rates:

For results comparison, we consider throughout this paper the following values of the parameters (20) [35]:

6. Numerical Simulation

In this section we will apply the LPDTM described in the previous sections to find approximate analytical solutions for three pollution models to illustrate the accuracy and effectiveness of the method. To simulate the pollution in the lakes we coded the LPDTM in Maple17.

6.1. Periodic Input Model

This input model is used when the pollutant is introduced into the lake 1 periodically. As an example we take , where is the average input of concentration of pollutant, is the amplitude of fluctuations, and is the frequency of fluctuations. Taking and the parameters values given in (23), then system (20) becomes

with the initial conditions

For the solution procedure with LPDTM, we take the differential transform of (25) and system (24), respectively, to get

and the recursion system

Using (26) and recursion system (27), we compute the first few terms

Using (9) and (28), we obtain the seventh and eight order solution approximations

The solutions series obtained from the DTM may have limited regions of convergence, even if we take more terms. We can increase accuracy by applying the Laplace-Padé posttreatment described in the previous sections. First, we apply -Laplace transforms to (29). Then, we substitute by and apply -Padé approximants to the transformed series. Finally, we substitute by and apply the inverse Laplace -transforms to the resulting expressions to obtain the approximate solution.

Applying Laplace transform to (29) yields

For the sake of simplicity we let in (30) to obtain

From (31) we compute the -Padé approximants , , and of , and , respectively, to get

Now since , we obtain , , and in terms of as follows:

Finally, applying the inverse -Laplace transforms to the Padé approximants (33), we obtain the following approximate solutions for pollution problem (24)-(25):

6.2. Exponentially Decaying Input Model

In this example, we assume the pollutant input to have the form . This corresponds to a heavy dumping of the pollutant. If we take , , and the parameters values given in (23), then system (20) becomes

with the initial conditions

For the solution procedure with LPDTM, we take the differential transform of (36) and system (35), respectively, to get

and the recursion system

Computing , for , and using (9) we obtain the seventh and fourth order solution approximations

Applying -Laplace transforms to (39). Then, substituting by and taking the -Padé approximants , , and of and , respectively. Finally, substituting by and applying the inverse Laplace -transforms to the resulting expressions, we obtain the following approximate solution for the pollution (35)-(36):

6.3. Linear Input Model

This input model is used when the pollutant is introduced into the first lake with a linear concentration; that is , where is a positive constant. If we take and the parameters values given in (23), then system (20) becomes

with the initial conditions

For the solution procedure with LPDTM, we take the differential transform of (42) and system (41), respectively, to get

and the recursion system

Computing , for , and using (9) we obtain the sixth order solution approximation

Applying -Laplace transforms to (45). Then, substituting by and applying Padé approximants to the transformed series, and finally, substituting by and applying the inverse Laplace -transforms to the resulting expressions, we obtain the approximate solution for (41)-(42)

7. Discussion

In this work, we presented the differential transform method (DTM) as a useful analytical tool to solve a pollution model for a system of three lakes. Three input models were successfully solved. For each of the three cases solved here, the DTM transformed the dynamic model into a recursion system for the coefficients of the power series solution. To improve the convergence of the DTM solution, a Laplace-Padé posttreatment is applied to the DTM’s truncated series leading to the approximate solution. Additionally, the solution procedure does not involve unnecessary computation like that related to noise terms [14]. This property of the DTM greatly reduces the volume of computation and improves the efficiency of the proposed method. It should be noticed that these problems were effectively handled by LPDTM method due to the malleability of DTM and resummation capability of Laplace-Padé.

For comparison purposes, the Fehlberg fourth-fifth order Runge-Kutta method with degree four interpolant (RKF45) [42, 43] build-in in Maple CAS software was used to obtain the exact solution of the pollution problems. The routine was configured to use an absolute error (A.E.) of . Figures 2, 3, and 4 show the comparison among exact solutions and the LPDTM approximations for the input models: periodic, exponential, and linear. In general terms, the absolute error is low for all approximations. If required, the error can be reduced by increasing the order of the DTM approximations in combination with a posttreatment of Laplace-Padé resummation method of higher order.

(a) Pollution for lake 1

(b) Pollution for lake 2

(c) Pollution for lake 3

(d) A.E. for pollution of lake 1

(e) A.E. for pollution of lake 2

(f) A.E. for pollution of lake 3

Figure 2: Periodic input signal and absolute error (A.E.) of approximations. The time is expressed in years.

(a) Pollution for lake 1

(b) Pollution for lake 2

(c) Pollution for lake 3

(d) A.E. for pollution of lake 1

(e) A.E. for pollution of lake 2

(f) A.E. for pollution of lake 3

Figure 3: Exponential input signal and absolute error (A.E.) of approximations. The time is expressed in years.

(a) Pollution for lake 1

(b) Pollution for lake 2

(c) Pollution for lake 3

(d) A.E. for pollution of lake 1

(e) A.E. for pollution of lake 2

(f) A.E. for pollution of lake 3

Figure 4: Linear input signal and absolute error (A.E.) of approximations. The time is expressed in years.

On the one hand, semianalytical methods like HPM, HAM, and VIM, among others, require an initial approximation for the solutions sought and the computation of one or several adjustment parameters. If the initial approximation is properly chosen, the results can be highly accurate. Nonetheless, no general methods are available to choose such initial approximation. This issue motivates the use of adjustment parameters obtained by minimizing the least-squares error with respect to the numerical solution.

On the other hand, DTM or LPDTM methods do not require any trial equation as a requisite for starting the method. Moreover, DTM obtains its coefficients using an easily computable straightforward procedure that can be implemented into programs like Maple or Mathematica.

8. Conclusion

This work presents LPDTM as a combination of DTM and a resummation method based on the Laplace transforms and Padé approximant. Firstly, the solutions of a pollution model of a system of three lakes are obtained in convergent series forms using DTM. Next, in order to enlarge the domain of convergence of the truncated power series, a posttreatment combining Laplace transform and Padé approximant is applied. This technique greatly improves the DTM’s truncated series solutions in convergence rate. Additionally, DTM is an attractive tool, because it does not require a perturbation parameter to work and it does not generate secular terms (noise terms) as other semianalytical methods like HPM, HAM, or VIM. The proposed method (LPDTM) is based on a straightforward procedure, suitable for engineers. Finally, further research should be performed to solve other highly nonlinear dynamic models.

Conflict of Interests

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

Acknowledgment

Hector Vazquez-Leal gratefully acknowledges the financial support provided by the National Council for Science and Technology of Mexico (CONACyT) through Grant CB-2010-01 no. 157024.

References

J. K. Zhou, Differential Transformation and Its Applications for Electrical Circuits, Huazhong University Press, Wuhan, China, 1986 (Chinese).

A.-M. Batiha and B. Batiha, “Differential transformation method for a reliable treatment of the nonlinear biochemical reaction model,” Advanced Studies in Biology, vol. 3, pp. 355–360, 2011.

Y. Khan, Z. Svoboda, and Z. Šmarda, “Solving certain classes of Lane-Emden type equations using the differential transformation method,” Advances in Difference Equations, vol. 2012, article 174, 2012.

Z. Šmarda, J. Diblík, and Y. Khan, “Extension of the differential transformation method to nonlinear differential and integro-differential equations with proportional delays,” Advances in Difference Equations, vol. 2013, article 69, 12 pages, 2013.

Y. Keskin and G. Oturanç, “Reduced differential transform method for partial differential equations,” International Journal of Nonlinear Sciences and Numerical Simulation, vol. 10, no. 6, pp. 741–749, 2009.

Y. Keskin and G. Oturanc, “The reduced differential transform method: a new approach to factional partial differential equations,” Nonlinear Science Letters, vol. 1, no. 2, pp. 207–217, 2010.

Y. Keskin and G. Oturanc, “Reduced differential transform method for generalized KdV equations,” Mathematical & Computational Applications, vol. 15, no. 3, pp. 382–393, 2010.

M. M. Rashidi, M. T. Rastegari, M. Asadi, and O. A. Bég, “A study of non-newtonian flow and heat transfer over a non-isothermal wedge using the homotopy analysis method,” Chemical Engineering Communications, vol. 199, no. 2, pp. 231–256, 2012.

M. M. Rashidi, S. A. M. Pour, T. Hayat, and S. Obaidat, “Analytic approximate solutions for steady flow over a rotating disk in porous medium with heat transfer by homotopy analysis method,” Computers and Fluids, vol. 54, no. 1, pp. 1–9, 2012.

J. Biazar and B. Ghanbari, “The homotopy perturbation method for solving neutral functional-differential equations with proportional delays,” Journal of King Saud University—Science, vol. 24, no. 1, pp. 33–37, 2012.

J. Biazar and M. Eslami, “A new homotopy perturbation method for solving systems of partial differential equations,” Computers and Mathematics with Applications, vol. 62, no. 1, pp. 225–234, 2011.

J. He, “Homotopy perturbation technique,” Computer Methods in Applied Mechanics and Engineering, vol. 178, no. 3-4, pp. 257–262, 1999.

H. Vazquez-Leal, Y. Khan, G. Fernández-Anaya, A. Herrera-May, and U. Filobello-Nino, “A general solution for Troesch's problem,” Mathematical Problems in Engineering, vol. 2012, Article ID 208375, 14 pages, 2012.

F. Soltanian, M. Dehghan, and S. M. Karbassi, “Solution of the differential algebraic equations via homotopy perturbation method and their engineering applications,” International Journal of Computer Mathematics, vol. 87, no. 9, pp. 1950–1974, 2010.

M. A. Asadi, F. Salehi, and M. M. Hosseini, “Modification of the homotopy perturbation method for nonlinear system of second-order BVPs,” Journal of Computer Science & Computational Mathematics, vol. 2, no. 5, 2012.

F. Salehi, M. A. Asadi, and M. M. Hosseini, “Solving system of DAEs by modified homotopy perturbation method,” Journal of Computer Science & Computational Mathematics, vol. 2, no. 6, pp. 1–5, 2012.

F. Guerrero, F. J. Santonja, and R. J. Villanueva, “Solving a model for the evolution of smoking habit in Spain with homotopy analysis method,” Nonlinear Analysis: Real World Applications, vol. 14, no. 1, pp. 549–558, 2013.

Y. Khan, H. Vazquez-Leal, L. Hernandez-Martinez, and N. Faraz, “Variational iteration algorithm-II for solving linear and non-linear ODEs,” International Journal of the Physical Sciences, vol. 7, no. 25, pp. 3099–4002, 2012.

H. Vazquez-Leal, “Generalized homotopy method for solving nonlinear differential equations,” Computational & Applied Mathematics, vol. 33, no. 1, pp. 275–288, 2014.

M. Yigider and E. Celik, “The numerical solution of partial differential-algebraic equations,” Advances in Difference Equations, vol. 2013, article 10, 2013.

Y. Keskin and G. Oturanc, “Differential transform method for solving linear and nonlinear wave equations,” Iranian Journal of Science & Technology, Transaction A, vol. 34, no. 2, pp. 113–122, 2010.

C. L. Chen and Y. C. Liu, “Solution of two-point boundary-value problems using the differential transformation method,” Journal of Optimization Theory and Applications, vol. 99, no. 1, pp. 23–35, 1998.

F. Ayaz, “Applications of differential transform method to differential-algebraic equations,” Applied Mathematics and Computation, vol. 152, no. 3, pp. 649–657, 2004.

F. Kangalgil and F. Ayaz, “Solitary wave solutions for the KDV and mKDV Equations by differential transform method,” Chaos, Solitons and Fractals, vol. 41, no. 1, pp. 464–472, 2009.

A. S. V. Ravi Kanth and K. Aruna, “Two-dimensional differential transform method for solving linear and non-linear Schrödinger equations,” Chaos, Solitons and Fractals, vol. 41, no. 5, pp. 2277–2281, 2009.

A. Arikoglu and I. Ozkol, “Solution of fractional differential equations by using differential transform method,” Chaos, Solitons and Fractals, vol. 34, no. 5, pp. 1473–1481, 2007.

J. Biazar and M. Eslami, “Differential transform method for quadratic riccati differential equation,” International Journal of Nonlinear Science, vol. 9, no. 4, pp. 444–447, 2010.

G. A. Baker, Essentials of Padé Approximants, Academic Press, London, UK, 1975.

Y. Yamamoto, C. Dang, Y. Hao, and Y. C. Jiao, “An aftertreatment technique for improving the accuracy of Adomian's decomposition method,” Computers and Mathematics with Applications, vol. 43, no. 6-7, pp. 783–798, 2002.

N. H. Sweilam and M. M. Khader, “Exact solutions of some coupled nonlinear partial differential equations using the homotopy perturbation method,” Computers and Mathematics with Applications, vol. 58, no. 11-12, pp. 2134–2141, 2009.

S. Momani, G. H. Erjaee, and M. H. Alnasr, “The modified homotopy perturbation method for solving strongly nonlinear oscillators,” Computers and Mathematics with Applications, vol. 58, no. 11-12, pp. 2209–2220, 2009.

S. Momani and V. S. Ertürk, “Solutions of non-linear oscillators by the modified differential transform method,” Computers and Mathematics with Applications, vol. 55, no. 4, pp. 833–842, 2008.

P. Tsai and C. Chen, “An approximate analytic solution of the nonlinear Riccati differential equation,” Journal of the Franklin Institute. Engineering and Applied Mathematics, vol. 347, no. 10, pp. 1850–1862, 2010.

A. E. Ebaid, “A reliable aftertreatment for improving the differential transformation method and its application to nonlinear oscillators with fractional nonlinearities,” Communications in Nonlinear Science and Numerical Simulation, vol. 16, no. 1, pp. 528–536, 2011.

J. Biazar, M. Shahbala, and H. Ebrahimi, “VIM for solving the pollution problem of a system of lakes,” Journal of Control Science and Engineering, vol. 2010, Article ID 829152, 6 pages, 2010.

J. Biazar, L. Farrokhi, and M. R. Islam, “Modeling the pollution of a system of lakes,” Applied Mathematics and Computation, vol. 178, no. 2, pp. 423–430, 2006.

Ş. Yüzbaşı, N. Şahin, and M. Sezer, “A collocation approach to solving the model of pollution for a system of lakes,” Mathematical and Computer Modelling, vol. 55, no. 3-4, pp. 330–341, 2012.

M. Merdan, “Homotopy perturbation method for solving modelling the pollution of a system of lakes,” Fen Dergisi, vol. 4, no. 1, pp. 99–111, 2009.

M. Merdan, “He's variational iteration method for solving modelling the pollution of a system of lakes,” Fen Bilimleri Dergisi, vol. 18, pp. 59–70, 2009.

M. Merdan, “A new application of modified differential transformation method for modelling the pollution of a system of lakes,” Selçuk Journal of Applied Mathematics, vol. 11, no. 2, pp. 27–40, 2010.

H. John, “Lake Pollution Modelling,” Virginia Tech

W. H. Enright, K. R. Jackson, S. P. Norsett, and P. G. Thomsen, “Interpolants for Runge-Kutta formulas,” ACM Transactions on Mathematical Software, vol. 12, pp. 193–218, 1986.

E. Fehlberg, “Klassische Runge-Kutta-Formeln vierter und niedrigerer Ordnung mit Schrittweiten-KONtrolle und ihre Anwendung auf Waermeleitungsprobleme,” Computing, vol. 6, pp. 61–71, 1970.