{'fraction_non_alphanumeric': 0.029649793996468513, 'fraction_numerical': 0.011477339611536198, 'mean_word_length': 4.7842553191489365, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Global Optimization on an Evolving Energy Landscape

Locating the global minimum of a complex potential energy surface is facilitated by considering a homotopy, namely a family of surfaces that interpolate continuously from an arbitrary initial potential to the system under consideration. Different strategies can be used to follow the evolving minima. It is possible to enhance the probability of locating the global minimum through a heuristic choice of interpolation schemes and parameters, and the continuously evolving potential landscape reduces the probability of trapping in local minima. In application to a model problem, finding the ground--state configuration and energy of rare--gas (Lennard--Jones) atomic clusters, we demonstrate the utility and efficacy of this method.

Introduction: Global optimization problems [1] can often be formulated in terms of finding the minimum (or maximum) of a multidimensional potential energy surface (PES). Such problems, which occur in a variety of areas, are of considerable practical and theoretical interest [2]. The "energy landscape" [3] paradigm is particularly useful when the potential energy function is continuously varying with the physical configurations relevant to the problem . An example of such a situation is the proteinfolding problem [4], namely determining the native configuration of complex molecules given their atomic composition. A simpler variant is the determination of the ground state configuration of atomic or molecular clusters [5].
In this Letter we propose a new homotopy method to study such problems by a controlled deformation of the potential energy surface. If V f is the potential energy hypersurface under consideration, we study the landscapes with α a parameter. Given a choice of initial potential, V i , this is a 1-parameter family of potential energy surfaces which smoothly evolves from V i into V f as α varies from 0 → 1. The minima of the landscapes continuously change with α, and in order to track them, one of two strategies are possible. Varying the interpolation parameter α in a finite number of steps, a standard technique such as conjugate gradient (CG) minimization [6] can be employed at each α. On the other hand, one can consider α as a time-dependent function such that the PES evolves according to where h(t) is suitably chosen with h(0) = 0, and lim t→T h(t) → 1. Over a timescale T , therefore, the potential deforms from the initial to the desired potential energy surface, and the evolving minima can be tracked, for example, by following the damped dynamics in this potential via molecular dynamics (MD) simulation.
In the present work we follow both these strategies, and show how homotopic deformation facilitates location of the global minimum in a model problem. Similar (so-called "continuation") homotopic methods have frequently been employed in related situations, as for example in finding roots of polynomial equations in several variables [7] or in the mean-field dynamics in attractor neural-networks [8].
Different global optimization methods frequently find optimal solutions by elimination, by seeking lower and lower minima. Trapping in local minima-and escape from these minima-is a major practical issue. A number of different strategies have been suggested in order to engineer escape from local minima. These include both techniques to allow for large excursions in the phase space by the use of temperature or similar auxiliary parameters (such as simulated annealing [10] and its variants [11,12]) as well as methods that deform the potential energy surface. The diffusion equation method [13] and the distance scaling method [14] fall in this latter class. Other methods utilize both strategies, as for example the stochastic tunneling method [15] where simulated annealing is performed on a surface where the barriers are exponentially reduced so as to facilitate escape from local minima, the landscape paving technique [16], or the basin hopping technique [17] which replaces the potential surface by a set of piecewise flat regions.
The present technique is in the class of optimization methods that exploit potential surface deformation to avoid trapping in local minima. The interpolation parameter α, or the switching functions h(t) smoothly convert one PES into another. The intermediate potentials are qualitatively not very different from the asymptotic potential in terms of the number of minima and maxima, although the relative depths and curvatures are quite different. As we discuss below, this feature contributes to efficiency of the present technique in locating minima. The lowest energy achieved when an ensemble of suitably compact initial configurations is evolved is taken as the ground state prediction of this method.
Application: The problem of minimum energy configuration determination for N particle atomic clusters is computationally hard, and the validity of a global solution cannot, typically, be verified. Existing data for global minima [19] are usually the "lowest minima as yet located" in all but the simplest cases. A variety of global optimization techniques have been applied to this problem [17,18] with differing degrees of success.
For the most extensively studied such systems, namely model rare-gas clusters, the potential energy surface (PES) is an additive pairwise Lennard-Jones interaction, where r ij is the distance between particles i and j, and ǫ, σ are the standard Lennard-Jones parameters. The potential energy landscape varies greatly with cluster size.
Notable difficult optimization problems in this regard are, for example, 38, 75, or 98 atom clusters, where the potential energy surface has the so-called multiple funnel structure [20,21,22].
In the implementation of the MD approach we proceed as follows. V i is taken to be a pairwise sum of harmonic terms V (r ij ) = (r ij − 2 1/6 σ) 2 /2. We perform molecular dynamics simulations [23] of the N particle system, with an additional damping term for each particle, where m is the mass of the particle and γ is the damping coefficient. The internal timescale of interparticle vibrations depends on the parameters m, σ and ǫ. For a given switching function h(t) (we have explored a variety of such functions listed in Table I) the adiabatic time scale is set by the parameter ζ; the entire system dynamics thus has two external time scales ζ −1 and mγ −1 . In the limit γ → ∞, our procedure reduces to a steepest descent minimization on the evolving potential. The dynamics of the system is followed until a stationary configuration is reached.
In order to quantitatively assess the efficiency of this procedure, we define the measure P g = Number of ground state configurations Total number of condensates , a condensate being a configuration such that all atoms are within a single cluster. This is clearly a function of γ and ζ. For the ground state energy, comparison is made to the existing benchmark calculations already available for Lennard-Jones clusters [19].
In the CG approach, V i is taken to be β N j=1 ( r j − r 0 j ) 2 , r 0 j being the (random) initial position for the jth atom. This choice of V i ensures that the initial configuration is the exact global minimum for the potential energy surface, Eq. (1) with α = 0; β is a constant that tunes the curvature of the PES. The parameter α is then varied from 0 to 1 in N s discrete steps; the result of the CG minimization (we follow the standard method [6]) at each step is taken to be the starting configuration for the CG minimization at the next value of α. In this latter approach, therefore, the attempt is to allow the global minimum itself to evolve homotopically. Results: The present application is intended to be illustrative rather than exhaustive. We have systematically studied different cluster sizes up to N = 40 and in all cases the calculated ground-state energy and configuration matches existing results exactly. This includes the difficult case of the 38-atom cluster which is an interesting and important test for any optimization method [20]. The number of minima increases exponentially with cluster size [24]; for LJ 7 there are 4 minima, while for LJ 55 the number exceeds 10 10 [2]. Detailed results, which clarify some aspects of the present technique are presented for the cases of N =19,22 (MD) and N =38 (CG).
In the MD version of the present technique, in the absence of switching, namely in the sudden limit V (t) = V f , the system quickly settles into the nearest available minimum based on the level of damping introduced. By starting from an ensemble of initial conditions, a variety of different minima are reached but the probability of finding the true ground state is essentially zero for large clusters.
With an adiabatic switch [9], the results are dramatically different. The continuous evolution of the potential energy landscape is a key factor in permitting escape from local minima. Only asymptotically does the system come to rest, but until then, there is always residual kinetic energy due to which the system avoids being trapped by small barriers. Shown in Fig. 1 is the typical variation of potential energy (in units of ǫ), which is nonmonotonic once the adiabatic switching is incorporated. Regardless of the actual form of the switching, more than 85% of all initially random configurations condense, except in the case where the switching is applied to the repulsive term of the potential. Representative data is given in Table I. As emphasized, the adiabatic optimization proposed here is heuristic. The optimal choice for the parameters γ, ζ for a given cluster size depend on a number of features such as the interaction potential parameters and the inherent time-scales. By scanning over reasonable values of the parameters, it is possible to determine regions in parameter space with a higher than average probability of reaching the ground state. It also appears that adiabaticity is crucial since the probability of reaching the ground state increases substantially with decreasing ζ: P g is shown versus ζ for the 19-atom case in Fig. 2.
In the CG method of following minima during homotopy, the probability of reaching the global minimum is enhanced through the modification of the PES curvature. Since V i adds a uniform positive curvature at the intermediate stages it effectively suppresses or eliminates many barrier. To perform some benchmarking of the advantage this gives, we present, in Table II, data pertaining to finding the global minimum for LJ 38 comparing the present method and the basin-hopping technique. The three lowest minima are at energies -173.928, -173.252 and -173.134 respectively. In either method, all particles are initially placed randomly inside a sphere of radius ( N 4 ) 1/3 . In V i , the parameter β= 100. In our implementation of the basin-hopping algorithm, coordinate displacements are random in the interval [-0.3,0.3] and the temperature is taken to be 2. An overall confining potential of the form V c = i exp(20(r i − a)), a = 1 + ( N 4 ) 1/3 was added to prevent dissociation, and a standard Polak-Ribiere algorithm was used for the conjugate-gradient minimization [6] with tolerance set between 10 −5 and 10 −7 . The average computational effort required is a product of the number of trials needed in order to get to the ground state on average and the number of function and derivative calls per trial. In our implementation of the algorithms, we find that the reduction in computational effort in locating the global minimum through the homotopy method is about 40%. The relative efficiencies can, however, vary depending on the actual choice of the various adjustable parameters in the two techniques. In either the MD or the CG version, configurations that do not reach the global minimum still typically tend to find the lowest energy states, so that a by-product of this methodology is a considerably detailed map of the low excitation regime of the cluster. This feature, however, is not unique to the present method.
Summary: We have presented here a method for global optimization which relies on the guided evolution of the underlying landscape. The methodology for finding minima on this surface can vary, and in the examples presented here, we have used both the conjugate gradient technique as well as damped molecular dynamics. (Dynamics in the landscape has been incorporated in other techniques, for example in genetic algorithms [25].) As in other methods, apart from the global minimum, we also obtain a detailed picture of the excitation spectrum.
Within the context of cluster geometry determination itself, several issues need to be addressed. The adiabatic method can be shown to locate ground states even when there are bifurcations along the deformation pathway [26]. Is it possible to design more efficient homotopic deformations? What is the role of V i in controlling the efficiency?
The application here, though in some ways a model problem, has all the complications that arise in more general optimization problems. The success of this simple technique is therefore encouraging.  Table II: Comparative analysis of the homotopy method and basin-hopping. For each method, N r initial configurations are evolved to find the global minimum in 100 instances for the LJ 38 cluster. N j , j = 0, 1, 2 are the number of times the lowest three minima are found in the two methods; the number of function and derivative calls needed (per initial condition) are also indicated to give an estimate of the computational effort involved.