FrictionQPotFEM 0.23.3
Loading...
Searching...
No Matches
VectorPartitioned.h
Go to the documentation of this file.
1
11#ifndef GOOSEFEM_VECTORPARTITIONED_H
12#define GOOSEFEM_VECTORPARTITIONED_H
13
14#include "Mesh.h"
15#include "Vector.h"
16#include "assertions.h"
17#include "config.h"
18
19namespace GooseFEM {
20
45class VectorPartitioned : public Vector {
46protected:
49 size_t m_nnu;
50 size_t m_nnp;
51
61
62public:
63 VectorPartitioned() = default;
64
77 {
79
80 m_iiu = xt::setdiff1d(m_dofs, m_iip);
81 m_nnp = m_iip.size();
82 m_nnu = m_iiu.size();
84
85 GOOSEFEM_ASSERT(xt::amax(m_iip)() <= xt::amax(m_dofs)());
86 }
87
91 size_t nnu() const
92 {
93 return m_nnu;
94 }
95
99 size_t nnp() const
100 {
101 return m_nnp;
102 }
103
108 {
109 return m_iiu;
110 }
111
116 {
117 return m_iip;
118 }
119
126 {
127 array_type::tensor<bool, 2> ret = xt::zeros<bool>(this->shape_nodevec());
128
129#pragma omp parallel for
130 for (size_t m = 0; m < m_nnode; ++m) {
131 for (size_t i = 0; i < m_ndim; ++i) {
132 if (m_part(m, i) < m_nnu) {
133 ret(m, i) = true;
134 }
135 }
136 }
137
138 return ret;
139 }
140
147 {
148 array_type::tensor<bool, 2> ret = xt::zeros<bool>(this->shape_nodevec());
149
150#pragma omp parallel for
151 for (size_t m = 0; m < m_nnode; ++m) {
152 for (size_t i = 0; i < m_ndim; ++i) {
153 if (m_part(m, i) >= m_nnu) {
154 ret(m, i) = true;
155 }
156 }
157 }
158
159 return ret;
160 }
161
176 const array_type::tensor<double, 2>& nodevec_src,
177 const array_type::tensor<double, 2>& nodevec_dest) const
178 {
179 array_type::tensor<double, 2> ret = nodevec_dest;
180 this->copy_u(nodevec_src, ret);
181 return ret;
182 }
183
196 void copy_u(
197 const array_type::tensor<double, 2>& nodevec_src,
198 array_type::tensor<double, 2>& nodevec_dest) const
199 {
200 GOOSEFEM_ASSERT(xt::has_shape(nodevec_src, {m_nnode, m_ndim}));
201 GOOSEFEM_ASSERT(xt::has_shape(nodevec_dest, {m_nnode, m_ndim}));
202
203#pragma omp parallel for
204 for (size_t m = 0; m < m_nnode; ++m) {
205 for (size_t i = 0; i < m_ndim; ++i) {
206 if (m_part(m, i) < m_nnu) {
207 nodevec_dest(m, i) = nodevec_src(m, i);
208 }
209 }
210 }
211 }
212
227 const array_type::tensor<double, 2>& nodevec_src,
228 const array_type::tensor<double, 2>& nodevec_dest) const
229 {
230 array_type::tensor<double, 2> ret = nodevec_dest;
231 this->copy_p(nodevec_src, ret);
232 return ret;
233 }
234
247 void copy_p(
248 const array_type::tensor<double, 2>& nodevec_src,
249 array_type::tensor<double, 2>& nodevec_dest) const
250 {
251 GOOSEFEM_ASSERT(xt::has_shape(nodevec_src, {m_nnode, m_ndim}));
252 GOOSEFEM_ASSERT(xt::has_shape(nodevec_dest, {m_nnode, m_ndim}));
253
254#pragma omp parallel for
255 for (size_t m = 0; m < m_nnode; ++m) {
256 for (size_t i = 0; i < m_ndim; ++i) {
257 if (m_part(m, i) >= m_nnu) {
258 nodevec_dest(m, i) = nodevec_src(m, i);
259 }
260 }
261 }
262 }
263
272 const array_type::tensor<double, 1>& dofval_u,
273 const array_type::tensor<double, 1>& dofval_p) const
274 {
275 array_type::tensor<double, 1> dofval = xt::empty<double>({m_ndof});
276 this->dofsFromParitioned(dofval_u, dofval_p, dofval);
277 return dofval;
278 }
279
288 const array_type::tensor<double, 1>& dofval_u,
289 const array_type::tensor<double, 1>& dofval_p,
290 array_type::tensor<double, 1>& dofval) const
291 {
292 GOOSEFEM_ASSERT(dofval_u.size() == m_nnu);
293 GOOSEFEM_ASSERT(dofval_p.size() == m_nnp);
294 GOOSEFEM_ASSERT(dofval.size() == m_ndof);
295
296 dofval.fill(0.0);
297
298#pragma omp parallel for
299 for (size_t d = 0; d < m_nnu; ++d) {
300 dofval(m_iiu(d)) = dofval_u(d);
301 }
302
303#pragma omp parallel for
304 for (size_t d = 0; d < m_nnp; ++d) {
305 dofval(m_iip(d)) = dofval_p(d);
306 }
307 }
308
319 const array_type::tensor<double, 1>& dofval_u,
320 const array_type::tensor<double, 1>& dofval_p) const
321 {
322 array_type::tensor<double, 2> nodevec = xt::empty<double>({m_nnode, m_ndim});
323 this->nodeFromPartitioned(dofval_u, dofval_p, nodevec);
324 return nodevec;
325 }
326
337 const array_type::tensor<double, 1>& dofval_u,
338 const array_type::tensor<double, 1>& dofval_p,
339 array_type::tensor<double, 2>& nodevec) const
340 {
341 GOOSEFEM_ASSERT(dofval_u.size() == m_nnu);
342 GOOSEFEM_ASSERT(dofval_p.size() == m_nnp);
343 GOOSEFEM_ASSERT(xt::has_shape(nodevec, {m_nnode, m_ndim}));
344
345#pragma omp parallel for
346 for (size_t m = 0; m < m_nnode; ++m) {
347 for (size_t i = 0; i < m_ndim; ++i) {
348 if (m_part(m, i) < m_nnu) {
349 nodevec(m, i) = dofval_u(m_part(m, i));
350 }
351 else {
352 nodevec(m, i) = dofval_p(m_part(m, i) - m_nnu);
353 }
354 }
355 }
356 }
357
368 const array_type::tensor<double, 1>& dofval_u,
369 const array_type::tensor<double, 1>& dofval_p) const
370 {
371 array_type::tensor<double, 3> elemvec = xt::empty<double>({m_nelem, m_nne, m_ndim});
372 this->elementFromPartitioned(dofval_u, dofval_p, elemvec);
373 return elemvec;
374 }
375
386 const array_type::tensor<double, 1>& dofval_u,
387 const array_type::tensor<double, 1>& dofval_p,
388 array_type::tensor<double, 3>& elemvec) const
389 {
390 GOOSEFEM_ASSERT(dofval_u.size() == m_nnu);
391 GOOSEFEM_ASSERT(dofval_p.size() == m_nnp);
392 GOOSEFEM_ASSERT(xt::has_shape(elemvec, {m_nelem, m_nne, m_ndim}));
393
394#pragma omp parallel for
395 for (size_t e = 0; e < m_nelem; ++e) {
396 for (size_t m = 0; m < m_nne; ++m) {
397 for (size_t i = 0; i < m_ndim; ++i) {
398 if (m_part(m_conn(e, m), i) < m_nnu) {
399 elemvec(e, m, i) = dofval_u(m_part(m_conn(e, m), i));
400 }
401 else {
402 elemvec(e, m, i) = dofval_p(m_part(m_conn(e, m), i) - m_nnu);
403 }
404 }
405 }
406 }
407 }
408
418 {
419 array_type::tensor<double, 1> dofval_u = xt::empty<double>({m_nnu});
420 this->asDofs_u(dofval, dofval_u);
421 return dofval_u;
422 }
423
433 const array_type::tensor<double, 1>& dofval,
434 array_type::tensor<double, 1>& dofval_u) const
435 {
436 GOOSEFEM_ASSERT(dofval.size() == m_ndof);
437 GOOSEFEM_ASSERT(dofval_u.size() == m_nnu);
438
439#pragma omp parallel for
440 for (size_t d = 0; d < m_nnu; ++d) {
441 dofval_u(d) = dofval(m_iiu(d));
442 }
443 }
444
453 {
454 array_type::tensor<double, 1> dofval_u = xt::empty<double>({m_nnu});
455 this->asDofs_u(nodevec, dofval_u);
456 return dofval_u;
457 }
458
467 const array_type::tensor<double, 2>& nodevec,
468 array_type::tensor<double, 1>& dofval_u) const
469 {
470 GOOSEFEM_ASSERT(xt::has_shape(nodevec, {m_nnode, m_ndim}));
471 GOOSEFEM_ASSERT(dofval_u.size() == m_nnu);
472
473 dofval_u.fill(0.0);
474
475#pragma omp parallel for
476 for (size_t m = 0; m < m_nnode; ++m) {
477 for (size_t i = 0; i < m_ndim; ++i) {
478 if (m_part(m, i) < m_nnu) {
479 dofval_u(m_part(m, i)) = nodevec(m, i);
480 }
481 }
482 }
483 }
484
493 {
494 array_type::tensor<double, 1> dofval_u = xt::empty<double>({m_nnu});
495 this->asDofs_u(elemvec, dofval_u);
496 return dofval_u;
497 }
498
507 const array_type::tensor<double, 3>& elemvec,
508 array_type::tensor<double, 1>& dofval_u) const
509 {
510 GOOSEFEM_ASSERT(xt::has_shape(elemvec, {m_nelem, m_nne, m_ndim}));
511 GOOSEFEM_ASSERT(dofval_u.size() == m_nnu);
512
513 dofval_u.fill(0.0);
514
515#pragma omp parallel for
516 for (size_t e = 0; e < m_nelem; ++e) {
517 for (size_t m = 0; m < m_nne; ++m) {
518 for (size_t i = 0; i < m_ndim; ++i) {
519 if (m_part(m_conn(e, m), i) < m_nnu) {
520 dofval_u(m_part(m_conn(e, m), i)) = elemvec(e, m, i);
521 }
522 }
523 }
524 }
525 }
526
536 {
537 array_type::tensor<double, 1> dofval_p = xt::empty<double>({m_nnp});
538 this->asDofs_p(dofval, dofval_p);
539 return dofval_p;
540 }
541
551 const array_type::tensor<double, 1>& dofval,
552 array_type::tensor<double, 1>& dofval_p) const
553 {
554 GOOSEFEM_ASSERT(dofval.size() == m_ndof);
555 GOOSEFEM_ASSERT(dofval_p.size() == m_nnp);
556
557#pragma omp parallel for
558 for (size_t d = 0; d < m_nnp; ++d) {
559 dofval_p(d) = dofval(m_iip(d));
560 }
561 }
562
571 {
572 array_type::tensor<double, 1> dofval_p = xt::empty<double>({m_nnp});
573 this->asDofs_p(nodevec, dofval_p);
574 return dofval_p;
575 }
576
585 const array_type::tensor<double, 2>& nodevec,
586 array_type::tensor<double, 1>& dofval_p) const
587 {
588 GOOSEFEM_ASSERT(xt::has_shape(nodevec, {m_nnode, m_ndim}));
589 GOOSEFEM_ASSERT(dofval_p.size() == m_nnp);
590
591 dofval_p.fill(0.0);
592
593#pragma omp parallel for
594 for (size_t m = 0; m < m_nnode; ++m) {
595 for (size_t i = 0; i < m_ndim; ++i) {
596 if (m_part(m, i) >= m_nnu) {
597 dofval_p(m_part(m, i) - m_nnu) = nodevec(m, i);
598 }
599 }
600 }
601 }
602
611 {
612 array_type::tensor<double, 1> dofval_p = xt::empty<double>({m_nnp});
613 this->asDofs_p(elemvec, dofval_p);
614 return dofval_p;
615 }
616
625 const array_type::tensor<double, 3>& elemvec,
626 array_type::tensor<double, 1>& dofval_p) const
627 {
628 GOOSEFEM_ASSERT(xt::has_shape(elemvec, {m_nelem, m_nne, m_ndim}));
629 GOOSEFEM_ASSERT(dofval_p.size() == m_nnp);
630
631 dofval_p.fill(0.0);
632
633#pragma omp parallel for
634 for (size_t e = 0; e < m_nelem; ++e) {
635 for (size_t m = 0; m < m_nne; ++m) {
636 for (size_t i = 0; i < m_ndim; ++i) {
637 if (m_part(m_conn(e, m), i) >= m_nnu) {
638 dofval_p(m_part(m_conn(e, m), i) - m_nnu) = elemvec(e, m, i);
639 }
640 }
641 }
642 }
643 }
644};
645
646} // namespace GooseFEM
647
648#endif
Generic mesh operations.
Methods to switch between storage types based on a mesh.
Reorder to lowest possible index, in specific order.
Definition: Mesh.h:2387
Class to switch between storage types, based on a mesh and DOFs that are partitioned in:
const array_type::tensor< size_t, 1 > & iip() const
array_type::tensor< double, 1 > AsDofs_u(const array_type::tensor< double, 3 > &elemvec) const
Convert "elemvec" to "dofval" (overwrite entries that occur more than once) and extract the unknown "...
void asDofs_u(const array_type::tensor< double, 1 > &dofval, array_type::tensor< double, 1 > &dofval_u) const
Extract the unknown "dofval":
void asDofs_p(const array_type::tensor< double, 3 > &elemvec, array_type::tensor< double, 1 > &dofval_p) const
Convert "elemvec" to "dofval" (overwrite entries that occur more than once) and extract the prescribe...
array_type::tensor< double, 1 > AsDofs_u(const array_type::tensor< double, 1 > &dofval) const
Extract the unknown "dofval":
void asDofs_p(const array_type::tensor< double, 1 > &dofval, array_type::tensor< double, 1 > &dofval_p) const
Extract the prescribed "dofval":
array_type::tensor< double, 1 > DofsFromParitioned(const array_type::tensor< double, 1 > &dofval_u, const array_type::tensor< double, 1 > &dofval_p) const
Combine unknown and prescribed "dofval" into a single "dofval" list.
void elementFromPartitioned(const array_type::tensor< double, 1 > &dofval_u, const array_type::tensor< double, 1 > &dofval_p, array_type::tensor< double, 3 > &elemvec) const
Combine unknown and prescribed "dofval" into a single "dofval" list and directly convert to "elemvec"...
array_type::tensor< bool, 2 > dofs_is_u() const
Per DOF (see Vector::dofs()) list if unknown ("u").
void asDofs_p(const array_type::tensor< double, 2 > &nodevec, array_type::tensor< double, 1 > &dofval_p) const
Convert "nodevec" to "dofval" (overwrite entries that occur more than once) and extract the prescribe...
array_type::tensor< double, 1 > AsDofs_p(const array_type::tensor< double, 3 > &elemvec) const
Convert "elemvec" to "dofval" (overwrite entries that occur more than once) and extract the prescribe...
VectorPartitioned(const array_type::tensor< size_t, 2 > &conn, const array_type::tensor< size_t, 2 > &dofs, const array_type::tensor< size_t, 1 > &iip)
Constructor.
void copy_p(const array_type::tensor< double, 2 > &nodevec_src, array_type::tensor< double, 2 > &nodevec_dest) const
Copy prescribed DOFs from "nodevec" to another "nodvec":
void dofsFromParitioned(const array_type::tensor< double, 1 > &dofval_u, const array_type::tensor< double, 1 > &dofval_p, array_type::tensor< double, 1 > &dofval) const
Combine unknown and prescribed "dofval" into a single "dofval" list.
void copy_u(const array_type::tensor< double, 2 > &nodevec_src, array_type::tensor< double, 2 > &nodevec_dest) const
Copy unknown DOFs from "nodevec" to another "nodvec":
const array_type::tensor< size_t, 1 > & iiu() const
array_type::tensor< size_t, 1 > m_iiu
See iiu()
array_type::tensor< size_t, 2 > m_part
Renumbered DOFs per node, such that.
array_type::tensor< double, 2 > Copy_u(const array_type::tensor< double, 2 > &nodevec_src, const array_type::tensor< double, 2 > &nodevec_dest) const
Copy unknown DOFs from "nodevec" to another "nodvec":
array_type::tensor< double, 2 > NodeFromPartitioned(const array_type::tensor< double, 1 > &dofval_u, const array_type::tensor< double, 1 > &dofval_p) const
Combine unknown and prescribed "dofval" into a single "dofval" list and directly convert to "nodeval"...
array_type::tensor< size_t, 1 > m_iip
See iip()
void asDofs_u(const array_type::tensor< double, 3 > &elemvec, array_type::tensor< double, 1 > &dofval_u) const
Convert "elemvec" to "dofval" (overwrite entries that occur more than once) and extract the unknown "...
array_type::tensor< double, 1 > AsDofs_p(const array_type::tensor< double, 1 > &dofval) const
Extract the prescribed "dofval":
void asDofs_u(const array_type::tensor< double, 2 > &nodevec, array_type::tensor< double, 1 > &dofval_u) const
Convert "nodevec" to "dofval" (overwrite entries that occur more than once) and extract the unknown "...
array_type::tensor< double, 2 > Copy_p(const array_type::tensor< double, 2 > &nodevec_src, const array_type::tensor< double, 2 > &nodevec_dest) const
Copy prescribed DOFs from "nodevec" to another "nodvec":
void nodeFromPartitioned(const array_type::tensor< double, 1 > &dofval_u, const array_type::tensor< double, 1 > &dofval_p, array_type::tensor< double, 2 > &nodevec) const
Combine unknown and prescribed "dofval" into a single "dofval" list and directly convert to "nodeval"...
array_type::tensor< bool, 2 > dofs_is_p() const
Per DOF (see Vector::dofs()) list if prescribed ("p").
array_type::tensor< double, 3 > ElementFromPartitioned(const array_type::tensor< double, 1 > &dofval_u, const array_type::tensor< double, 1 > &dofval_p) const
Combine unknown and prescribed "dofval" into a single "dofval" list and directly convert to "elemvec"...
array_type::tensor< double, 1 > AsDofs_u(const array_type::tensor< double, 2 > &nodevec) const
Convert "nodevec" to "dofval" (overwrite entries that occur more than once) and extract the unknown "...
array_type::tensor< double, 1 > AsDofs_p(const array_type::tensor< double, 2 > &nodevec) const
Convert "nodevec" to "dofval" (overwrite entries that occur more than once) and extract the prescribe...
Class to switch between storage types.
Definition: Vector.h:23
array_type::tensor< size_t, 2 > m_dofs
See dofs()
Definition: Vector.h:752
const array_type::tensor< size_t, 2 > & conn() const
Definition: Vector.h:92
std::array< size_t, 2 > shape_nodevec() const
Shape of "nodevec".
Definition: Vector.h:280
size_t m_ndof
See ndof.
Definition: Vector.h:757
size_t m_nelem
See nelem.
Definition: Vector.h:753
array_type::tensor< size_t, 2 > m_conn
See conn()
Definition: Vector.h:751
size_t m_ndim
See ndim.
Definition: Vector.h:756
const array_type::tensor< size_t, 2 > & dofs() const
Definition: Vector.h:100
size_t m_nne
See nne.
Definition: Vector.h:754
size_t m_nnode
See nnode.
Definition: Vector.h:755
#define GOOSEFEM_ASSERT(expr)
All assertions are implementation as::
Definition: config.h:96
xt::xtensor< T, N > tensor
Fixed (static) rank array.
Definition: config.h:176
Toolbox to perform finite element computations.
Definition: Allocate.h:14
bool is_unique(const T &arg)
Returns true is a list is unique (has not duplicate items).
Definition: assertions.h:20