11#ifndef GOOSEFEM_VECTORPARTITIONED_H
12#define GOOSEFEM_VECTORPARTITIONED_H
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) {
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) {
180 this->
copy_u(nodevec_src, ret);
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) {
207 nodevec_dest(m, i) = nodevec_src(m, i);
231 this->
copy_p(nodevec_src, ret);
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) {
258 nodevec_dest(m, i) = nodevec_src(m, i);
298#pragma omp parallel for
299 for (
size_t d = 0; d <
m_nnu; ++d) {
300 dofval(
m_iiu(d)) = dofval_u(d);
303#pragma omp parallel for
304 for (
size_t d = 0; d <
m_nnp; ++d) {
305 dofval(
m_iip(d)) = dofval_p(d);
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) {
349 nodevec(m, i) = dofval_u(
m_part(m, i));
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) {
439#pragma omp parallel for
440 for (
size_t d = 0; d <
m_nnu; ++d) {
441 dofval_u(d) = dofval(
m_iiu(d));
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) {
479 dofval_u(
m_part(m, i)) = nodevec(m, i);
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) {
557#pragma omp parallel for
558 for (
size_t d = 0; d <
m_nnp; ++d) {
559 dofval_p(d) = dofval(
m_iip(d));
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) {
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) {
Methods to switch between storage types based on a mesh.
Reorder to lowest possible index, in specific order.
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.
array_type::tensor< size_t, 2 > m_dofs
See dofs()
const array_type::tensor< size_t, 2 > & conn() const
std::array< size_t, 2 > shape_nodevec() const
Shape of "nodevec".
array_type::tensor< size_t, 2 > m_conn
See conn()
const array_type::tensor< size_t, 2 > & dofs() const
#define GOOSEFEM_ASSERT(expr)
All assertions are implementation as::
xt::xtensor< T, N > tensor
Fixed (static) rank array.
Toolbox to perform finite element computations.
bool is_unique(const T &arg)
Returns true is a list is unique (has not duplicate items).