Key factors affecting the parallel efficiency of archetypical quantum-chemical calculations are discussed. Effective load balancing schemes are proposed. Introduction of the memory affinity to the balancing process is shown to result in super-linear scaling

Brela, Mateusz

Makowski, Marcin

Mazur, Grzegorz

English

Computing and Informatics (E-Journal - Institute of Informatics, SAS, Bratislava)

Computing and Informatics, Vol. 30, 2011, 761–771EFFECTIVE RESOURCE ALLOCATION IN PARALLELQUANTUM-CHEMICAL CALCULATIONSGrzegorz MazurDepartment of Computational Methods in ChemistryJagiellonian UniversityIngardena 3, 30-060 Kraków, Polande-mail: mazur@chemia.uj.edu.plMarcin MakowskiDepartment of Theoretical ChemistryJagiellonian UniversityIngardena 3, 30-060 Kraków, Polande-mail: makowskm@chemia.uj.edu.plMateusz BrelaFaculty of ChemistryJagiellonian UniversityIngardena 3, 30-060 Kraków, Polande-mail: mbrela@student.chemia.uj.edu.plCommunicated by Jacek KitowskiAbstract. Key factors affecting the parallel efficiency of archetypical quantum-chemical calculations are discussed. Effective load balancing schemes are proposed.Introduction of the memory affinity to the balancing process is shown to result insuper-linear scaling.Keywords: Computational chemistry, parallelization, optimization, super-linearscaling762 G. Mazur, M. Makowski, M. Brela1 INTRODUCTIONSignificantly growing demand for higher quality computational chemistry results forlarge molecular systems can be recently observed. As a result, the chemistry commu-nity needs for computational power are constantly growing. While this trend is ac-companied by quick development of high performance hardware solutions, the typicalarchitecture of modern supercomputers and growing popularity of computer clustersrequire algorithmic changes in software to efficiently use the accessible resources.The key issue is that the cluster architecture is inherently highly parallel. Thesame holds for virtually all modern supercomputers. This necessitates the switchfrom traditional sequential algorithms to those that are able to exploit the pa-rallelism of the hardware architecture. Fortunately, most of the typical quantumchemistry problems are relatively easy to reformulate in such a way.Most of the existing procedures used for commonly executed quantum-chemicalcalculations are parallelized within the Single Program Multiple Data (SPMD)model. This approach is in general not optimal from the efficiency point of view. Itis, however, justified by its conceptual simplicity and near-optimal performance fortypical quantum-chemical applications. The SPMD model requires a load-balancingalgorithm to partition the problem into tasks performed by the computational nodes.The partitioning applies to all resources, which in case of typical quantum chemicalcalculations means the CPU time and memory.In the paper we review archetypical quantum-chemical algorithms, analyzinghow their key features affect their parallel performance. This allows us to formulateeffective load-balancing schemes.The proposed load-balancing procedures were implemented in development ver-sion of Niedoida [1], a general-purpose computational chemistry package built asa set of libraries implementing quantum-chemical and microelectrostatic calcula-tions. Niedoida is a parallel program designed according to the SPMD model. Mes-sage Passing Interface (MPI) [2] is used as the low-level parallelization framework.The main reasons for using Niedoida in this work is that it provides a plug-in load-balancing interface and a wide range of production quality quantum-chemical pro-cedures, allowing for benchmarking different algorithms in a realistic and consistentenvironment.The paper is structured as follows. The core computational quantum chemistryalgorithms are introduced in Section 2. In the next section the parallel Hartree-Fock calculations are analyzed. In Section 4 the parallel Moller-Pleset PerturbationTheory implementation is discussed. The paper is concluded with short summaryin Section 5.2 CORE ALGORITHMS OF COMPUTATIONALQUANTUM CHEMISTRYComputational quantum chemistry can be considered to be just a glorified namefor simple tensor algebra. Specifically, these are various transformations of theEffective Resource Allocation 763two-electron integrals (electron repulsion integrals, ERI) tensor which consitute therate-determining step of the calculations. The ERI are defined as(µν|κλ) =∫∫dr1dr2χµ(r1)χν(r1)1|r1 − r2|χκ(r2)χλ(r2) (1)where χ denote the atomic orbitals. It has to be noted that the sheer size of theERI tensor and its sparsity pattern prevent standard algebraic treatment.While the actual resource requirements depend strongly on the type of tensorcontraction performed by a specific algorithm, we can distinguish two main types ofthem, the two-index contractionGµν =∑κλPκλ [2(µν|κλ)− (µλ|κν)] (2)and the four-index transformation(ia|jb) =∑µνλσCµiCνaCκjCλb(µν|κλ). (3)The former consitutes the core of the Hartree-Fock [3, 4] and hybrid Kohn-Sham [5] methods. The latter is crucial in the post-Hartree-Fock methods for whichwe used the Moller-Pleset second order perturbation theory (MP2) [4] as a repre-sentative example.In the next sections we analyze how the structure of the key transformationsaffects the performance of the various load-balancing schemes and propose effictientsolutions.3 PARALLEL HARTREE-FOCK CALCULATIONSThe bottleneck of the two-index contraction of Equation (2) is the two-electron in-tegral generation step, being responsible for roughly 80% of the whole calculationtime. This is caused by both the sheer number of the integrals and by the rela-tively high computational cost of calculating an integral. The latter holds even formodern quasi-optimal algorithms of integral generation. To considerable extent, thesame applies to the Kohn-Sham (KS) calculations with hybrid exchange-correlationfunctionals.The structure of the problem results in the natural partitioning scheme whichsplits the ERI generation into batches and scatters them across the nodes.3.1 CPU-Time Driven Load-BalancingWe start the analysis with the simplest possible load-balancing algorithm, the staticone. It splits the computational problem into equally sized tasks which are uniformlydistributed between nodes. While conceptually easy, the approach is bound to be764 G. Mazur, M. Makowski, M. Brelavery inefficient. Even assuming highly uniform computational enviroment, the timespent by the processors on performing assigned tasks will be significantly different.This is because the time necessary to calculate integrals involving orbitals withvarious angular momentum differs significantly.Recently we proposed an efficient load-balancing algorithm, which dynamicallyadapts to the inhomogeneity of either the computing environment or the data beingcalculated [6]. The proposed algorithm works as follows. The computational prob-lem is divided into tasks of different sizes. The tasks are stored in a task queue.A node gets the next task from the queue as soon as it completes the previous one.The splitting of the problem into tasks is organized as follows. A fraction 1/f of theoriginal problem is divided into n tasks, where f is the splitting factor and n standsfor the number of nodes. Then the procedure is repeated recursively for the remain-ing part of the problem. The recurrence is stopped when the size of the remainingpart is smaller than the threshold t. The algorithm is parametrized by f and t.Test calculations were performed for sexithiophene (see Figure 1) at the Hartree-Fock level of theory in the 6-31G** basis. The size of the model system is bestcharacterized by the number of orbitals. In this case the number of orbitals is 514,which is of the magnitude typical for commonly performed quantum-chemical cal-culations. The scaling of the proposed load-balancing algorithms with the numberof CPUs is shown in Figure 2. The dynamic load-balancing algorithms performsconsistently better than the static one.Fig. 1. The sexithiophene moleculeThe relative time required to perform calculations using n CPUs may be de-scribed byTn =αnβ+ γ (4)where α represents the fraction of the parallelized part and γ stands for the (ef-fectively) serial fraction of the calculations. β describes the deviation from linearspeedup. The scaling function described above was fitted to the datapoints obtainedfor the model system. The fitted parameter values are presented in Table 1.Algorithm α β γStatic 0.92 0.69 0.07Dynamic 0.96 1.00 0.03Table 1. Fitted parameter values for the load-balancing algorithms. See text for details.Effective Resource Allocation 765Fig. 2. Time of the SCF calculations for the model system as a function of the number ofCPUs. The time is relative to the 1 CPU case. See text for details.The non-zero value of the γ parameter stems from both non-parallelized partsof the code (mainly the Fock matrix diagonalization) and, to a lesser extent, fromthe communication overhead. The difference between the values for the static anddynamic dispatcher seems to be mainly a consequence of the relatively poor fit ofthe γ parameter. Still, in both cases the value is small, which shows that the bulkof the code is executed in the parallel mode.The main difference between the algorithms is the value of the β parameter,which decides how much faster works the parallelized part of the code with increasingnumber of the computational nodes. While the scaling for the static case seems tobe poor, for the dynamic case it has reached the theoretical limit of 1. However,the performance of the dynamic load-balancing algorithm can still be improved bytaking advantage of the hyper-cache effect.3.2 Integral CacheTo further improve the efficiency of the integral contraction, some of the quantum-chemical programs use integral caching. The cache mechanism relies on storingalready computed integrals in fast memory. As the same integrals are used in eachiteration of SCF cycle there are significant profits to be realized if their values arecached in the fast memory efficiently. However, for real systems the total number766 G. Mazur, M. Makowski, M. Brelaof integrals is too large to allow for keeping them all in the fast memory. Therefore,only a fraction of the integrals can be stored in the cache. The amount of cachedintegrals depends on the accessible memory.The caching strategy used in this paper is a very simple one. The integralsare stored in the order they were requested by the Fock matrix generation proce-dure until the assigned memory is exhausted. Other, prospectively more effectiveapproaches are possible, like giving preference to the integrals involving orbitals ofhigh angular momentum.Analogous technique can be applied to other computationally intensive proce-dures, like three-electron integrals used for the density-fitting and atomic orbitalvalues used for exchange-correlation potential generation. However, in this paperwe concentrate on the ERI cache.3.3 Cache AffinityIn the case of single processor, caching is conceptually simple but limited due tousually small memory to use for the purpose. n-node calculations offer n timesextended cache memory space and potentially considerable time savings are possibleif the load-balancing algorithm allows to use it in an efficient way. The efficiency isdetermined by cache affinity, which decides how large is the probability that if anintegral was processed by given node, it will be processed by the same node in thenext SCF cycle.From the point of view of the cache affinity, the static load-balancing algorithmis optimal. This is because the integrals in the consecutive SCF cycles are alwaysprocessed by the same nodes as in the previous cycles. On the other hand, the dy-namic load-balancing algorithm represents the case of minimal cache affinity, withthe probability of an integral being processed by the same node being equal to 1/n.However, despite the low cache affinity, the dynamic algorithm performs much bet-ter than the static one. Therefore, we introduce cache-awareness to the dynamicalgorithm to improve its cache affinity without loosing its adaptivity.The cache-aware load-balancing algorithm works as follows. The computationalproblem is divided into tasks of different sizes. The tasks are bundled into groupsof size n where n stands for the number of nodes. The task groups are stored in thetask groups queue. A node gets the next task from the current task group as soonas it completes the previous one. The task to be handed out to the node i is takenfrom the i-th slot if available or chosen randomly otherwise. When no more tasksare available in the current group, the group is retired and the next one is takenfrom the queue.The splitting of the problem into tasks is organized as follows. A fraction 1/fof the original problem is divided into n tasks, where f is the splitting factor. Thetasks are grouped together, and the group is put to the queue. Then the procedureis repeated recursively for the remaining part of the problem. The recurrence isstopped when the size of the remaining part is smaller than the threshold t. Thealgorithm is parametrized by the values of f and t.Effective Resource Allocation 767The algorithm was implemented in Niedoida. Test calculations were performedfor the same model system as in Section 3.1.Fig. 3. Cache hit ratio as a function of the number of CPUs. Sexithiophene molecule,HF/6-31G** level of theory, 256 MB of cache memory per CPU.The cache affinity of a load-balancing algorithm can be characterized by thecache hit ratio. Cache hit ratio is the fraction of the integrals which were retrievedfrom the cache to the total number of integrals required to perform the whole SCFprocess. The dependence of the cache hit ratio on the number of computationalnodes for the cache-aware algorithm, while worse than for the static one, is signifi-cantly better than for the dynamic scheme (see Figure 3).Comparison of the speed of the analyzed algorithms shows that for small cachesizes the cache-aware algorithms performs on par with the dynamic one. With thegrowing cache size, the dynamic algorithm performs slightly worse and the cache-aware algorithm performs better. For large cache size (larger than 256 MB per node)the cache-aware algorithm is the fastest one.To analyze the scaling of the cache-aware algorithm in more detail we performedcalculations for the model system using large cache size (512 MB per node). Therelative time of the calculations is depicted in Figure 4. Fitting the parameters ofEquation (4) yields α = 1.00, β = 1.07 and γ = 0.02. Comparison with Table 1shows that the cache-aware algorithm scales better than the other analyzed algo-rithms. The value of β exceeding the theoretical limit of 1 suggests presence of thehyper-cache effect.768 G. Mazur, M. Makowski, M. BrelaFig. 4. Time of the SCF calculations for the model system as a function of the number ofCPUs. The time is relative to the 1 CPU case. See text for details.4 PARALLEL MP2 CALCULATIONSNaively, the time complexity of the transformation of Equation (3) reaches O(N 8).It can be reduced to O(N 5) using the four-stage algorithmfor all i, a, j, b(iν|κλ) =∑µCµi(µν|κλ)(ia|κλ) =∑νCνa(iν|κλ)(ia|jλ) =∑κCκj(ia|κλ)(ia|jb) =∑λCλb(ia|jλ).However, the reduction of the time complexity comes at the cost of memory com-plexity reaching O(N 4) as the tensors of partially-transformed integrals have to bekept in memory.Effective Resource Allocation 769Given the memory available in currently used computer architectures and thesize of typical calculations, the four-stage transformation has to be split into severalpasses Pfor all Pfor i ∈ Pfor all a, j, b(iν|κλ) =∑µCµi(µν|κλ)(ia|κλ) =∑νCνa(iν|κλ)(ia|jλ) =∑κCκj(ia|κλ)(ia|jb) =∑λCλb(ia|jλ).This way the memory consumption is reduced to O(N 3). Such division of thewhole transformation into passes maps directly to the natural partitioning of theproblem into parallel tasks.The characteristic feature of the partitioning is that large part of the two-electronintegrals tensor is calculated in every pass. Therefore, contrary to the Hartree-Fockcase, CPU-time for integral generation averages out effectively, leaving the time ofthe single pass to depend solely on the number of orbitals treated. Hence, assuminghomogeneous execution environment, the static load-balancing algorithm should beoptimal.Test calculations were performed for linear chain of 15 water molecules in theSTO-3G basis. The scaling of the proposed load-balancing algorithms with thenumber of CPUs is shown in Figure 5.Fitting the parameters of Equation (4) yields α = 0.94, β = 1.06 and γ = 0.05.The modest super-linear scaling is attributed to the hyper-cache effect.5 SUMMARYAnalysis of the key aspects affecting parallel performance of the core computationalquantum chemistry algorithms was performed. It has been shown that despite thesuperficial similarity between the considered algorithms, their parallel form requiresdifferent load balancing schemes to yield an efficient solution. The analysis allowedus to design close to optimal load balancing procedures for the studied families ofquantum-chemical calculations.It has been shown that explicit introduction of the memory affinity to the loadbalancing procedures allows for superlinear scaling. We were able to exploit thehypercache effect in all computational quantum-chemistry algorithms considered inthis work.770 G. Mazur, M. Makowski, M. BrelaFig. 5. Time of the MP2 calculations for the model system as a function of the number ofCPUs. The time is relative to the 1 CPU case. See text for details.AcknowledgmentsSupport from ACK Cyfronet (grant MEiN/SGI3700/UJ/077/2006) is gratefully ac-knowledged.REFERENCES[1] Mazur, G.—Makowski, M.—Piskorz, W.—Ćwiklik,  L.—Sterzel, M.—Radoń, M.—Kulig, W.—Jagoda-Ćwiklik, B.—B lażewicz, D.: Niedoida 0.3.2007.[2] MPI-2: Extensions to the message-passing interface. http://www-unix.mcs.anl.gov/mpi/mpi-standard/mpi-report-2.0/mpi2-report.htm.[3] Roothaan, C. C. J.: New Developments in Molecular Orbital Theory. Rev. Mod.Phys., Vol. 23, 1951, p. 69.[4] Helgaker, T.—Jorgensen P.—Olsen, J.: Molecular Electronic Structure The-ory. John Wiley and Sons, Ltd., 2000.[5] Becke, A. D.: Density-Functional Exchange-Energy Approximation With CorrectAsymptotic Behaviour. Phys. Rev. A, Vol. 38, 1988, p. 3098.Effective Resource Allocation 771[6] Mazur, G.—Makowski, M.: Development and Optimization of ComputationalChemistry Algorithms. Computing and Informatics, Vol. 28, 2009, p. 115.Grzegorz Mazur is an Assistant Professor at the Departmentof Computational Methods in Chemistry of Jagiellonian Univer-sity in Cracow (Poland). He attained his Ph.D. in chemistry in2001 at the same university. His current research interests in-clude theoretical description of the excited states properties andoptimization of quantum chemistry algorithms.Marcin Makowski is an Assistant Professor at the Depart-ment of Theoretical Chemistry of Jagiellonian University in Cra-cow (Poland). He attained his M. Sc. degree in 2001 and hisPh.D. in 2004, both in chemistry, at the same university, andhis B. Sc. degree in computer science at University of Mining andMetallurgy in Cracow in 2004. His current research interests in-clude theoretical molecular spectroscopy, linear scaling methodsin quantum chemistry and optimization of quantum chemistryalgorithms. He participates in several research projects relatedwith the development of theoretical chemistry formalisms andnumerical methodologies that allow to efficiently calculate electronic structure of largemolecular systems.Mateusz Brela is working towards his M. Sc. Degree in theo-retical molecular spectroscopy from the Faculty of chemistry ofJagiellonian University in Cracow (Poland). His current researchinterests include optimization of quantum chemistry algorithms,linear scaling methods in quantum chemistry, theoretical mole-cular spectroscopy and molecular dynamics.

Effective Resource Allocation in Parallel Quantum-Chemical Calculations

http://www.cai.sk/ojs/index.php/cai/article/download/195/165

Effective Resource Allocation in Parallel Quantum-Chemical Calculations

Abstract

Similar works

Full text

Available Versions

Computing and Informatics (E-Journal - Institute of Informatics, SAS, Bratislava)