Abstract

A kinetic model, which involves emigration, immigration, birth, death, and fluctuation terms, is utilized to investigate the evolution of city-size distribution. The Boltzmann-type equation and the corresponding Fokker–Planck equation are obtained to analyze the urban size distribution. Three different population variable functions, containing both birth and death rates, are considered. For each population variable function, the closed form of its stationary solution is derived. It is found that the size of urban population distribution follows a power law. Numerical simulation illustrates the correctness of the results.

1. Introduction

A number of works have shown that city-size distribution follows power laws. Benguigui and Blumenfeld-Liebertha [1] create a rank-size logarithmic plot’s equation and several models to describe different types of city-size distribution. Their results illustrate that the distribution of city-size obeys a power law distribution or the concave distribution on the logarithmic graph of the rank size. Based on the city-size data of the southeast and southwest in the United States, Garmestani et al. [2] test the size category and power law behavior of the urban size distribution. According to the size distribution of France communes, Calderın-Ojeda [3] obtains a conclusion that the upper quartile of the commune size data in the early period (1962–1975) strictly follows the power law. On account of the Chinese city-size data of 22 provincial units, Li and Zhang [4] verify that the urban population conforms to the power law distribution through statistical methods. Taking China as an example, the empirical research on the power law distribution of urban rank size is presented in the study of Wu and Yang [5].

Another classical explanation of the city-size distribution is Zipf’s law [6], which is a power law with an index equal to 1. Marsili and Zhang [7] establish a master equation according to the exchange rate of the increase or decrease of a city’s population size. Ghosh et al. [8] introduce a resource utilization model to study the city-size distribution. On account of the data from India’s population census and China’s population census, Gangopadhyay and Basu [9] prove that the urban size distribution of India and China satisfies Zipf’s law. Eeckout [10, 11] concludes that the distribution of urban population is in line with the lognormal distribution. González-Val [12] analyzes the size distribution of Spanish cities from a new perspective, focusing on the role of distance. Levy [13], Zhang et al. [14], Xu and Zhu [15], and Giesen et al. [16] investigated the city-size distribution from different angles.

Gabaix [17] creates models of random growth with Gibrat’s law in the upper tail to explain Zipf’s law for cities. Zanette and Manrubia [18] utilize a stochastic model that incorporates the essential mechanisms with govern city formation and find that population distribution follows a power law. Based on the Indian population census data in the year 2011, Devadoss et al. [19] use statistical approaches to conclude that lower tail cities in India obey a power law. Utilizing a kinetic model, we find that the city-size distribution obeys a power law.

Combining the immigration and emigration factors, Gualandi and Toscani [20] utilize a kinetic model to investigate the features of urban agglomerations. The number of urban population is represented by the variable . The interaction equation for the change in urban population reads as follows:where is a random variable. The functions and describe the rates of people who move out and move in the city. The variable denotes the quantity of population that can migrate toward a city from out of the city. When the urban population alters, it goes from to .

Considering a certain migration function and using the Boltzmann-type equation and corresponding Fokker–Planck equation, Gualandi and Toscani [20] find that the equilibrium density follows a power law for large cities and a lognormal law for middle and low cities.

Inspired by the work of Gualandi and Toscani [20], we consider the interaction rule equation, which includes immigration, emigration, birth, death, and fluctuation terms. A kind of Boltzmann equation is employed to describe the time evolution of population density by using the rarefied gas dynamics theory. We obtain the equation that represents the time–behavior of the average urban population. The Fokker–Planck equation relating to the Boltzmann-type equation is derived by using the asymptotic limit method. Because the birth and death are unfixed, we consider a population variable function that involves birth and death. We consider three different population variable functions, which correspond to different Fokker–Planck equations, from which we obtain an explicit form of the stationary solution. It is found that the size of urban population distribution follows a power law distribution.

Our work is different from that of Gualandi and Toscani [20]. First, the interaction rule in the study of Gualandi and Toscani [20] contains immigration, emigration, and random fluctuation factors, while the interaction rule in our work includes immigration, emigration, birth, death, and random fluctuation factors. Second, due to different interaction rules, the Boltzmann-type equation is not the same as that of Gualandi and Toscani [20]. Third, we consider three different population variable functions about birth and death rates, which are not investigated by Gualandi and Toscani [20]. Finally, the corresponding Fokker–Planck equations are different, and the steady-state solutions are distinct from those of Gualandi and Toscani [20]. Here, we state that our obtained results are similar to those presented by Gualandi and Toscani [20]. Namely, the city-size distribution still obeys a power law.

The structure of this work is as follows: Section 2 provides the kinetic model and its basic parameter settings. The model describing the evolution of the urban population is considered in Section 3. A Fokker–Planck equation is discussed in Section 4. Section 5 presents the numerical analysis of our results. Section 6 is the conclusion.

2. Kinetic Model

Based on the methods in the study of Gualandi and Toscani [20], we assume that the city population change is related to the factors about immigration, emigration, birth, death, and random fluctuation. Although the number of urban population is obviously a natural number, we assume that in the rest of the paper. The interaction equation for the change of urban population evolves as follows:

