27 #ifndef _SF_CONTAINER_H 28 #define _SF_CONTAINER_H 74 return {arr[0], arr[1], arr[2]};
80 arr[0] = p.
x, arr[1] = p.
y, arr[2] = p.
z;
92 return a.
x*b.
x + a.
y*b.
y + a.
z*b.
z;
95 inline void outer_prod(
const Point & a,
const Point & b,
const double s,
double* buff,
const bool add =
false)
98 memset(buff, 0, 9*
sizeof(
double));
100 buff[0] += a.
x*b.
x*s, buff[1] += a.
x*b.
y*s, buff[2] += a.
x*b.
z*s;
101 buff[3] += a.
y*b.
x*s, buff[4] += a.
y*b.
y*s, buff[5] += a.
y*b.
z*s;
102 buff[6] += a.
z*b.
x*s, buff[7] += a.
z*b.
y*s, buff[8] += a.
z*b.
z*s;
118 double norm =
mag(vect);
119 return {vect.
x / norm, vect.
y / norm, vect.
z / norm};
124 return {a.
x - b.
x, a.
y - b.
y, a.
z - b.
z};
129 return {a.
x + b.
x, a.
y + b.
y, a.
z + b.
z};
134 return {a.
x * s, a.
y * s, a.
z * s};
159 if ( !strcmp( eletype,
"Tt" ) ) {
161 }
else if ( !strcmp( eletype,
"Hx" ) ) {
163 }
else if ( !strcmp( eletype,
"Oc" ) ) {
165 }
else if ( !strcmp( eletype,
"Py" ) ) {
167 }
else if ( !strcmp( eletype,
"Pr" ) ) {
169 }
else if ( !strcmp( eletype,
"Qd" ) ) {
171 }
else if ( !strcmp( eletype,
"Tr" ) ) {
179 template<
class S,
class POINT>
inline 183 for(
size_t i=0; i<pts.
size(); i++)
184 pts[i] = {arr[i*3+0], arr[i*3+1], arr[i*3+2]};
243 for(
size_t i=0; i<a.
size(); i++)
245 T aidx = a[i], bidx = b[i];
282 size_t lsize = idx.
size();
296 T val = it != _fwd.
end() ? it->second : T(-1);
298 while( (c < lsize) && (cc == idx[c]) )
318 size_t lsize = idx.
size();
332 T val = it != _bwd.
end() ? *it : T(-1);
334 while( (c < lsize) && (cc == idx[c]) )
351 return _fwd.
count(idx) == 1;
355 return _bwd.
count(idx) == 1;
382 template<
class T,
class S>
414 std::map< SF_nbr, vector<T> >
nbr;
421 meshdata(): g_numelem(0), l_numelem(0), g_numpts(0), l_numpts(0),
434 if(nbr.count(nbr_type) == 0)
437 return nbr[nbr_type];
454 typename std::map< SF_nbr, vector<T> >::iterator it = nbr.find(nbr_type);
455 assert(it != nbr.end());
472 typename std::map< SF_nbr, vector<T> >::const_iterator it = nbr.find(nbr_type);
473 assert(it != nbr.end());
487 vector<T> & nod = this->register_numbering(nbr_type);
490 for(
const T & c : this->con) g2l[c] = 0;
497 for(
auto it = g2l.
begin(); it != g2l.
end(); ++it, widx++) {
498 nod[widx] = it->first;
503 this->l_numpts = nod.
size();
516 vector<T> & nod = this->get_numbering(nbr_type);
517 for(T & c : con) c = nod[c];
570 template<
class T,
class S>
595 template<
class T,
class S>
602 for(
size_t nidx=0; nidx < n2n_cnt.
size(); nidx++)
603 if(nmax < n2n_cnt[nidx]) nmax = n2n_cnt[nidx];
605 MPI_Allreduce(MPI_IN_PLACE, &nmax, 1, MPI_INT, MPI_MAX, mesh.
comm);
610 template<
class T,
class S>
613 const vector<T> & alg_nod = mesh.
pl.algebraic_nodes();
616 for(
const T & n : alg_nod) alg_mask[n] =
true;
643 for(
size_t i=0; i<scnt.
size(); i++)
644 scnt[i] *= fac, rcnt[i] *= fac;
646 for(
size_t i=0; i<sdsp.
size(); i++)
647 sdsp[i] *= fac, rdsp[i] *= fac;
654 scnt = rcnt, sdsp = rdsp;
655 rcnt = tcnt, rdsp = tdsp;
670 MPI_Comm_size(comm, &size); MPI_Comm_rank(comm, &rank);
675 MPI_Alltoall(scnt.
data(),
sizeof(T), MPI_BYTE, rcnt.
data(),
sizeof(T), MPI_BYTE, comm);
692 MPI_Comm_size(comm, &size); MPI_Comm_rank(comm, &rank);
697 for(
size_t i=0; i<dest.
size(); i++) scnt[dest[i]] += cnt[i];
699 MPI_Alltoall(scnt.
data(),
sizeof(T), MPI_BYTE, rcnt.
data(),
sizeof(T), MPI_BYTE, comm);
716 size_t size = rcnt.
size();
719 for(
size_t i=0, widx=0; i<size; i++)
720 for(T j=0; j<rcnt[i]; j++, widx++)
731 template<
class T,
class S>
745 template<
class T,
class S,
class V>
761 bool operator<(const struct tuple<T> & lhs,
const struct tuple<T> & rhs)
764 return lhs.
v1 < rhs.v1;
766 return lhs.v2 < rhs.v2;
768 template<
class T,
class S>
769 bool operator<(const struct mixed_tuple<T, S> & lhs,
const struct mixed_tuple<T, S> & rhs)
772 return lhs.
v1 < rhs.v1;
774 return lhs.v2 < rhs.v2;
778 bool operator<(const struct triple<T> & lhs,
const struct triple<T> & rhs)
781 return lhs.
v1 < rhs.v1;
782 else if(lhs.v2 != rhs.v2)
783 return lhs.v2 < rhs.v2;
785 return lhs.v3 < rhs.v3;
789 bool operator<(const struct quadruple<T> & lhs,
const struct quadruple<T> & rhs)
792 return lhs.
v1 < rhs.v1;
793 else if(lhs.v2 != rhs.v2)
794 return lhs.v2 < rhs.v2;
795 else if(lhs.v3 != rhs.v3)
796 return lhs.v3 < rhs.v3;
798 return lhs.v4 < rhs.v4;
812 return a.
v1 == b.
v1 && a.
v2 == b.
v2;
881 void sort_triple(
const T in1,
const T in2,
const T in3, T & out1, T & out2, T & out3)
883 bool t12 = in1 < in2;
884 bool t13 = in1 < in3;
885 bool t23 = in2 < in3;
891 out1=in1; out2=in2; out3=in3;
897 out1=in1; out2=in3; out3=in2;
901 out1=in3; out2=in1; out3=in2;
911 out1=in2; out2=in1; out3=in3;
915 out1=in2; out2=in3; out3=in1;
920 out1=in3; out2=in2; out3=in1;
void generate_par_layout()
Set up the parallel layout.
The mesh storage class. It contains both element and vertex data.
hm_uint mkhash(hm_uint a, hm_uint b)
The vector class and related algorithms.
non_overlapping_layout< T > epl
element parallel layout
vector< T > & get_numbering(SF_nbr nbr_type)
Get the vector defining a certain numbering.
Point normalize(const Point &vect)
void sort_triple(const T in1, const T in2, const T in3, T &out1, T &out2, T &out3)
sort the "in" triple into the "out" triple
Point arr_to_point(V *arr)
void multiply_connectivities(const vector< T > &a_cnt, const vector< T > &a_con, const vector< T > &b_cnt, const vector< T > &b_con, vector< T > &c_cnt, vector< T > &c_con)
unsigned long int hm_uint
static hm_uint hash(SF::tuple< T > a)
MPI_Comm comm
the parallel mesh is defined on a MPI world
Dense matrix class and associated funcs.
vector< T > tag
element tag
static hm_uint hash(SF::triple< T > a)
PETSc numbering of nodes.
vector< T > dsp
connectivity starting index of each element
Point orthogonalize(const Point &a, const Point &b)
void interval(vector< T > &vec, size_t start, size_t end)
Create an integer interval between start and end.
T backward_map(T idx) const
Map one index from b to a.
overlapping_layout< T > pl
nodal parallel layout
vector< S > fib
fiber direction
void resize(size_t size)
Resize all vectors to size.
void array_to_points(const vector< S > &arr, vector< POINT > &pts)
The nodal numbering of the reference mesh (the one stored on HD).
std::map< SF_nbr, vector< T > > nbr
container for different numberings
void count(const vector< T > &data, vector< S > &cnt)
Count number of occurrences of indices.
double inner_prod(const Point &a, const Point &b)
index_mapping()
empty constructor.
vector< T > rdsp
Displacements w.r.t. rcnt.
dmat< S > operator-(const dmat< S > &a, const dmat< S > &b)
Sequential linear algebra kernels.
vector< S > xyz
node cooridnates
T * data()
Pointer to the vector's start.
void point_to_arr(Point &p, V *arr)
const vector< T > & get_numbering(SF_nbr nbr_type) const
Get the vector defining a certain numbering.
SF_nbr
Enumeration encoding the different supported numberings.
void transpose_connectivity(const vector< T > &a_cnt, const vector< T > &a_con, vector< T > &b_cnt, vector< T > &b_con)
Transpose CRS matrix graph A into B.
#define SF_COMM
the default SlimFem MPI communicator
The element numbering of the reference mesh (the one stored on HD).
const hashmap::unordered_map< T, T > & get_bwd_map() const
Submesh nodal numbering: The globally ascending sorted reference indices are reindexed.
void get_alg_mask(const meshdata< T, S > &mesh, vector< bool > &alg_mask)
int max_nodal_edgecount(const meshdata< T, S > &mesh)
Compute the maximum number of node-to-node edges for a mesh.
The class holds the communication graph for a MPI_Exchange() call.
Point project(const Point &a, const Point &b)
project b onto a
double distance(const Point &a, const Point &b)
void sort(Compare comp=Compare())
Sort data entries.
void cnt_from_dsp(const vector< T > &dsp, vector< T > &cnt)
Compute counts from displacements.
void configure(const vector< V > &dest, const vector< T > &cnt, MPI_Comm comm)
Set up the communication graph.
void clear_data()
Clear the mesh data from memory.
Submesh element numbering: The globally ascending sorted reference indices are reindexed.
hashmap::unordered_set< int > extr_tag
the element tags based on which the mesh has been extracted
elem_t getElemTypeID(char *eletype)
Generate element type enum from string.
bool in_b(const T idx)
return whether idx is in set B
static bool cmp(SF::tuple< T > a, SF::tuple< T > b)
void configure(const vector< V > &dest, MPI_Comm comm)
Set up the communication graph.
T sum(const vector< T > &vec)
Compute sum of a vector's entries.
dmat< S > operator+(const dmat< S > &a, const dmat< S > &b)
meshdata()
construct empty mesh
static hm_uint hash(SF::quadruple< T > a)
dmat< S > operator*(const dmat< S > &a, const dmat< S > &b)
static bool cmp(SF::quadruple< T > a, SF::quadruple< T > b)
void binary_sort_copy(vector< T > &_V, vector< S > &_W)
Index mapping class. This is a bijective mapping.
void global_to_local(const vector< T > &glob, vector< T > &data, bool sortedData, bool doWarn)
hashmap::unordered_map< T, T > & get_bwd_map()
void scale(V fac)
scale comm graph layout data
void outer_prod(const Point &a, const Point &b, const double s, double *buff, const bool add=false)
static bool cmp(SF::triple< T > a, SF::triple< T > b)
void localize(SF_nbr nbr_type)
Localize the connectivity data w.r.t. a given numbering.
std::string name
the mesh name
hm_int count(const K &key) const
Check if key exists.
size_t size() const
The current size of the vector.
vector< T > scnt
Number of elements sent to each rank.
void nodal_connectivity_graph(const meshdata< T, S > &mesh, vector< T > &n2n_cnt, vector< T > &n2n_con)
Compute the node-to-node connectivity.
The parallel layout of non overlapping indices.
vector< T > sdsp
Displacements w.r.t. scnt.
size_t num_global_idx() const
Retrieve the global number of indices.
index_mapping(const vector< T > &a, const vector< T > &b)
Constructor that uses assign() to set up the index mapping.
const hm_uint mkhash_init
void assign(const vector< T > &a, const vector< T > &b)
Set up the index mapping between a and b.
size_t l_numpts
local number of points
hashmap::unordered_map< T, T > & get_fwd_map()
void assign(const vector< T > &idx, MPI_Comm comm)
Initialization function.
size_t l_numelem
local number of elements
size_t g_numpts
global number of points
void forward_map(vector< T > &idx) const
Map a whole array of indices in a to indices in b.
void backward_map(vector< T > &idx) const
Map a whole array of indices in b to indices in a.
void assign(InputIterator s, InputIterator e)
Assign a memory range.
void source_ranks(vector< V > &source)
For every received data element, get the rank indices it was receive from.
iterator find(const K &key)
Search for key. Return iterator.
T forward_map(T idx) const
Map one index from a to b.
vector< S > she
sheet direction
bool in_a(const T idx)
return whether idx is in set A
void transpose()
transpose comm graph (receive becomes send, and vice versa)
size_t g_numelem
global number of elements
const hashmap::unordered_map< T, T > & get_fwd_map() const
vector< T > rcnt
Number of elements received from each rank.
vector< T > & register_numbering(SF_nbr nbr_type)
Register a new numbering to the mesh and return the associated index vector.
void assign(const vector< T > &ref_eidx, MPI_Comm comm)
Generate the layout.
size_t size() const
number of entries in bijective map
Point cross(const Point &a, const Point &b)
cross product
void resize(size_t n)
Resize a vector.
vector< elem_t > type
element type
Classes similar to unordered_set and unordered_map, but with better performance.
The overlapping_layout class contains the algorithms related to managing overlapping parallel index s...
void dsp_from_cnt(const vector< T > &cnt, vector< T > &dsp)
Compute displacements from counts.
void globalize(SF_nbr nbr_type)
Localize the connectivity data w.r.t. a given numbering.
double mag(const Point &vect)
vector magnitude