9#ifndef GOOSEFEM_TYINGSPERIODIC_H
10#define GOOSEFEM_TYINGSPERIODIC_H
16#include <Eigen/Sparse>
62 template <
class C,
class D,
class S,
class T>
82 template <
class C,
class D,
class S,
class T,
class U>
86 const S& control_dofs,
92 m_ndim = m_coor.shape(1);
93 m_nties = m_tyings.shape(0);
106 xt::view(
dofs, xt::keep(dependent), xt::all());
107 U
iid = xt::flatten(dependent_dofs);
118 m_dofs = reorder.apply(
dofs);
119 m_control = reorder.apply(control_dofs);
190 return xt::arange<size_t>(m_nni, m_nni + m_nnd);
200 return xt::arange<size_t>(m_nni);
210 return xt::arange<size_t>(m_nnu);
220 return xt::arange<size_t>(m_nnp) + m_nnu;
237 Eigen::SparseMatrix<double>
Cdi()
const
239 std::vector<Eigen::Triplet<double>> data;
241 data.reserve(m_nties * m_ndim * (m_ndim + 1));
243 for (
size_t i = 0; i < m_nties; ++i) {
244 for (
size_t j = 0; j < m_ndim; ++j) {
246 size_t ni = m_tyings(i, 0);
247 size_t nd = m_tyings(i, 1);
249 data.push_back(Eigen::Triplet<double>(i * m_ndim + j, m_dofs(ni, j), +1.0));
251 for (
size_t k = 0; k < m_ndim; ++k) {
252 data.push_back(Eigen::Triplet<double>(
253 i * m_ndim + j, m_control(j, k), m_coor(nd, k) - m_coor(ni, k)));
258 Eigen::SparseMatrix<double>
Cdi;
259 Cdi.resize(m_nnd, m_nni);
260 Cdi.setFromTriplets(data.begin(), data.end());
270 Eigen::SparseMatrix<double>
Cdu()
const
272 std::vector<Eigen::Triplet<double>> data;
274 data.reserve(m_nties * m_ndim * (m_ndim + 1));
276 for (
size_t i = 0; i < m_nties; ++i) {
277 for (
size_t j = 0; j < m_ndim; ++j) {
279 size_t ni = m_tyings(i, 0);
280 size_t nd = m_tyings(i, 1);
282 if (m_dofs(ni, j) < m_nnu) {
283 data.push_back(Eigen::Triplet<double>(i * m_ndim + j, m_dofs(ni, j), +1.0));
286 for (
size_t k = 0; k < m_ndim; ++k) {
287 if (m_control(j, k) < m_nnu) {
288 data.push_back(Eigen::Triplet<double>(
289 i * m_ndim + j, m_control(j, k), m_coor(nd, k) - m_coor(ni, k)));
295 Eigen::SparseMatrix<double>
Cdu;
296 Cdu.resize(m_nnd, m_nnu);
297 Cdu.setFromTriplets(data.begin(), data.end());
307 Eigen::SparseMatrix<double>
Cdp()
const
309 std::vector<Eigen::Triplet<double>> data;
311 data.reserve(m_nties * m_ndim * (m_ndim + 1));
313 for (
size_t i = 0; i < m_nties; ++i) {
314 for (
size_t j = 0; j < m_ndim; ++j) {
316 size_t ni = m_tyings(i, 0);
317 size_t nd = m_tyings(i, 1);
319 if (m_dofs(ni, j) >= m_nnu) {
321 Eigen::Triplet<double>(i * m_ndim + j, m_dofs(ni, j) - m_nnu, +1.0));
324 for (
size_t k = 0; k < m_ndim; ++k) {
325 if (m_control(j, k) >= m_nnu) {
326 data.push_back(Eigen::Triplet<double>(
328 m_control(j, k) - m_nnu,
329 m_coor(nd, k) - m_coor(ni, k)));
335 Eigen::SparseMatrix<double>
Cdp;
336 Cdp.resize(m_nnd, m_nnp);
337 Cdp.setFromTriplets(data.begin(), data.end());
370 template <
class C,
class N>
379 size_t nnode =
coor.shape(0);
380 size_t ndim =
coor.shape(1);
382 m_control_dofs = xt::arange<size_t>(ndim * ndim).reshape({ndim, ndim});
383 m_control_dofs += xt::amax(
dofs)() + 1;
385 m_control_nodes = nnode + xt::arange<size_t>(ndim);
387 m_coor = xt::concatenate(xt::xtuple(
coor, xt::zeros<double>({ndim, ndim})));
388 m_dofs = xt::concatenate(xt::xtuple(
dofs, m_control_dofs));
418 return m_control_dofs;
428 return m_control_nodes;
Reorder to lowest possible index, in specific order.
Add control nodes to an existing system.
const array_type::tensor< double, 2 > & coor() const
Nodal coordinates, for the system with control nodes added to it.
const array_type::tensor< size_t, 2 > & dofs() const
DOF-numbers per node, for the system with control nodes added to it.
const array_type::tensor< size_t, 1 > & controlNodes() const
Node-numbers of the control nodes.
Control(const C &coor, const N &dofs)
Constructor.
const array_type::tensor< size_t, 2 > & controlDofs() const
DOF-numbers of each control node.
Nodal tyings per periodic boundary conditions.
Periodic(const C &coor, const D &dofs, const S &control_dofs, const T &nodal_tyings)
Constructor.
Periodic(const C &coor, const D &dofs, const S &control_dofs, const T &nodal_tyings, const U &iip)
Constructor.
Eigen::SparseMatrix< double > Cdp() const
Prescribed part of the partitioned tying matrix, see Cdi().
Eigen::SparseMatrix< double > Cdu() const
Unknown part of the partitioned tying matrix, see Cdi().
array_type::tensor< size_t, 1 > iip() const
Independent prescribed DOFs.
array_type::tensor< size_t, 1 > iiu() const
Independent unknown DOFs.
array_type::tensor< size_t, 1 > iid() const
Dependent DOFs.
array_type::tensor< size_t, 1 > iii() const
Independent DOFs.
const array_type::tensor< size_t, 2 > & dofs() const
const array_type::tensor< size_t, 2 > & control() const
Eigen::SparseMatrix< double > Cdi() const
Return tying matrix such as to get the dependent DOFs from the independent DOFs as follows.
const array_type::tensor< size_t, 2 > & nodal_tyings() const
Return the applied nodal tyings.
#define GOOSEFEM_ASSERT(expr)
All assertions are implementation as::
xt::xtensor< T, N > tensor
Fixed (static) rank array.
Toolbox to perform finite element computations.