#### Abstract

We present a fourth order finite difference method for doubly singular boundary value problem with boundary conditions and , where , , and are finite constants. Here and is allowed to be discontinuous at the singular point . The method is based on uniform mesh. The accuracy of the method is established under quite general conditions and also corroborated through one numerical example.

#### 1. Introduction

Consider the following class of singular two point boundary value problems: with boundary conditions where , and , are finite constants. Here and is unbounded near . The condition says that the problem is singular and as is unbounded near , so the problem is doubly singular [1]. Let , , and satisfy the following conditions:â€‰(i) in ,(ii) ,(iii) ,,, and(iv) and exists on ,â€‰â€‰where on .â€‰(i) in , is unbounded near ,(ii) ,(iii) , with ,(iv) and exists on , and(v) . â€‰ â€‰ is continuous on , exists, and is continuous and for all and for all real .

Thomas [2] and Fermi [3] independently derived a boundary value problem for determining the electrical potential in an atom. The analysis leads to the nonlinear singular second order problem with a set of boundary conditions. The following are of our interest:(i)the neutral atom with Bohr radius given by ;(ii)the ionized atom given by .Furthermore, Chan and Hon [4] have considered the generalized Thomas-Fermi equation for parameter values ,,, and .

Such singular problems have been the concern of several researchers [5â€“8]. The existence-uniqueness of the solution of the boundary value problem (1.1) with boundary condition (1.2) or (1.3) is established in [1, 9â€“13]. Bobisud [1] has mentioned that in case , the condition is quite severe, that is, it is sufficient but not necessary for forcing the solution to be differentiable at . In fact if , then Thus if either or . But if has discontinuity at ; thus it is natural to consider the weaker boundary condition .

There is a considerable literature on numerical methods for but to the best of our knowledge very few numerical methods are available to tackle doubly singular boundary value problems. Reddien [14] has considered the linear form of (1.1) and derived numerical methods for which is stronger assumption than .

Some second order methods (Chawla and Katti [15], Pandey and Singh [16, 17]) as well as fourth order methods (Chawla et al. [18â€“21], Pandey and Singh [22â€“24]) have been developed for . Most of the researchers have developed methods for the function , and the boundary conditions and .

Chawla [19] has given fourth order method for the problem (1.1)-(1.2) with and . The method is extended by Pandey and Singh [22] to a class of function satisfying with analytic in the neighborhood of the singular point.

In this work we extend the fourth order accuracy method developed in [22] to doubly singular boundary value problems (1.1)-(1.2) and (1.1), (1.3) both. That is, we allow the data function to be discontinuous at the singular point . Further, we do not require analyticity of the function . The convergence of the method is established under quite general conditions on the functions , and . Fourth order convergence of the method is corroborated through one example and maximum absolute errors are displayed in Table 1.

The work is organized in the following manner. In Section 2, first the method is described and then its construction is explained. In Section 3, convergence of the method is established and in Section 4, the order of the method is corroborated through one example.

#### 2. Finite Difference Method

This section is divided in to two parts: (i) description of the method and (ii) derivation of the method. The coefficients not specified explicitly in this section are specified in the appendices.

##### 2.1. Description of the Method

In this section we first state the method; detailed derivation is given in Section 2.2.

For positive integer , we consider the mesh over :â€‰, . For uniform mess , where is mentioned in condition . Let â€‰â€‰ is the solution, and so forth. Now we approximate the differential equation (1.1) with boundary condition (1.2) on the grid by the following difference equations: for boundary value problem (1.1)-(1.2). In the case of boundary condition (1.3) we use the following difference equation: Here denotes the approximate solution,,, and Note that the functions and are defined in assumption . Thus the method for boundary value problem (1.1)-(1.2) is given by difference equations (2.1)-(2.2). Further, for boundary value problem (1.1) and (1.3) the method is given by difference equations (2.1) (with ) and (2.2)-(2.3).