The left-hand side denotes the population after the change. The first item on the right-hand side indicates the initial population of the city. The functions and denote the rates of people who move out and move in the city. The variable represents persons who migrate toward a city from out of the city. And the forth item shows the change in population involving birth and death rates, where is a population variable function. Assume that is a random variable with mean zero and variance .

Birth and death can be described by using statistical data or stochastic processes [21]. Here, we assume that birth and death are only related to the size of the urban population. Suppose that takes the following form:where and are functions describing the birth and death rates, respectively.

The total number of people moving out cannot be greater than the original population of a city. Thus we assume that , where and . Similarly, we suppose that the upper and lower bounds of and are , , where and .

Next, suppose that the population distribution in the external environment is known. is the probability density function of variable . We suppose that is bounded variable and satisfies as follows:

In Equation (4), indicates the average population which migrates toward a city from out of the city.

We use the dynamic model of particle collision to describe the changing characteristics of urban population distribution. Let denote the probability density function of city-size, which is related to time and variable . It is generally assumed that satisfies as follows:

Note that the change of obeys a weak form of Boltzmann-type equation, namely, for all smooth functions , as follows:where the bracket represents the mathematical expectation, and the left side of Equation (6) reflects the population distribution evolution of the city caused by the interaction at time . Choosing , the right-hand side is zero, which is verified that Equation (6) satisfies the conservation property. The average number of urban population changes over time, namely,

Choosing in Equation (6), we get the evolution of the mean value as follows:

On the basis of the previous discussion, using the interaction rule (Equation (2)) derives as follows:

Substituting Equation (9) into Equation (8) yields as follows:

According to the assumption , , and . Then the mean satisfies the following inequality:

Following the idea in the study of Tchorbadjieff and Mayster [22], we suppose that the average size of a city’s population does not exceed the maximum value , where is defined as follows:where the constants and satisfy .

If , , and are constants, namely , , and . Equation (10) becomes the following:

Solving Equation (13), we acquire the expression of , which reads as follows:

If , from Equation (14), when , the city’s average population converges to .

On the other hand, we consider that and are constants. In the Logistic growth model, , where is the inherent growth rate, and , where is the maximum population capacity of the city. Then , Equation (10) becomes the following:

According to Jensen’s inequalityfrom which we have the following:

Because the right-hand side of inequality (Equation (17)) contains the term , we obtain that the right-hand formula holds if and only if does not cross the bounded value as follows:

Finally, we obtain the following:

Note that the city’s average population is bounded.

3. The Fokker–Plank Equation

In the ordinary way, it is difficult to find the analytical solution of Equation (6). A feasible approach is to investigate its large-time behavior, which is found in the study of Zhou et al. [23] and Gualandi and Toscani [24], to consider the kinetic model. We use the asymptotic limit method to derive the Fokker–Planck model relating to the Boltzmann-type equation. We employ the scaling in the study of Stanley et al. [25] as follows:where and are positive constants, and .

The scaling causes a small change in the average urban population in Equation (13), but this change is independent of the value of the exponent . According to the proportional setting in Equation (20), the average urban population satisfies as follows:

Considering , we obtain the following:

Note that the scale parameter is not included in Equation (22). Therefore, the variables in Equation (8) can be simplified by scaling. The properties of Equation (22) remain unchanged when the interaction rule (Equation (2)) changes slightly.

Using to represent the change of urban population as follows:then the interaction rule (Equation (2)) satisfies as follows:and

Assuming . After scaling, we obtain . Then, we have the following:and

Using Taylor expansion of around , we write as follows:where , .

Using Equations (26) and (27), we obtain the following:where

Substituting Equation (29) into Equation (6), after the time scaling , we obtain the following:where

Considering the remainder term , we get the following:

Thus, we obtain the following:

Applying integration by parts and the boundary conditions vanishing at infinity, we derive the following Fokker–Planck equation:

4. Steady-State Distribution of City-Size

Generally speaking, it is difficult to obtain the explicit solution of the Fokker–Planck equation (Equation (35)). However, we recognize that the Fokker–Planck equation exists as a steady-state solution, which keeps the main feature of the Boltzmann-type equation. As usual in the kinetic model, it is convenient to study certain asymptotics, which leads to simplified models of the Fokker–Planck equation. By means of this way, it is to investigate stationary state solution which retains the same information on the microscopic interaction at a macroscopic level [26]. Thus, we consider the steady-state solution of the Fokker–Planck equation (Equation (35)).

Assume that , , and are constants, namely . Equation (35) becomes the following:

Integrating both sides of the above, we obtain the following:

Setting , we have the following:

Solving Equation (38), we obtain the explicit stationary solution as follows:

In Equation (39), the positive constant is to ensure that the probability density integral in steady state is one, namely . Considering , we have and . From Equation (39), we obtain the following:

Using yields as follows:from which we get the following:

The rank of in a city is as follows:where is a constant. In Zipf’s law, the rank .

