The difficulties as well as possible remedies in solving boundary value problems and eigenvalue problems on unbounded domains are thoroughly considered in the well known monograph of Boyd JPB and recent monograph of Shen, Tang and Wang STW . There exists a vast literature on domain truncation type methods, the most used methods for such problem, but considerable less rigorous estimations of the domain truncation errors. Domain truncation was coupled with different techniques such as shooting, finite differences as well as Galerkin schemes. Galerkin or collocation methods, based on the Laguerre functions, and the Chebyshev polynomials with various mappings, are the other two alternatives. The application of these methods restricts to various first and second order differential problems (see the above mentioned works and also BRB and Shen1 ).

Thus, the main aim of this paper is to accurately solve higher order boundary layer type problems, including generalized eigenvalue ones, by a collocation method based on generalized Laguerre functions. The method is general, fairly efficient in imposing the boundary conditions in origin and at infinity, easy implementable and require no domain truncation. Actually, the boundary conditions at infinity are exactly satisfied by Laguerre functions and the conditions at zero are introduced by a removing technique of independent boundary conditions due to J. Hoepffner JH . Once the differentiation matrices are accurate computed the method is directly implementable. Throughout this paper the Laguerre differentiation matrices are those provided by the paper of Weideman and Reddy WR (see also the monograph of Funaro F ).

We will analyze two distinct, difficult, and representative boundary layer type problems. The first one is the general Falkner-Skan two-point boundary value problem. This is a nonlinear third-order one, with not known closed-form solution (see for instance Schlichting HS , Ockendon & Ockendon OO or Acheson Ach ) . It is one of the most challenging problems in viscous fluid dynamics. It appears in the study of the viscous boundary layer on a flat plate.

As a very short history, we observe that in the years ’30 and ’40 of the 20th century important scientist such as Falkner, Skan, Hartree and Weyl tried approximate analytical solution to this problem. Around ’70, Keller KHB and Cebeci et al. CC proposed finite difference solutions, based mainly on the shooting method, for this problem. Riley and Weidman RW , Fang and Zhang FZ and Sachdev et al. SKB were concerned, more recently, with analytical (series solutions) to Falkner-Skan type problems. In spite of these considerable efforts, reliable numerical values of the unknown function and its first two derivatives, on a grid which concentrates its nodes toward the fixed boundary, are not available.

The aim of this work is to show that Laguerre collocation provides numerical solutions that fill this gap. In solving numerically, the famous Falkner-Skan equation, one encounter a supplementary difficulty, namely two genuinely nonlinear terms. To cope with this difficulty the so called homotopy perturbation method was used in the recent years. The papers of Motsa and Sibanda MS (see also the numerous literature quoted) is illustrative in this respect. The authors truncate the domain to a finite one, discretize the equation using Chebyshev collocation and then solve the nonlinear discrete system by a homotopy technique. Unfortunately, no indication about the convergence of this iterative procedure is supplied and the computation is carried out on a rough regular grid.

Up to our knowledge, only in the papers PDP and PDP1 the authors use sinc function representation in order to avoid the truncation of the half-line. The sinc method, in its original form, is intended for solving problems on the real line. With the aid of a mapping function the method may be used on the half line. The authors solve the particular Blasius problem in the first paper, and Falkner-Skan problem in the second paper quoted above, on equidistant nodes. They use a homotopy method to cope with nonlinearity and unfortunately the resolution of their results is poor.

In our strategy, the nonlinear algebraic system produced by Laguerre discretization is solved by a Newton type method. This strategy depends only on a single positive spatial scaling parameter which maps half-line to itself. This parameter is used in order to improve the resolution of the method inside the boundary layer. In this way, the solution of the problem, its first derivative, i.e., tangential velocity as a function of normal distance and its second derivative, i.e., the drag inside the boundary layer are accurately estimated.

Two boundary value problems, one linear and another nonlinear, both supplied with Falkner-Skan boundary conditions are solved as non trivial tests for the accuracy of our algorithm. Using some exact solutions, we observe that the linear problem can be solved with spectral accuracy. Unfortunately, the nonlinear terms in Falkner-Skan problem severely reduce this accuracy.