##### 2.2. Construction of the Method

In this section we describe the derivation of the method. Local truncation errors are mentioned without proof.

With the differential equation (1.1) becomes . Now integrating twice, first from to and then from to and changing the order of integration we get where , and . In an analogous manner, we get Eliminating from (2.5) and (2.6) we obtain the Chawla's identity where Via Taylor expansion of , and about in and taking the following approximations for and we get the approximation for the smooth solution as: where is local truncation error given by The coefficients , and so forth, are specified in the appendices and the functions , are defined in the assumptions , .

###### 2.2.1. Discretization of the Boundary Condition at

We write (2.6) for

Now, we use the boundary condition at , Taylor expansion of and following approximations for and in (2.12) and get the discretization at as follows: where the local truncation error is given by The coefficients are specified in the appendices.

The discretization (2.14) involves unknown , which we approximate in the following way: where Now, let then replacing by in (2.14) we obtain the following discretization for a smooth solution at where

###### 2.2.2. Discretization of the Boundary Condition at

Integrating the differential equation (1.1) twice, first from to ; then from to and interchanging the order of integration we get Via Taylor expansion of ,, and about we get the following discretization for : where the local truncation error is given by The coefficients are specified in the appendices and the functions , are defined in the assumptions and .

The discretization (2.22) involves ; so for we use the following approximation: where Let ; then in (2.21) replacing by we obtain the following discretization for : where , and coefficients , are specified in appendices.

#### 3. Convergence of the Method

In this section we show that under suitable conditions the method (2.1)-(2.2) for the boundary value problems (1.1)-(1.2) described in Section 2 is of fourth order accuracy.

Let ,,, and , then the finite difference method given by (2.1)-(2.2) for the boundary value problems (1.1)-(1.2) can be expressed in matrix form as and corresponding to smooth solution satisfies the perturbed system where and are tridiagonal matrices.

Let , then from Mean Value Theorem ,, then from (3.1) and (3.2) we get the error equation as

We recall that by the notation we mean that all components of the vector satisfy . Similarly by we mean all the elements of the matrix satisfy .

From Corollary of Theorems 7.2 and 7.4 in [25] an (irreducible) tridiagonal matrix is monotone if Furthermore, from Theorem 7.3 in [25], exists and the elements of are nonnegative. Now, if is monotone and , then from Theorem 7.5 in [25] Let vector norm and matrix norm be defined by then from (3.3) we get where .

Furthermore, if and is also monotone matrix, then from (3.5) we get Now from [22], of matrix is given by where Let , for , and be bounded on . Now, as is fixed in and is chosen such that , for sufficiently small we get for suitable constants and .

Now from (3.8)â€“(3.11) we get where . It is easy to see that From (3.12)-(3.13) we obtain where In view of the following inequalities for , and from (3.14) it is easy to establish that Thus we have established the following result.

Theorem 3.1. *Assume that , , and satisfy assumptions given in ,, and , respectively. Let , then the method is of fourth order accuracy for sufficiently small mesh size provided, , and are bounded on . *

*Remark 3.2. * We have developed the numerical method for the boundary condition but could not establish the fourth order convergence although the order of the accuracy is verified through one example in the next section.

#### 4. Numerical Illustrations

To illustrate the convergence of the method and to corroborate their order of the accuracy, we apply the method to following example. The maximum absolute errors are displayed in Table 1.

*Example 4.1. *Consider

with exact solution .

The method is applied on the example for . Maximum absolute errors and order of accuracy for Example 4.1 are displayed in Table 1.

We use the following approximation of in our numerical program:

##### 4.1. Numerical Results for Uniform Mesh Case

*Remark 4.2. *In the discretization (2.18), may take any value in but from the Table 1, maximum absolute errors are less for the value of close to one. So it is better to take value of close to one.

#### Appendices

The coefficients involved in the method and their approximations are mentioned below.

#### A. Coefficients for Method

where