According to Equation (39), it is found that the steady-state density is an inverse gamma function. The equilibrium density has a polynomial rate of decay about , namely,

The rank of in this situation is as follows:

When , namely the classical Zipf’s law is obtained for . However, Zipf’s law generally does not hold in the study of Verbavatz and Barthelemy [27].

If and are constants and is variable. In general, the number of births and deaths is not fixed. The parameters considering birth and death should vary. Ignoring other influencing factors, we discuss that the birth rate and death rate of a city are only related to the number of population.

Case 1. According to the Logistic growth model,  and . Then , Equation (35) is written as follows:Solving Equation (46), the steady-state solution takes the following form:In Equation (47), the positive constant has been chosen to fix . From Equation (47), it is a polynomial rate decay with respect to at infinity. The steady density function (Equation (47)) has a power law tail. Ifit follows that the term in Equation (47) decays more slowly than in Equation (39).

Case 2. (a)We know that the population growth gradually flatten out when becomes very big. Thus, we assume the following:where is the inherent growth rate. We write Equation (35) as follows:Usingfrom Equation (50), we have the following:In Equation (52), the constant is to ensure that the probability density integral in the steady-state is 1. From Equation (52), it is found that the steady-state density function is a polynomial rate decay with respect to at infinity and follows a power ().(b)Consider that , namelyEquation (35) becomes the following:from which we obtain the following:In Equation (55), the constant assures that the integral of the equilibrium density is 1. From Equation (55), we find that the equilibrium density function is a polynomial rate decay with respect to at infinity and follows a power (). Whether is positive or negative, the steady-state solution has the same term ().

Case 3. (a)Suppose that is represented as follows:Equation (35) is expressed in the following form:Combiningwe obtain the explicit steady-state density function as follows:In Equation (59), the constant is chosen to satisfy that the integral of the equilibrium density is one. It is polynomial rate decay with respect to at infinity. From Equation (59), we know that steady-state distribution follows power law distribution.(b)Suppose that , namelyEquation (35) becomes the following:from which we get the following:In Equation (62), the constant guarantees that the integral of the steady-state density is one. From Equation (62), we conclude that the equilibrium density obeys a power law. The equilibrium density (Equation (62)) is a polynomial rate decay with respect to at infinity. Equations (59) and (62) have the same term ().

Remark 4. From Equations (39), (47), (52), (55), (59), and (62), we conclude that as goes to infinity, the city-size density function is a polynomial rate decay and obeys a power law. We also say that the city-size distribution follows a power law (see [5, 8, 20]).

5. Numerical Simulation

The empirical study about the Chinese city rank-size obeys a power law distribution by using statistical methods [5]. Thus, theoretically speaking, in this part, we take different values of the parameters of the stationary city-size density and provide several numerical results. When the population growth rate is constant or variable, the steady-state density function of the urban population is a polynomial rate decay with respect to at infinity. Therefore, Equation (39) is used to describe the asymptotic behavior at steady state. According to Equations (39) and (42), we only need to consider the values of , , , , and .

First, we consider different stationary city-size densities of urban population growth rate. When the urban population growth rate ranges from −0.29 to 0.2, the probability density of the steady state is shown in Figure 1. In Figure 1(a), corresponds to (blue line), (black line), and (red line). In Figure 1(b), the values of are 0.01, 0.1, and 0.19, corresponds to (blue line), (black line), and (red line). In Figure 1(c), the stationary profile corresponds to (blue line), (black line), and (red line). Figure 1(d) contains (blue line), (black line), and (red line). From Figure 1, we obtain the conclusion that the increase of , depicting the urban population growth rates, results in stationary density possessing the lower vertex and fatter tails.

Second, we consider the probability density of the urban population size corresponding to different urban population emigration rates when the urban population growth rate . When the emigration rate of the urban population is 0.02, 0.04, 0.06, 0.2, 0.4, and 0.6. The probability density of the steady-state is shown in Figure 2, which shows that the lower the emigration rate of the urban population is, the higher the vertex of probability density is. Lower emigration has fatter tails.

Third, we consider the steady-state density of urban population distribution corresponding to the difference of variance when the urban population growth rate and emigration . When the variance is 0.3, 0.4, 0.5, 0.6, 0.8, and 1, the probability density of the steady urban population is shown in Figure 3. We conclude from Figure 3 that the larger the variance is, the higher the vertex of probability density is. Meanwhile, lower variance possesses fatter tails.

6. Conclusion

The evolution of the urban population is related to many factors, including urban immigration, emigration, birth, and death. We use a dynamic model to study the city-size distribution involving urban immigration, emigration, birth, and death. We describe the interaction rule of urban population change and give the Boltzmann-type equation, from which we obtain the corresponding Fokker–Planck equation. We discuss three different population variable functions and obtain stationary city-size density functions. It is concluded that when the population growth rate is a constant or a variable function, the equilibrium city-size distribution follows a power law. Numerical simulation illustrates the results.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This research is supported by the National Natural Science Foundation of China (no. 11471263).