1 /* Ergo, version 3.8, a program for linear scaling electronic structure
2  * calculations.
3  * Copyright (C) 2019 Elias Rudberg, Emanuel H. Rubensson, Pawel Salek,
4  * and Anastasia Kruchinina.
5  *
6  * This program is free software: you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation, either version 3 of the License, or
9  * (at your option) any later version.
10  *
11  * This program is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with this program.  If not, see <http://www.gnu.org/licenses/>.
18  *
19  * Primary academic reference:
20  * Ergo: An open-source program for linear-scaling electronic structure
21  * calculations,
22  * Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, and Anastasia
23  * Kruchinina,
24  * SoftwareX 7, 107 (2018),
25  * <http://dx.doi.org/10.1016/j.softx.2018.03.005>
26  *
27  * For further information about Ergo, see <http://www.ergoscf.org>.
28  */
29 
30 /** @file integrals_1el_potential_prep.h
31 
32     @brief Code for 1-electron integrals, preparatory work for
33     computation of electron-nuclear potential energy matrix V.
34 
35     @author: Elias Rudberg <em>responsible</em>
36 */
37 
38 #ifndef INTEGRALS_1EL_POTENTIAL_PREP_HEADER
39 #define INTEGRALS_1EL_POTENTIAL_PREP_HEADER
40 
41 #include "basisinfo.h"
42 #include <algorithm>    // std::sort
43 
44 struct DistributionSpecStructWithIndexes2 {
45   DistributionSpecStruct distr;
46   int basisFuncIdx1;
47   int basisFuncIdx2;
48 };
49 
50 struct group_struct {
51   int startIndex;
52   int count;
53   int maxNoOfMoments;
54   int maxDegree;
55   ergo_real maxExtent;
56 };
57 
58 struct maxMomentVectorNormStruct {
59   ergo_real maxMomentVectorNormList[MAX_MULTIPOLE_DEGREE_BASIC+1];
60 };
61 
62 struct SetOfDistrsForVInfo {
63   ergo_real maxExtentForAll;
64   maxMomentVectorNormStruct maxMomentVectorNormForAll;
65   ergo_real boundingCubeCenterCoords[3];
66   ergo_real boundingCubeWidth;
67 };
68 
69 struct SetOfDistrsForV {
70   std::vector<DistributionSpecStructWithIndexes2> distrList;
71   std::vector<multipole_struct_small> multipoleList; // same size as distrList
72   std::vector<group_struct> groupList;
73   std::vector<maxMomentVectorNormStruct> maxMomentVectorNormList; // size same as groupList
74   SetOfDistrsForVInfo info;
75   // Stuff needed for Chunks&Tasks usage
76   SetOfDistrsForV();
77   SetOfDistrsForV(const SetOfDistrsForV & other);
78   void write_to_buffer ( char * dataBuffer, size_t const bufferSize ) const;
79   size_t get_size() const;
80   void assign_from_buffer ( char const * dataBuffer, size_t const bufferSize);
81 };
82 
83 void
84 organize_distrs_for_V(const IntegralInfo & integralInfo,
85 		      SetOfDistrsForV & setOfDistrsForV,
86 		      const std::vector<DistributionSpecStructWithIndexes2> & inputList,
87 		      ergo_real threshold,
88 		      ergo_real maxCharge);
89 
90 template <typename DistributionSpecStructType>
91 int
compare_distrs(const void * p1,const void * p2)92 compare_distrs(const void* p1, const void* p2) {
93   DistributionSpecStructType* d1 = (DistributionSpecStructType*)p1;
94   DistributionSpecStructType* d2 = (DistributionSpecStructType*)p2;
95   /* FIXME: Not nice to have these two hardcoded values here.  */
96   const ergo_real tolernance_dist = 1e-10;
97   const ergo_real tolernance_exponent = 1e-11;
98   ergo_real dx = d1->distr.centerCoords[0] - d2->distr.centerCoords[0];
99   if(dx > tolernance_dist)
100     return 1;
101   if(dx < -tolernance_dist)
102     return -1;
103   ergo_real dy = d1->distr.centerCoords[1] - d2->distr.centerCoords[1];
104   if(dy > tolernance_dist)
105     return 1;
106   if(dy < -tolernance_dist)
107     return -1;
108   ergo_real dz = d1->distr.centerCoords[2] - d2->distr.centerCoords[2];
109   if(dz > tolernance_dist)
110     return 1;
111   if(dz < -tolernance_dist)
112     return -1;
113   ergo_real de = d1->distr.exponent - d2->distr.exponent;
114   if(de > tolernance_exponent)
115     return 1;
116   if(de < -tolernance_exponent)
117     return -1;
118   return 0;
119 }
120 
121 template <typename DistributionSpecStructType>
122 bool
compare_distrs_bool(const DistributionSpecStructType & p1,const DistributionSpecStructType & p2)123 compare_distrs_bool(const DistributionSpecStructType & p1, const DistributionSpecStructType & p2) {
124   int i = compare_distrs<DistributionSpecStructType>(&p1, &p2);
125   return (i == 1);
126 }
127 
128 template <typename DistributionSpecStructType>
129 int
sort_distr_list(DistributionSpecStructType * list,int n)130 sort_distr_list(DistributionSpecStructType* list, int n) {
131   std::sort(&list[0], &list[n], compare_distrs_bool<DistributionSpecStructType>);
132   return 0;
133 }
134 
135 #endif
136