We consider boundary value problems with an essential singularity (singularity of the second kind),
tαz'(t)=f(t,z(t)), t∈ (0,1], (1)
z ∈ C[0,1]. (3)
where α is greater than 1. Here, f, g are smooth functions of dimension n and p, respectively. In general p is smaller than n, and condition (3) provides the additional n-p relations that guarantee the well-posedness of the problem. Boundary value problems with an essential singularity are frequently encountered in applications. In particular, problems posed on infinite intervals belong to the above problem class when they are suitably transformed to a finite interval. Flow problems (Blasius equation, von Karman swirling flow) and the classical electromagnetic self-interaction problem, are sources for the models we are interested in. The purpose of the paper is to investigate a numerical approach which may be successfully applied to obtain high-order solutions for boundary value problems of the type (1)-(3). The results are original and have not been published elsewhere. Especially, we examine the empirical convergence order of collocation methods at either equidistant or Gaussian points. The motivation to apply these methods was their satisfactory performance when solving boundary value problems with a singularity of the first kind, where α = 1 in (1). Therefore collocation has been implemented in our MATLAB code sbvp designed for the latter class of problems, together with an a posteriori global error estimate based on the defect correction principle. This error estimate is asymptotically correct when a singularity of the first kind is present, but unfortunately, it does not work for problems with an essential singularity. In this case however, a strategy based on mesh halving seems to be a promising candidate which can provide an asymptotically correct error estimate for the collocation solution of (1)-(3).