Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-2023 The plumed team
3 : (see the PEOPLE file at the root of the distribution for a list of names)
4 :
5 : See http://www.plumed.org for more information.
6 :
7 : This file is part of plumed, version 2.
8 :
9 : plumed is free software: you can redistribute it and/or modify
10 : it under the terms of the GNU Lesser General Public License as published by
11 : the Free Software Foundation, either version 3 of the License, or
12 : (at your option) any later version.
13 :
14 : plumed is distributed in the hope that it will be useful,
15 : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 : GNU Lesser General Public License for more details.
18 :
19 : You should have received a copy of the GNU Lesser General Public License
20 : along with plumed. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 : #ifndef __PLUMED_tools_Tools_h
23 : #define __PLUMED_tools_Tools_h
24 :
25 : #include "AtomNumber.h"
26 : #include "Vector.h"
27 : #include "Tensor.h"
28 : #include "small_vector/small_vector.h"
29 : #include <vector>
30 : #include <string>
31 : #include <cctype>
32 : #include <cstdio>
33 : #include <cmath>
34 : #include <limits>
35 : #include <algorithm>
36 : #include <sstream>
37 : #include <memory>
38 : #include <cstddef>
39 : #include <queue>
40 : #include <mutex>
41 : #include <filesystem>
42 : #include <utility>
43 : #include <unordered_map>
44 : #include <map>
45 : #include <condition_variable>
46 :
47 : namespace PLMD {
48 :
49 : class IFile;
50 :
51 : /// \ingroup TOOLBOX
52 : /// Very small non-zero number
53 : constexpr double epsilon(std::numeric_limits<double>::epsilon());
54 :
55 : /// \ingroup TOOLBOX
56 : /// Boltzman constant in kj/K
57 : constexpr double kBoltzmann(0.0083144621);
58 :
59 : /// \ingroup TOOLBOX
60 : /// PI
61 : constexpr double pi(3.141592653589793238462643383279502884197169399375105820974944592307);
62 :
63 : constexpr double dp2cutoff(6.25);
64 :
65 : constexpr double dp2cutoffA=1.00193418799744762399; // 1.0/(1-std::exp(-dp2cutoff));
66 : constexpr double dp2cutoffB=-.00193418799744762399; // -std::exp(-dp2cutoff)/(1-std::exp(-dp2cutoff));
67 :
68 5719 : inline static bool dp2cutoffNoStretch() {
69 5719 : static const auto* res=std::getenv("PLUMED_DP2CUTOFF_NOSTRETCH");
70 5719 : return res;
71 : }
72 :
73 : /// \ingroup TOOLBOX
74 : /// Empty class which just contains several (static) tools
75 : class Tools {
76 : /// class to convert a string to a generic type T
77 : template<class T>
78 : static bool convertToAny(const std::string & str,T &t);
79 : /// class to convert a string to a real type T.
80 : /// T should be either float, double, or long double
81 : template<class T>
82 : static bool convertToReal(const std::string & str,T &t);
83 : /// class to convert a string to a int type T
84 : template<class T>
85 : static bool convertToInt(const std::string & str,T &t);
86 : /// @brief the recursive part of the template fastpow implementation
87 : template <int exp, typename T=double, std::enable_if_t< (exp >=0), bool> = true>
88 : static inline /*consteval*/ T fastpow_rec(T base, T result);
89 : public:
90 : /// Split the line in words using separators.
91 : /// It also take into account parenthesis. Outer parenthesis found are removed from
92 : /// output, and the text between them is considered as a single word. Only the
93 : /// outer parenthesis are processed, to allow nesting them.
94 : /// parlevel, if not NULL, is increased or decreased according to the number of opened/closed parenthesis
95 : static std::vector<std::string> getWords(std::string_view line,const char* sep=NULL,int* parlevel=NULL,const char* parenthesis="{", const bool& delete_parenthesis=true);
96 : /// Faster version
97 : /// This version does not parse parenthesis and operates on a preallocated small_vector of string_view's
98 : static void getWordsSimple(gch::small_vector<std::string_view> & words,std::string_view line);
99 : /// Get a line from the file pointer ifile
100 : static bool getline(FILE*,std::string & line);
101 : /// Get a parsed line from the file pointer ifile
102 : /// This function already takes care of joining continued lines and splitting the
103 : /// resulting line into an array of words
104 : static bool getParsedLine(IFile&ifile,std::vector<std::string> & line, const bool trimcomments=true);
105 : /// compare two string in a case insensitive manner
106 : static bool caseInSensStringCompare(const std::string & str1, const std::string &str2);
107 : /// Convert a string to a double, reading it
108 : static bool convertNoexcept(const std::string & str,double & t);
109 : /// Convert a string to a long double, reading it
110 : static bool convertNoexcept(const std::string & str,long double & t);
111 : /// Convert a string to a float, reading it
112 : static bool convertNoexcept(const std::string & str,float & t);
113 : /// Convert a string to a int, reading it
114 : static bool convertNoexcept(const std::string & str,int & t);
115 : /// Convert a string to a long int, reading it
116 : static bool convertNoexcept(const std::string & str,long int & t);
117 : /// Convert a string to a long long int, reading it
118 : static bool convertNoexcept(const std::string & str,long long int & t);
119 : /// Convert a string to an unsigned int, reading it
120 : static bool convertNoexcept(const std::string & str,unsigned & t);
121 : /// Convert a string to a long unsigned int, reading it
122 : static bool convertNoexcept(const std::string & str,long unsigned & t);
123 : /// Convert a string to a long long unsigned int, reading it
124 : static bool convertNoexcept(const std::string & str,long long unsigned & t);
125 : /// Convert a string to a atom number, reading it
126 : static bool convertNoexcept(const std::string & str,AtomNumber & t);
127 : /// Convert a string to a string (i.e. copy)
128 : static bool convertNoexcept(const std::string & str,std::string & t);
129 : /// Convert anything into a string
130 : template<typename T>
131 : static bool convertNoexcept(T i,std::string & str);
132 : /// Convert anything into anything, throwing an exception in case there is an error
133 : /// Remove trailing blanks
134 : static void trim(std::string & s);
135 : /// Remove leading blanks
136 : static void ltrim(std::string & s);
137 : /// Remove trailing comments
138 : static void trimComments(std::string & s);
139 : /// Apply pbc for a unitary cell
140 : static double pbc(double);
141 : /// Retrieve a key from a vector of options.
142 : /// It finds a key starting with "key=" or equal to "key" and copy the
143 : /// part after the = on s. E.g.:
144 : /// line.push_back("aa=xx");
145 : /// getKey(line,"aa",s);
146 : /// will set s="xx"
147 : static bool getKey(std::vector<std::string>& line,const std::string & key,std::string & s,int rep=-1);
148 : /// Find a keyword on the input line, eventually deleting it, and saving its value to val
149 : template <typename T,typename U>
150 24810670 : static void convert(const T & t,U & u) {
151 24811734 : plumed_assert(convertNoexcept(t,u)) <<"Error converting "<<t;
152 24810669 : }
153 : template <typename T>
154 : static bool parse(std::vector<std::string>&line,const std::string&key,T&val,int rep=-1);
155 : /// Find a keyword on the input line, eventually deleting it, and saving its value to a vector
156 : template <class T>
157 : static bool parseVector(std::vector<std::string>&line,const std::string&key,std::vector<T>&val,int rep=-1);
158 : /// Find a keyword without arguments on the input line
159 : static bool parseFlag(std::vector<std::string>&line,const std::string&key,bool&val);
160 : /// Find a keyword on the input line, just reporting if it exists or not
161 : static bool findKeyword(const std::vector<std::string>&line,const std::string&key);
162 : /// Interpret atom ranges
163 : static void interpretRanges(std::vector<std::string>&);
164 : /// Remove duplicates from a vector of type T
165 : template <typename T>
166 : static void removeDuplicates(std::vector<T>& vec);
167 : /// interpret ":" syntax for labels
168 : static void interpretLabel(std::vector<std::string>&s);
169 : /// list files in a directory
170 : static std::vector<std::string> ls(const std::string&);
171 : /// removes leading and trailing blanks from a string
172 : static void stripLeadingAndTrailingBlanks( std::string& str );
173 : /// Extract the extensions from a file name.
174 : /// E.g.: extension("pippo.xyz")="xyz".
175 : /// It only returns extensions with a length between 1 and 4
176 : /// E.g.: extension("pippo.12345")="" whereas extenion("pippo.1234")="1234";
177 : /// It is also smart enough to detect "/", so that
178 : /// extension("pippo/.t")="" whereas extension("pippo/a.t")="t"
179 : static std::string extension(const std::string&);
180 : /// Fast int power
181 : static double fastpow(double base,int exp);
182 : /// Fast int power for power known at compile time
183 : template <int exp, typename T=double>
184 : static inline /*consteval*/ T fastpow(T base);
185 : /// Modified 0th-order Bessel function of the first kind
186 : static double bessel0(const double& val);
187 : /// Check if a string full starts with string start.
188 : /// Same as full.find(start)==0
189 : static bool startWith(const std::string & full,const std::string &start);
190 : /**
191 : Tool to create a vector of raw pointers from a vector of unique_pointers (const version).
192 : Returning a vector is fast in C++11. It can be used in order to feed a vector<unique_ptr<T>>
193 : to a function that takes a vector<T*>.
194 : \verbatim
195 : // some function that takes a vec
196 : void func(std::vector<Data*> & vec);
197 : std::vector<std::unique_ptr<Data>> vec;
198 : // func(vec); // does not compile
199 : func(Tools::unique2raw(vec)); // compiles
200 : \endverbatim
201 : Notice that the conversion is fast but takes
202 : some time to allocate the new vector and copy the pointers. In case the function
203 : acting on the vector<T*> is very fast and we do not want to add significant overhead,
204 : it might be convenient to store a separate set of raw pointers.
205 : \verbatim
206 : // some function that takes a vec
207 : void func(std::vector<Data*> & vec);
208 : std::vector<std::unique_ptr<Data>> vec;
209 :
210 : // conversion done only once:
211 : auto vec_ptr=Tools::unique2raw(vec);
212 :
213 : for(int i=0;i<1000;i++){
214 : func(vec_ptr);
215 : }
216 : \endverbatim
217 : */
218 : template <typename T>
219 : static std::vector<T*> unique2raw(const std::vector<std::unique_ptr<T>>&);
220 : /// Tool to create a vector of raw pointers from a vector of unique_pointers.
221 : /// See the non const version.
222 : template <typename T>
223 : static std::vector<const T*> unique2raw(const std::vector<std::unique_ptr<const T>>&);
224 : /// Tiny class that changes directory and comes back when going out of scope.
225 : /// In case system calls to change dir are not available it throws an exception.
226 : /// \warning By construction, changing directory breaks thread safety! Use with care.
227 : class DirectoryChanger {
228 : const std::filesystem::path path;
229 : public:
230 : explicit DirectoryChanger(const char*path);
231 : ~DirectoryChanger();
232 : };
233 :
234 : template<class T, class... Args>
235 : static auto make_unique(Args&&... args) {
236 154288 : return std::make_unique<T>(std::forward<Args>(args)...);
237 : }
238 :
239 : static void set_to_zero(double*ptr,unsigned n) {
240 1276640 : for(unsigned i=0; i<n; i++) {
241 1247901 : ptr[i]=0.0;
242 : }
243 : }
244 :
245 : template<unsigned n>
246 28739 : static void set_to_zero(std::vector<VectorGeneric<n>> & vec) {
247 28739 : unsigned s=vec.size();
248 28739 : if(s==0) {
249 : return;
250 : }
251 28739 : set_to_zero(&vec[0][0],s*n);
252 : }
253 :
254 : template<unsigned n,unsigned m>
255 : static void set_to_zero(std::vector<TensorGeneric<n,m>> & vec) {
256 : unsigned s=vec.size();
257 : if(s==0) {
258 : return;
259 : }
260 : set_to_zero(&vec[0](0,0),s*n*m);
261 : }
262 :
263 : static std::unique_ptr<std::lock_guard<std::mutex>> molfile_lock();
264 : /// Build a concatenated exception message.
265 : /// Should be called with an in-flight exception.
266 : static std::string concatenateExceptionMessages();
267 :
268 :
269 : /// Tiny class implementing faster std::string_view access to an unordered_map
270 : /// It exposes a limited number of methods of std::unordered_map. Others could be added.
271 : /// Importantly, when it is accessed via a std::string_view, the access does not
272 : /// require constructing a std::string and is thus faster.
273 : /// Deletion would be slower instead. It's not even implemented yet.
274 : template<class T>
275 806844 : class FastStringUnorderedMap {
276 : std::unordered_map<std::string_view,T> map;
277 : std::vector<std::unique_ptr<const char[]>> keys;
278 :
279 : // see https://stackoverflow.com/questions/34596768/stdunordered-mapfind-using-a-type-different-than-the-key-type
280 507956 : std::unique_ptr<const char[]> conv(std::string_view str) {
281 507956 : auto p=std::make_unique<char[]>(str.size()+1);
282 : std::memcpy(p.get(), str.data(), str.size()+1);
283 507956 : return p;
284 : }
285 :
286 : public:
287 :
288 : FastStringUnorderedMap() = default;
289 16385 : FastStringUnorderedMap(std::initializer_list<std::pair<const std::string_view,T>> init) {
290 510335 : for(const auto & c : init) {
291 493950 : (*this)[c.first]=c.second;
292 : }
293 16385 : }
294 :
295 3557350 : T& operator[]( const std::string_view & key ) {
296 : auto f=map.find(key);
297 3557350 : if(f!=map.end()) {
298 3049394 : return f->second;
299 : }
300 1015912 : keys.push_back(conv(key));
301 507956 : return map[keys.back().get()];
302 : }
303 :
304 : auto begin() {
305 : return map.begin();
306 : }
307 : auto end() {
308 : return map.end();
309 : }
310 : auto begin() const {
311 : return map.begin();
312 : }
313 : auto end() const {
314 : return map.end();
315 : }
316 : auto find(const std::string_view & key) {
317 : return map.find(key);
318 : }
319 : auto find(const std::string_view & key) const {
320 : return map.find(key);
321 : }
322 : };
323 :
324 : /// Utility to create named critical sections
325 : /// Key should be usable in a std::map
326 : template<class Key>
327 : class CriticalSectionWithKey {
328 : std::mutex mutex;
329 : std::condition_variable notify;
330 : std::map<Key, int> in_progress;
331 : public:
332 23 : void start(const Key & key) {
333 23 : std::unique_lock<std::mutex> lock(mutex);
334 128 : while (in_progress[key] > 0) {
335 : // Wait if this command is already in progress.
336 105 : notify.wait(lock);
337 : }
338 : // Mark this command as in progress.
339 23 : in_progress[key]++;
340 23 : }
341 23 : void stop(const Key & key) {
342 23 : std::unique_lock<std::mutex> lock(mutex);
343 : // Mark this command as completed.
344 23 : in_progress[key]--;
345 : // Notify other threads that may be waiting for this command to complete.
346 23 : notify.notify_all();
347 23 : }
348 : class Handler {
349 : CriticalSectionWithKey* section{nullptr};
350 : Key key;
351 23 : Handler(CriticalSectionWithKey* section,const Key& key):
352 23 : section(section),
353 23 : key(key) {
354 23 : section->start(key);
355 23 : }
356 : friend class CriticalSectionWithKey;
357 : public:
358 : /// Default constructor
359 : Handler() = default;
360 : /// Default copy constructor is deleted (not copyable)
361 : Handler(const Handler & handler) = delete;
362 : /// Default copy assignment is deleted (not copyable)
363 : Handler & operator=(const Handler & handler) = delete;
364 : /// Move constructor.
365 : Handler(Handler && handler) noexcept :
366 : section(handler.section),
367 : key(std::move(handler.key)) {
368 : handler.section=nullptr;
369 : };
370 : /// Move assignment.
371 : Handler & operator=(Handler && handler) noexcept {
372 : if(this!=&handler) {
373 : if(section) {
374 : section->stop(key);
375 : }
376 : section=handler.section;
377 : key=std::move(handler.key);
378 : }
379 : handler.section=nullptr;
380 : return *this;
381 : }
382 : /// Destructor
383 23 : ~Handler() {
384 23 : if(section) {
385 23 : section->stop(key);
386 : }
387 23 : }
388 : };
389 :
390 : Handler startStop(const Key & key) {
391 23 : return Handler(this,key);
392 : }
393 :
394 : };
395 :
396 : };
397 :
398 : template <class T>
399 386534 : bool Tools::parse(std::vector<std::string>&line,const std::string&key,T&val,int rep) {
400 : std::string s;
401 773068 : if(!getKey(line,key+"=",s,rep)) {
402 : return false;
403 : }
404 311779 : if(s.length()>0 && !convertNoexcept(s,val)) {
405 3 : return false;
406 : }
407 : return true;
408 : }
409 :
410 : template <class T>
411 283702 : bool Tools::parseVector(std::vector<std::string>&line,const std::string&key,std::vector<T>&val,int rep) {
412 : std::string s;
413 567404 : if(!getKey(line,key+"=",s,rep)) {
414 : return false;
415 : }
416 : val.clear();
417 261605 : std::vector<std::string> words=getWords(s,"\t\n ,");
418 5525012 : for(unsigned i=0; i<words.size(); ++i) {
419 : T v;
420 5263407 : std::string s=words[i];
421 5263407 : const std::string multi("@replicas:");
422 5263407 : if(rep>=0 && startWith(s,multi)) {
423 6 : s=s.substr(multi.length(),s.length());
424 6 : std::vector<std::string> words=getWords(s,"\t\n ,");
425 6 : plumed_assert(rep<static_cast<int>(words.size()));
426 6 : s=words[rep];
427 6 : }
428 5263407 : if(!convertNoexcept(s,v)) {
429 : return false;
430 : }
431 5263407 : val.push_back(v);
432 : }
433 : return true;
434 261605 : }
435 :
436 : template<typename T>
437 91907 : void Tools::removeDuplicates(std::vector<T>& vec) {
438 91907 : std::sort(vec.begin(), vec.end());
439 91907 : vec.erase(std::unique(vec.begin(), vec.end()), vec.end());
440 91907 : }
441 :
442 : inline
443 98713 : bool Tools::parseFlag(std::vector<std::string>&line,const std::string&key,bool&val) {
444 791283 : for(auto p=line.begin(); p!=line.end(); ++p) {
445 704755 : if(key==*p) {
446 12185 : val=true;
447 : line.erase(p);
448 : return true;
449 : }
450 : }
451 : return false;
452 : }
453 :
454 : /// beware: this brings any number into a pbc that ranges from -0.5 to 0.5
455 : inline
456 1786274180 : double Tools::pbc(double x) {
457 : #ifdef __PLUMED_PBC_WHILE
458 : while (x>0.5) {
459 : x-=1.0;
460 : }
461 : while (x<-0.5) {
462 : x+=1.0;
463 : }
464 : return x;
465 : #else
466 : if constexpr (std::numeric_limits<int>::round_style == std::round_toward_zero) {
467 : constexpr double offset=100.0;
468 1786274180 : const double y=x+offset;
469 1786274180 : if(y>=0) {
470 1786269602 : return y-int(y+0.5);
471 : } else {
472 4578 : return y-int(y-0.5);
473 : }
474 : } else if constexpr (std::numeric_limits<int>::round_style == std::round_to_nearest) {
475 : return x-int(x);
476 : } else {
477 : return x-floor(x+0.5);
478 : }
479 : #endif
480 : }
481 :
482 : template<typename T>
483 8081659 : bool Tools::convertNoexcept(T i,std::string & str) {
484 8081659 : std::ostringstream ostr;
485 810734 : ostr<<i;
486 8081659 : str=ostr.str();
487 8081659 : return true;
488 8081659 : }
489 :
490 : inline
491 56283115 : double Tools::fastpow(double base, int exp) {
492 56283115 : if(exp<0) {
493 0 : exp=-exp;
494 0 : base=1.0/base;
495 : }
496 : double result = 1.0;
497 289151934 : while (exp) {
498 199268627 : if (exp & 1) {
499 140987083 : result *= base;
500 : }
501 199268627 : exp >>= 1;
502 199268627 : base *= base;
503 : }
504 :
505 56283115 : return result;
506 : }
507 :
508 : template <int exp, typename T, std::enable_if_t< (exp >=0), bool>>
509 : inline T Tools::fastpow_rec(T const base, T result) {
510 : if constexpr (exp == 0) {
511 : return result;
512 : }
513 : if constexpr (exp & 1) {
514 19283644 : result *= base;
515 : }
516 11328670 : return fastpow_rec<(exp>>1),T> (base*base, result);
517 : }
518 :
519 : template <int exp, typename T>
520 : inline T Tools::fastpow(T const base) {
521 : if constexpr (exp<0) {
522 : return fastpow_rec<-exp,T>(1.0/base,1.0);
523 : } else {
524 : return fastpow_rec<exp,T>(base, 1.0);
525 : }
526 : }
527 :
528 : template<typename T>
529 2307 : std::vector<T*> Tools::unique2raw(const std::vector<std::unique_ptr<T>> & x) {
530 2307 : std::vector<T*> v(x.size());
531 3681 : for(unsigned i=0; i<x.size(); i++) {
532 1374 : v[i]=x[i].get();
533 : }
534 2307 : return v;
535 : }
536 :
537 : template<typename T>
538 : std::vector<const T*> Tools::unique2raw(const std::vector<std::unique_ptr<const T>> & x) {
539 : std::vector<const T*> v(x.size());
540 : for(unsigned i=0; i<x.size(); i++) {
541 : v[i]=x[i].get();
542 : }
543 : return v;
544 : }
545 :
546 : }
547 :
548 : #endif
|