The second problem considered is the hydrodynamic stability of the Ekman boundary layer adjacent to a rigid surface. This problem has been analyzed in the pioneer work of Lilly DKL , and then thoroughly reconsidered in AB , GM and MVM . The dynamical (physical) aspects of this problem are available in some details in the work of Allen and Bridges AB . From mathematical point of view the problem is a generalized eigenvalue one (GEP for short) containing derivatives up to fourth order. Greenberg and Marletta GM treat the problem in the general framework of the operator theory and then reduce that to a first order system of differential equations. The shooting method coupled with the domain truncation produces numerical results. Allen and Bridges advocate in their paper AB that ”the compound matrix method is clearly the best approach to hydrodynamic stability problems on infinite domains”. In fact they also couple the compound method with a domain truncation. However, we will also provide the leftmost eigenvalue and additionally the most unstable eigenmode but with a less important computational effort. Moreover, in our method the discretized algebraic (GEP) is accurately solved by Arnoldi methods and the leftmost part of the spectrum is comparatively found by Jacobi-Davidson type methods.

The outline of the paper is as follows. In Sect. 2 we quantify the performances of Laguerre collocation when a standard linear third order two-point boundary value problem on the half line is solved. The problem can be seen as a linear variant of the Falkner-Skan equation and we assume for this, solutions that have different typical decay at infinity. In Sect. 3 we obtain the Laguerre collocation equations for the general Falkner-Skan boundary value problem. In Sect. 4 we get the so called Blasius solutions and solve the general equation. Some comments on the accuracy of the method are also provided. In Sect. 5 we solve the Ekman problem and Sect. 6 concludes.

We introduce the method with the representative problem

where the coefficients $a\left(x\right),b\left(x\right)$ and $c\left(x\right),x\in (0,\mathrm{\infty }),$ are smooth enough in order to guarantee a correct application of the strong collocation technique.

The interpolation operators are built on the weighted Laguerre polynomials, i.e., on functions $e^{-x/2}{L}_n\left(x\right)$, $L_n\left(x\right)$ being the classical Laguerre polynomial of order $n$. The usual Laguerre interpolation in weighted Sobolev spaces was recently systematically analyzed in the monograph of Shen et al. STW . However, the best result remains that of Bernardi and Maday BM which reads

where $I_N$ is the usual interpolation operator, the weights $\omega$ are defined as ${\omega }_0:=e^{-x}$, and the norms have usual significance.

As our numerical experiments showed that the Laguerre differentiation matrices are seriously polluted by round off errors for $N>100$ (see also Fig. 3) this result seems to be not very encouraging. Consequently, we try to examine numerically the convergence behavior of some typical functions. Thus we consider the following four exact solutions to (1):

for some natural $k$.

These functions have typical decay properties. The errors corresponding to the first case ($inf$ and $L_2$ norms respectively) are reported in Fig. 1. It is fairly clear the the accuracy depends on the scaling factor $b$. Whenever this factor is properly chosen, even for moderate values of $N$, the method produces results with spectral accuracy.

The presence of an oscillation factor worse considerable the approximation and requires a larger interval. The situation corresponding to the second example is summarized in Fig. 2.

Two important conclusions can be formulated. First, in our boundary layer type problems we expect an exponential decay without oscillations of solutions. Thus the spectral accuracy can be altered only by nonlinearity of the problem. Second, for oscillatory decaying solutions, mapping techniques or even truncation domain methods can be a serious concurrent for Laguerre collocation. In Shen1 , and partially in Shen , the authors study the accuracy of Galerkin method applied to a self adjoint second order elliptic problem with similar conclusions.

We have introduced a general applicable Laguerre collocation type method in order to discretize high order differential equations on the half line. The method handles efficiently any type of boundary conditions and avoids the domain truncation. Moreover it has provided accurate solutions inside the boundary layer as well as in characterizing the linear hydrodynamic stability of Ekman layer. With a simple smoothing process the second derivative of the unknown in the Falkner-Skan equation is correctly recovered on a fine grid clustered toward the fixed boundary.

The method proved to be stable with respect to the order of discretization and scaling factor.

However, a systematic quantitative comparison between our numerical results and the results already obtained, mainly with shooting and truncation, is rather difficult to be carried out because our method is conceptually different.

Qualitatively, the result are similar with the older ones in both cases, a third order nonlinear boundary value problem and a fourth order eigenvalue one. However, the method is directly and easy implementable and, at least in the case of eigenvalue problems, confirm the sharpest theoretical results available, based on the operator theory. More than that, using the notion of pseudospectrum it indicates possible limitations of the linear stability analysis.