9#ifndef GOOSEFEM_MESHTRI3_H
10#define GOOSEFEM_MESHTRI3_H
46 m_nnode = (m_nelx + 1) * (m_nely + 1);
47 m_nelem = m_nelx * m_nely * 2;
54 size_t nelx_impl()
const
59 size_t nely_impl()
const
69 array_type::tensor<double, 2> coor_impl()
const
71 array_type::tensor<double, 2> ret = xt::empty<double>({m_nnode, m_ndim});
73 array_type::tensor<double, 1> x =
74 xt::linspace<double>(0.0, m_h *
static_cast<double>(m_nelx), m_nelx + 1);
75 array_type::tensor<double, 1> y =
76 xt::linspace<double>(0.0, m_h *
static_cast<double>(m_nely), m_nely + 1);
80 for (
size_t iy = 0; iy < m_nely + 1; ++iy) {
81 for (
size_t ix = 0; ix < m_nelx + 1; ++ix) {
82 ret(inode, 0) = x(ix);
83 ret(inode, 1) = y(iy);
91 array_type::tensor<size_t, 2> conn_impl()
const
93 array_type::tensor<size_t, 2> ret = xt::empty<size_t>({m_nelem, m_nne});
97 for (
size_t iy = 0; iy < m_nely; ++iy) {
98 for (
size_t ix = 0; ix < m_nelx; ++ix) {
99 ret(ielem, 0) = (iy) * (m_nelx + 1) + (ix);
100 ret(ielem, 1) = (iy) * (m_nelx + 1) + (ix + 1);
101 ret(ielem, 2) = (iy + 1) * (m_nelx + 1) + (ix);
103 ret(ielem, 0) = (iy) * (m_nelx + 1) + (ix + 1);
104 ret(ielem, 1) = (iy + 1) * (m_nelx + 1) + (ix + 1);
105 ret(ielem, 2) = (iy + 1) * (m_nelx + 1) + (ix);
113 array_type::tensor<size_t, 1> nodesBottomEdge_impl()
const
115 array_type::tensor<size_t, 1> ret = xt::empty<size_t>({m_nelx + 1});
117 for (
size_t ix = 0; ix < m_nelx + 1; ++ix) {
124 array_type::tensor<size_t, 1> nodesTopEdge_impl()
const
126 array_type::tensor<size_t, 1> ret = xt::empty<size_t>({m_nelx + 1});
128 for (
size_t ix = 0; ix < m_nelx + 1; ++ix) {
129 ret(ix) = ix + m_nely * (m_nelx + 1);
135 array_type::tensor<size_t, 1> nodesLeftEdge_impl()
const
137 array_type::tensor<size_t, 1> ret = xt::empty<size_t>({m_nely + 1});
139 for (
size_t iy = 0; iy < m_nely + 1; ++iy) {
140 ret(iy) = iy * (m_nelx + 1);
146 array_type::tensor<size_t, 1> nodesRightEdge_impl()
const
148 array_type::tensor<size_t, 1> ret = xt::empty<size_t>({m_nely + 1});
150 for (
size_t iy = 0; iy < m_nely + 1; ++iy) {
151 ret(iy) = iy * (m_nelx + 1) + m_nelx;
157 array_type::tensor<size_t, 1> nodesBottomOpenEdge_impl()
const
159 array_type::tensor<size_t, 1> ret = xt::empty<size_t>({m_nelx - 1});
161 for (
size_t ix = 1; ix < m_nelx; ++ix) {
168 array_type::tensor<size_t, 1> nodesTopOpenEdge_impl()
const
170 array_type::tensor<size_t, 1> ret = xt::empty<size_t>({m_nelx - 1});
172 for (
size_t ix = 1; ix < m_nelx; ++ix) {
173 ret(ix - 1) = ix + m_nely * (m_nelx + 1);
179 array_type::tensor<size_t, 1> nodesLeftOpenEdge_impl()
const
181 array_type::tensor<size_t, 1> ret = xt::empty<size_t>({m_nely - 1});
183 for (
size_t iy = 1; iy < m_nely; ++iy) {
184 ret(iy - 1) = iy * (m_nelx + 1);
190 array_type::tensor<size_t, 1> nodesRightOpenEdge_impl()
const
192 array_type::tensor<size_t, 1> ret = xt::empty<size_t>({m_nely - 1});
194 for (
size_t iy = 1; iy < m_nely; ++iy) {
195 ret(iy - 1) = iy * (m_nelx + 1) + m_nelx;
201 size_t nodesBottomLeftCorner_impl()
const
206 size_t nodesBottomRightCorner_impl()
const
211 size_t nodesTopLeftCorner_impl()
const
213 return m_nely * (m_nelx + 1);
216 size_t nodesTopRightCorner_impl()
const
218 return m_nely * (m_nelx + 1) + m_nelx;
239inline array_type::tensor<int, 1>
246 size_t nelem = conn.shape(0);
250 for (
size_t ielem = 0; ielem < nelem; ++ielem) {
252 xt::view(coor, conn(ielem, 0), xt::all()) - xt::view(coor, conn(ielem, 1), xt::all());
254 xt::view(coor, conn(ielem, 2), xt::all()) - xt::view(coor, conn(ielem, 1), xt::all());
256 k = v1(0) * v2(1) - v2(0) * v1(1);
283 int orientation = -1)
292 size_t nelem = conn.shape(0);
295 for (
size_t ielem = 0; ielem < nelem; ++ielem) {
296 if ((orientation == -1 && val(ielem) > 0) || (orientation == +1 && val(ielem) < 0)) {
297 std::swap(ret(ielem, 2), ret(ielem, 1));
315 int orientation = -1)
CRTP base class for regular meshes in 2d.
CRTP base class for regular meshes.
auto nely() const
Number of elements in y-direction == height of the mesh, in units of h,.
auto h() const
Linear edge size of one 'block'.
auto nelx() const
Number of elements in x-direction == width of the mesh in units of h.
Regular grid of squares, with each square cut into two triangular elements.
Regular(size_t nelx, size_t nely, double h=1.0)
Constructor.
#define GOOSEFEM_ASSERT(expr)
All assertions are implementation as::
array_type::tensor< size_t, 2 > setOrientation(const array_type::tensor< double, 2 > &coor, const array_type::tensor< size_t, 2 > &conn, const array_type::tensor< int, 1 > &val, int orientation=-1)
Set the orientation of a mesh of triangular elements of type ElementType::Tri3.
array_type::tensor< int, 1 > getOrientation(const array_type::tensor< double, 2 > &coor, const array_type::tensor< size_t, 2 > &conn)
Read the orientation of a mesh of triangular elements of type ElementType::Tri3.
ElementType
Enumerator for element-types.
@ Tri3
Triangle: 3-noded element in 2-d.
xt::xtensor< T, N > tensor
Fixed (static) rank array.
Toolbox to perform finite element computations.