openCARP
Doxygen code documentation for the open cardiac electrophysiology simulator openCARP
basics.h
Go to the documentation of this file.
1 // ----------------------------------------------------------------------------
2 // openCARP is an open cardiac electrophysiology simulator.
3 //
4 // Copyright (C) 2020 openCARP project
5 //
6 // This program is licensed under the openCARP Academic Public License (APL)
7 // v1.0: You can use and redistribute it and/or modify it in non-commercial
8 // academic environments under the terms of APL as published by the openCARP
9 // project v1.0, or (at your option) any later version. Commercial use requires
10 // a commercial license (info@opencarp.org).
11 //
12 // This program is distributed without any warranty; see the openCARP APL for
13 // more details.
14 //
15 // You should have received a copy of the openCARP APL along with this program
16 // and can find it online: http://www.opencarp.org/license
17 // ----------------------------------------------------------------------------
18 
27 #ifndef _BASICS_H
28 #define _BASICS_H
29 
30 #include <cstdio>
31 #include <cstdlib>
32 #include <cstdarg>
33 #include <cstring>
34 #include <cassert>
35 #include <string>
36 
37 #include <unistd.h>
38 #include <sys/stat.h>
39 #include <mpi.h>
40 #include <petscsys.h> // TODO: review this (needed for PETSC_COMM_WORLD)
41 
42 #include "vect.h"
43 #include "SF_vector.h"
44 
45 #ifndef CARP_PARAMS
46 #define CARP_PARAMS
47 #include "openCARP_p.h"
48 #include "openCARP_d.h"
49 #endif
50 
51 #define STRUCT_ZERO(S) memset(&(S), 0, sizeof(S))
52 #define SLIST_APPEND( S, P ) sltlst_append( S, &(P), sizeof(P) )
53 
54 namespace opencarp {
55 
57 struct Salt_list {
58  void *data;
59  int chunk;
60  int nitems;
61  int size;
62 };
63 
73 void sltlst_append(Salt_list* sl, void *p, int quantum);
74 
81 char* dupstr(const char* old_str);
82 
89 char *stringify(double r);
90 
91 std::string get_basename(const std::string & path);
92 
102 template<class STRVEC>
103 void split_string(const std::string & input, const char s, STRVEC & list)
104 {
105  size_t fnd = 0, ofnd = 0;
106  int num_parts = 1, idx = 0;
107 
108  fnd = input.find(s);
109  while (fnd != std::string::npos)
110  {
111  num_parts++;
112  ofnd = fnd+1;
113  fnd = input.find(s, ofnd);
114  }
115 
116  list.resize(num_parts);
117 
118  fnd = input.find(s); ofnd = 0;
119  while (fnd != std::string::npos)
120  {
121  list[idx++].assign(input.begin() + ofnd, input.begin() + fnd);
122  ofnd = fnd+1;
123  fnd = input.find(s, ofnd);
124  }
125 
126  if(ofnd < input.size())
127  list[idx].assign(input.begin() + ofnd, input.end());
128 }
129 
133 struct file_desc {
134  FILE* fd;
135  const char* name;
136 };
137 
139 
147 bool f_exist(const char *fname);
148 
157 FILE_SPEC f_open(const char *fname, const char *mode);
158 
164 void f_close(FILE_SPEC & f);
165 
175 void f_read_par(void *ptr, size_t size, size_t nmemb, FILE_SPEC stream, MPI_Comm comm = PETSC_COMM_WORLD);
176 
187 void f_write_par(void* ptr, size_t size, size_t nmemb, int source_pid, FILE_SPEC stream,
188  MPI_Comm comm = PETSC_COMM_WORLD);
189 
198 char *read_bin_string(FILE_SPEC in);
199 
210 char* read_bin_string_par(FILE_SPEC in);
211 
212 char* f_gets_par(char* s, int size, FILE_SPEC stream, MPI_Comm comm = PETSC_COMM_WORLD);
213 
219 void write_bin_string(FILE_SPEC out, const char *s);
220 
232 template<typename T> inline
233 T get_global(T in, MPI_Op OP, MPI_Comm comm = PETSC_COMM_WORLD)
234 {
235  double dat = in;
236  MPI_Allreduce(MPI_IN_PLACE, &dat, 1, MPI_DOUBLE, OP, comm);
237  return T(dat);
238 }
239 
240 template<> inline
241 int get_global<int>(int in, MPI_Op OP, MPI_Comm comm)
242 {
243  int dat = in;
244  MPI_Allreduce(MPI_IN_PLACE, &dat, 1, MPI_INT, OP, comm);
245  return dat;
246 }
247 
248 template<> inline
249 long int get_global<long int>(long int in, MPI_Op OP, MPI_Comm comm)
250 {
251  long int dat = in;
252  MPI_Allreduce(MPI_IN_PLACE, &dat, 1, MPI_LONG, OP, comm);
253  return dat;
254 }
255 
263 template<typename T> inline
264 void get_global(SF::vector<T> & vec, int sender = 0, MPI_Comm comm = PETSC_COMM_WORLD)
265 {
266  size_t sz = vec.size();
267  MPI_Bcast(&sz, sizeof(size_t), MPI_BYTE, sender, comm);
268  vec.resize(sz);
269  MPI_Bcast(vec.data(), sz*sizeof(T), MPI_BYTE, sender, comm);
270 }
271 
276 inline int get_rank(MPI_Comm comm = PETSC_COMM_WORLD)
277 {
278  int rank;
279  MPI_Comm_rank(comm, &rank);
280 
281  return rank;
282 }
283 
290 inline int get_size(MPI_Comm comm = PETSC_COMM_WORLD)
291 {
292  int size;
293 
294  MPI_Comm_size(comm, &size);
295 
296  return size;
297 }
298 
299 inline int get_remote_size(MPI_Comm intercomm)
300 {
301  int size;
302 
303  MPI_Comm_remote_size(intercomm, &size);
304  return size;
305 }
306 
307 // flags for log messages
308 #define ECHO 1 // screen output
309 #define LOCAL 2 // output all cores
310 #define SYNCED 4 // synchronized
311 #define FLUSH 8 // flush screen output
312 #define NONL 16 // no new line
313 
314 #define MAX_MESG_LEN 2048
315 #define MAX_LOG_LEVEL 5
316 
340 void log_msg(FILE_SPEC out, int level, unsigned char flag, const char *fmt, ...);
341 
343 void init_iterations_logger(FILE_SPEC* & logger, const char* filename);
344 
345 inline void remove_preceding_char(char* buff, const int buffsize, const char c)
346 {
347  int ridx = 0, widx=0;
348  while(ridx < buffsize && buff[ridx] == c) ridx++;
349 
350  while(ridx < (buffsize - 1))
351  buff[widx++] = buff[ridx++];
352 
353  buff[widx] = '\0';
354 }
355 
356 inline void remove_char(char* buff, const int buffsize, const char c)
357 {
358  int ridx = 0, widx=0;
359  while(ridx < buffsize) {
360  if(buff[ridx] != c)
361  buff[widx++] = buff[ridx];
362  ridx++;
363  }
364 
365  buff[widx] = '\0';
366 }
367 
368 inline bool has_char(char* buff, const int buffsize, const char c)
369 {
370  int ridx = 0;
371  while(ridx < buffsize) {
372  if(buff[ridx] == '\0') return false;
373  if(buff[ridx] == c) return true;
374  ridx++;
375  }
376 
377  return false;
378 }
379 
381 class geom_shape {
382  public:
383  enum shape_t {sphere = 1, block = 2, cylinder = 3};
387  double radius;
388 };
389 
391 bool point_in_shape(const Point & p, const geom_shape & shape);
392 
398 bool is_big_endian();
399 
407 bool file_can_be_opened(const char* file);
408 
409 
411 bool path_is_absolute(const char* path);
412 
413 /* fmemopen for OSX / BSD compilations */
414 #ifdef __APPLE__
415 #define USE_FMEM_WRAPPER 1
416 #endif
417 
418 #ifdef OS_FREEBSD
419 #define USE_FMEM_WRAPPER 1
420 #endif
421 
422 #ifdef USE_FMEM_WRAPPER
423 struct fmem {
424  size_t pos;
425  size_t size;
426  char *buffer;
427 };
428 typedef struct fmem fmem_t;
429 
430 FILE *fmemopen_(void *, size_t, const char *);
431 #else
432 /* Else use the normal fmemopen */
433 #define fmemopen_ fmemopen
434 #endif
435 
436 inline void get_time(double & tm)
437 {
438  tm = MPI_Wtime();
439 }
440 
441 inline double get_time()
442 {
443  double tm = MPI_Wtime();
444  return tm;
445 }
446 
447 template<typename V> inline
448 V timing(V & t2, const V & t1)
449 {
450  t2 = get_time();
451  return t2 - t1;
452 }
453 
454 } // namespace opencarp
455 
456 // TODO: PETSc currently support neither float nor int64_T, so it's desactivated
457 // template _macro(std::int64_t, double)
458 #define OPENCARP_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(_macro) \
459  template _macro(std::int32_t, double)
460 
461 #endif
long int get_global< long int >(long int in, MPI_Op OP, MPI_Comm comm)
Definition: basics.h:249
int get_global< int >(int in, MPI_Op OP, MPI_Comm comm)
Definition: basics.h:241
The vector class and related algorithms.
void sltlst_append(Salt_list *sl, void *p, int quantum)
Definition: basics.cc:33
class to store shape definitions
Definition: basics.h:381
char * dupstr(const char *old_str)
Definition: basics.cc:44
void remove_preceding_char(char *buff, const int buffsize, const char c)
Definition: basics.h:345
const char * name
Definition: basics.h:135
bool is_big_endian()
Definition: basics.cc:290
void split_string(const std::string &input, const char s, STRVEC &list)
Split a string holding a character-seperated list into a vector of strings.
Definition: basics.h:103
void * data
the buffer
Definition: basics.h:58
bool has_char(char *buff, const int buffsize, const char c)
Definition: basics.h:368
void f_read_par(void *ptr, size_t size, size_t nmemb, FILE_SPEC stream, MPI_Comm comm)
Parallel fread. Root reads, then broadcasts.
Definition: basics.cc:171
3D vector algebra.
T * data()
Pointer to the vector&#39;s start.
Definition: SF_vector.h:91
void write_bin_string(FILE_SPEC out, const char *s)
Definition: basics.cc:219
void remove_char(char *buff, const int buffsize, const char c)
Definition: basics.h:356
char * read_bin_string_par(FILE_SPEC in)
Definition: basics.cc:239
char * read_bin_string(FILE_SPEC in)
Definition: basics.cc:228
int chunk
allocate memory to hold this many items
Definition: basics.h:59
bool file_can_be_opened(const char *file)
Check wheterh a file can be opened for reading.
Definition: basics.cc:301
void get_time(double &tm)
Definition: basics.h:436
T get_global(T in, MPI_Op OP, MPI_Comm comm=PETSC_COMM_WORLD)
Do a global reduction on a variable.
Definition: basics.h:233
V timing(V &t2, const V &t1)
Definition: basics.h:448
bool f_exist(const char *fname)
Definition: basics.cc:130
int get_remote_size(MPI_Comm intercomm)
Definition: basics.h:299
File descriptor struct.
Definition: basics.h:133
std::string get_basename(const std::string &path)
Definition: basics.cc:61
size_t size() const
The current size of the vector.
Definition: SF_vector.h:104
void log_msg(FILE_SPEC out, int level, unsigned char flag, const char *fmt,...)
Definition: basics.cc:72
int nitems
number of items
Definition: basics.h:60
A vector storing arbitrary data.
Definition: SF_vector.h:42
#define fmemopen_
Definition: basics.h:433
void f_close(FILE_SPEC &f)
Close a FILE_SPEC.
Definition: basics.cc:162
saltatory list – memory is allocated in chunks
Definition: basics.h:57
void init_iterations_logger(FILE_SPEC *&logger, const char *filename)
init a logger for solver iterations
bool point_in_shape(const Point &p, const geom_shape &shape)
test if a point is inside a simple geometric shape
Definition: basics.cc:250
file_desc * FILE_SPEC
Definition: basics.h:138
void resize(size_t n)
Resize a vector.
Definition: SF_vector.h:209
int get_size(MPI_Comm comm=PETSC_COMM_WORLD)
Definition: basics.h:290
char * f_gets_par(char *s, int size, FILE_SPEC stream, MPI_Comm comm)
Definition: basics.cc:199
int size
size of list items
Definition: basics.h:61
int get_rank(MPI_Comm comm=PETSC_COMM_WORLD)
Definition: basics.h:276
FILE_SPEC f_open(const char *fname, const char *mode)
Open a FILE_SPEC.
Definition: basics.cc:135
bool path_is_absolute(const char *path)
check whether path is absolute
Definition: basics.cc:319
void f_write_par(void *ptr, size_t size, size_t nmemb, int source_pid, FILE_SPEC stream, MPI_Comm comm)
Write in parallel. Data comes from one rank, rank 0 writes.
Definition: basics.cc:179
char * stringify(double r)
Definition: basics.cc:54