@dharr 

Thank you. Sorry. That was my mistake. Sure, I need z(0)=0 instead of z(-infinity)=0. But nevetheless even for z(-15)=-1 I got the same error "Newton iteration is not converging". Did you use Maple 17?

In general case this equation has an analytical soltion in terms of Jacobi ellptic sine that under certain conditions can approxomated by the hyperbolic function like tanh(x). 

P.S. Update. I see that you used an approximation. Thanks!

@dharr 

Many thanks for the detailed repsonse.

This is exactly the result I expected with 1/6 for the positive x and y approaching zero from the positive direction and the positive parameters.

P.S. Sorry for the long reply, I recovered from COVID

It is not classified. Maple 17. The code is the simplest

First, I just take the indefinite integral

restart;
R0 := 1/(a-sqrt(b+c*cos(x)));
R1 := int(R0, x);

Then, bearing in mind the limits of integration, I omit the contribution of the term with the arctan function because of zero and try to simplify the expression without elliptic integrals:

R10 := 2*sqrt(-(2*c*cos((1/2)*x)^2+b-c)*(cos((1/2)*x)^2-1))*sqrt(-cos((1/2)*x)^2+1)*sqrt((2*c*cos((1/2)*x)^2+b-c)/(-c+b))/((a^2-b+c)*sqrt(-2*cos((1/2)*x)^4*c-cos((1/2)*x)^2*b+3*c*cos((1/2)*x)^2+b-c)*sin((1/2)*x)*sqrt(c*(2*cos((1/2)*x)^2-1)+b))

At the end I obtain the analytical result, which I specified in the question by means of complete integrals of the first and the third kind

@mmcdara 

Thanks for the rapid reply and the detailed explanation. But I worry if before I consider the indefinite inegral and obtain the analytical expression for it

2*a*arctan((a^2-b+c)*tan((1/2)*x)/sqrt((a^2-b-c)*(a^2-b+c)))/sqrt((a^2-b-c)*(a^2-b+c))+2*sqrt(-(2*c*cos((1/2)*x)^2+b-c)*(cos((1/2)*x)^2-1))*sqrt(-cos((1/2)*x)^2+1)*sqrt((2*c*cos((1/2)*x)^2+b-c)/(-c+b))*(EllipticF(cos((1/2)*x), sqrt(-2*c/(-c+b)))*a^2-a^2*EllipticPi(cos((1/2)*x), 2*c/(a^2-b+c), sqrt(-2*c/(-c+b)))-EllipticF(cos((1/2)*x), sqrt(-2*c/(-c+b)))*b+EllipticF(cos((1/2)*x), sqrt(-2*c/(-c+b)))*c)/((a^2-b+c)*sqrt(-2*cos((1/2)*x)^4*c-cos((1/2)*x)^2*b+3*c*cos((1/2)*x)^2+b-c)*sin((1/2)*x)*sqrt(c*(2*cos((1/2)*x)^2-1)+b))

then perform simplification and eventually substitute the integration limits 0 and 2*Pi. This is the procedure that gives me that discrepancy.

@mmcdara 

Many thanks for the rapid response and the detailed explanation. It looks clear but I worry that if I consider the indefinite integral before and the obtain analytical result 

2*a*arctan((a^2-b+c)*tan((1/2)*x)/sqrt((a^2-b-c)*(a^2-b+c)))/sqrt((a^2-b-c)*(a^2-b+c))+2*sqrt(-(2*c*cos((1/2)*x)^2+b-c)*(cos((1/2)*x)^2-1))*sqrt(-cos((1/2)*x)^2+1)*sqrt((2*c*cos((1/2)*x)^2+b-c)/(-c+b))*(EllipticF(cos((1/2)*x), sqrt(-2*c/(-c+b)))*a^2-a^2*EllipticPi(cos((1/2)*x), 2*c/(a^2-b+c), sqrt(-2*c/(-c+b)))-EllipticF(cos((1/2)*x), sqrt(-2*c/(-c+b)))*b+EllipticF(cos((1/2)*x), sqrt(-2*c/(-c+b)))*c)/((a^2-b+c)*sqrt(-2*cos((1/2)*x)^4*c-cos((1/2)*x)^2*b+3*c*cos((1/2)*x)^2+b-c)*sin((1/2)*x)*sqrt(c*(2*cos((1/2)*x)^2-1)+b))

then perform simplification and eventually I substitute the integration limits 0 and 2*Pi to evaluate the integral.

From your answer it is obvious that I missed some terms but the problem how to take them into account in this case?

@Rouben Rostamian  

The problem has a simple physical background. In the solid-state physics, it is possible to determine the real and imaginary parts of the dielectric permittivity \epsilon through the polarization function, depending on the frequency omega and a certain parameter q. In fact and roughly speaking, the integral is the same polarization operator, which is defined through the law of dispersion . 

Usually, if one knows the dispersion law (e.g., the most common example is quadratic one), one can predict the frequency behavior of the dielectric permittivity function. I want to solve the inverse problem.

Let us assume that we know the form of the dielectric permittivity. Let us imagine that it has the form of a slowly varying function (https://en.wikipedia.org/wiki/Slowly_varying_function). The task is to determine the corresponding dispersion law. 

I must confess that I chose the integration limits almost arbitrarily. In principle, there are no physical restrictions on other values of integration limits. I wanted to understand for a start how it is technically possible to numerically solve this equation, which can be called as "delay-integral equation". 

@Carl Love 

Sorry for this misunderstanding. I remember you from your helpful responses to my previous questions, so I just wanted to say thank you. I really didn't know that this type (DDE) of equation exists. 
And here I am just looking for how to numerically solve this equation, or at least get a couple of hints that will push me towards solving

@Carl Love 

As usual, thank you very much for your comment and clarification. I had assumed that this type of equation referred to integro-functional equations. 

Above, I provided some additional information about the values of the parameters and the unknown function.

@Rouben Rostamian  thanks for your response.

1. You are right absolutely.  At least, at the begininig I wanted to undertand qualitatve behaviour of this function .

2. I forgot to mention that definitely I was seeking the numerical solution. It would be great if the analytical solution exists even in the form of very cumbersome expression or even in the form of highly nonlinear differential equation. 

As for numerical values of parameters:

a=1, b=100

and  is positevely defined. 

Values of omega are also positive.

@Preben Alsholm 

Thank you fo the answer and for the suggestion. Here for this differential equation a parameter p is an ineteger number. At least for p=0 there is an analytical solution expressed in terms of periodic functions, namely sine and cosine fucntions. So obviously based on that periodic conditions I expect to obtain periodic solutuion also.

I have tried already to use initial function as sin(x)+cos(x) but now I experience another problem "Error, (in dsolve/numeric/bvp) initial Newton iteration is not converging".

@nm Many thank for the reply and for the improvment! It's just my unwanted error with the incorrect call of the type for the solution.

I have chosen n(x) as a function for the convinience to export later this expression to latex. 

@dharr 

Thank you for the comment and for the help. I found the solution (see my comment above).

@Preben Alsholm 

It's ok. I have resolved that issue and I have extracted the data.

Thanks a lot for the reply.

I know that the execution of the code is rather long but my problem is that I can't extract the data file from the plot due to this warning.

@vv  I apologize for the question again where did you find the asymptotic expansion of incomplete elliptic integrals in the case complex argument and modulus? I mean in what textbook?