513 2 * this->n_dofs_per_vertex() + this->n_dofs_per_line();
515 this->interface_constraints.TableBase<2, double>::reinit(
516 this->interface_constraints_size());
529 for (
unsigned int i = 0; i <
dofs_1d; ++i)
530 this->interface_constraints(0, i) =
dofs_subcell[0](1, i);
534 for (
unsigned int i = 0; i <
dofs_1d; ++i)
536 this->interface_constraints(1 + c * (this->degree - 1) +
j - 2,
546 this->interface_constraints(0, face_renumber[i]) =
551 this->interface_constraints(1, face_renumber[i]) =
554 this->interface_constraints(2, face_renumber[i]) =
557 this->interface_constraints(3, face_renumber[i]) =
560 this->interface_constraints(4, face_renumber[i]) =
565 for (
unsigned int j = 0;
j < (this->degree - 1); ++
j)
567 this->interface_constraints(5 +
j, face_renumber[i]) =
570 this->interface_constraints(5 + (this->degree - 1) +
j,
574 this->interface_constraints(5 + 2 * (this->degree - 1) +
j,
578 this->interface_constraints(5 + 3 * (this->degree - 1) +
j,
585 for (
unsigned int j = 0;
j < (this->degree - 1); ++
j)
588 this->interface_constraints(5 + 4 * (this->degree - 1) +
j,
592 this->interface_constraints(5 + 4 * (this->degree - 1) +
593 (this->degree - 1) +
j,
598 this->interface_constraints(5 + 4 * (this->degree - 1) +
599 2 * (this->degree - 1) +
j,
603 this->interface_constraints(5 + 4 * (this->degree - 1) +
604 3 * (this->degree - 1) +
j,
609 this->interface_constraints(5 + 4 * (this->degree - 1) +
610 4 * (this->degree - 1) +
j,
614 this->interface_constraints(5 + 4 * (this->degree - 1) +
615 5 * (this->degree - 1) +
j,
620 this->interface_constraints(5 + 4 * (this->degree - 1) +
621 6 * (this->degree - 1) +
j,
625 this->interface_constraints(5 + 4 * (this->degree - 1) +
626 7 * (this->degree - 1) +
j,
633 for (
unsigned int j = 0;
j < (this->degree - 1); ++
j)
634 for (
unsigned int k = 0;
k < (this->degree - 1); ++
k)
637 this->interface_constraints(5 + 12 * (this->degree - 1) +
638 j +
k * (this->degree - 1),
643 this->interface_constraints(5 + 12 * (this->degree - 1) +
644 j +
k * (this->degree - 1) +
651 this->interface_constraints(5 + 12 * (this->degree - 1) +
652 j +
k * (this->degree - 1) +
653 2 * (this->degree - 1) *
659 this->interface_constraints(5 + 12 * (this->degree - 1) +
660 j +
k * (this->degree - 1) +
661 3 * (this->degree - 1) *
687 2 * this->n_dofs_per_vertex() + this->n_dofs_per_line();
690 const std::vector<unsigned int> &renumber =
695 this->prolongation[
iso][c].reinit(this->n_dofs_per_cell(),
696 this->n_dofs_per_cell());
697 this->restriction[
iso][c].reinit(this->n_dofs_per_cell(),
698 this->n_dofs_per_cell());
724 for (
unsigned int j = 0;
j < this->n_dofs_per_cell(); ++
j)
725 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
730 for (
unsigned int c = 0;
735 const unsigned int c0 = c % 2;
737 const unsigned int c1 = c / 2;
739 this->prolongation[
iso][c](
j, i) =
744 (renumber[i] - (renumber[i] %
dofs_1d)) /
747 this->restriction[
iso][c](
j, i) =
752 (renumber[i] - (renumber[i] %
dofs_1d)) /
760 for (
unsigned int c = 0;
765 const unsigned int c0 = c % 2;
767 const unsigned int c1 = (c % 4) / 2;
769 const unsigned int c2 = c / 4;
771 this->prolongation[
iso][c](
j, i) =
789 this->restriction[
iso][c](
j, i) =
1020 const unsigned int subface,
1022 const unsigned int face_no)
const
1071 for (
unsigned int dof = 2;
1072 dof < this->n_dofs_per_face(
face_no);)
1080 for (
unsigned int i = 2; i < this->n_dofs_per_face(
face_no);
1088 for (
unsigned int j = i + 1;
1089 j < this->n_dofs_per_face(
face_no);
1095 if (((i +
j) & 1) != 0
u)
1115 for (
unsigned int dof = 2;
1116 dof < this->n_dofs_per_face(
face_no);)
1126 for (
unsigned int i = 2; i < this->n_dofs_per_face(
face_no);
1134 for (
unsigned int j = i + 1;
1135 j < this->n_dofs_per_face(
face_no);
1163 for (
unsigned int i = 0; i < this->degree - 1;)
1165 for (
unsigned int line = 0;
1169 i + line * (this->degree - 1) +
1172 for (
unsigned int j = 0;
j < this->degree - 1;)
1175 i + (
j + 4) * this->degree -
j) =
1186 for (
unsigned int vertex = 0;
1193 for (
unsigned int i = 2; i <= this->degree; ++i)
1198 i + 2 *
this->degree) = tmp;
1203 i +
this->degree + 1) = tmp;
1205 i + 2 *
this->degree) = tmp;
1207 i + 3 *
this->degree - 1) = tmp;
1210 for (
unsigned int j = 0;
j < this->degree - 1;)
1213 (i + 2) *
this->degree +
j + 2 -
1216 i + (
j + 4) *
this->degree -
j -
1223 for (
unsigned int j = 2;
j <= this->degree; ++
j)
1227 i + (
j + 2) *
this->degree -
j) =
1233 for (
unsigned int k =
j + 1;
k < this->degree; ++
k)
1238 if (((
j +
k) & 1) != 0
u)
1241 i + (
k + 2) *
this->degree -
k) =
1248 i + (
k + 2) *
this->degree -
k) =
1258 for (
unsigned int j = i + 1;
j <= this->degree; ++
j)
1263 if (((i +
j) & 1) != 0
u)
1269 j + 2 *
this->degree) = tmp;
1272 for (
unsigned int k = 2;
k <= this->degree; ++
k)
1276 j + (
k + 2) *
this->degree -
k) =
1282 for (
unsigned int l =
k + 1;
1289 if (((
k + l) & 1) != 0
u)
1292 j + (l + 2) *
this->degree - l) =
1300 j + (l + 2) *
this->degree - l) =
1311 j +
this->degree + 1) = tmp;
1313 j + 2 *
this->degree) = tmp;
1315 j + 3 *
this->degree - 1) =
1319 for (
unsigned int k = 0;
k < this->degree - 1;)
1322 (
j + 2) *
this->degree +
1326 j + (
k + 4) *
this->degree -
k - 2) = tmp;
1336 j + 2 *
this->degree) = tmp;
1339 for (
unsigned int k = 2;
k <= this->degree; ++
k)
1343 j + (
k + 2) *
this->degree -
k) =
1349 for (
unsigned int l =
k + 1;
1356 if (((
k + l) & 1) != 0
u)
1359 j + (l + 2) *
this->degree - l) =
1367 j + (l + 2) *
this->degree - l) =
1378 j +
this->degree + 1) = tmp;
1380 j + 2 *
this->degree) = tmp;
1382 j + 3 *
this->degree - 1) =
1386 for (
unsigned int k = 0;
k < this->degree - 1;)
1389 (
j + 2) *
this->degree +
1393 j + (
k + 4) *
this->degree -
k - 2) = tmp;
1411 for (
unsigned int i = 0; i < this->degree - 1;)
1413 for (
unsigned int line = 0;
1417 i + line * (this->degree - 1) +
1420 for (
unsigned int j = 0;
j < this->degree - 1;)
1423 i + (
j + 4) * this->degree -
j) =
1434 for (
unsigned int vertex = 0;
1441 for (
unsigned int i = 2; i <= this->degree; ++i)
1447 i + 2 *
this->degree) = tmp;
1449 i + 3 *
this->degree - 1) = tmp;
1452 for (
unsigned int j = 0;
j < this->degree - 1;)
1455 j + (i + 2) * this->degree + 2 -
1458 i + (
j + 4) *
this->degree -
j -
1465 i +
this->degree + 1) = tmp;
1467 i + 2 *
this->degree) = tmp;
1472 for (
unsigned int j = i + 1;
j <= this->degree; ++
j)
1479 j + 2 *
this->degree) = tmp;
1482 for (
unsigned int k = 2;
k <= this->degree; ++
k)
1486 j + (
k + 2) *
this->degree -
k) =
1492 for (
unsigned int l =
k + 1; l <= this->degree;
1498 if (((
k + l) & 1) != 0
u)
1501 j + (l + 2) *
this->degree - l) =
1509 j + (l + 2) *
this->degree - l) =
1518 for (
unsigned int k = 0;
k < this->degree - 1;)
1521 j + (
k + 4) *
this->degree -
1528 j + 2 *
this->degree) = tmp;
1530 j + 3 *
this->degree - 1) = tmp;
1532 if (((i +
j) & 1) != 0
u)
1539 j +
this->degree + 1) =
1543 for (
unsigned int k = 0;
k < this->degree - 1;)
1546 k + (
j + 2) * this->degree +
1554 for (
unsigned int j = 2;
j <= this->degree; ++
j)
1558 i + (
j + 2) *
this->degree -
j) =
1564 for (
unsigned int k =
j + 1;
k <= this->degree; ++
k)
1569 if (((
j +
k) & 1) != 0
u)
1572 i + (
k + 2) *
this->degree -
k) =
1579 i + (
k + 2) *
this->degree -
k) =
1597 for (
unsigned int i = 0; i < this->degree - 1;)
1599 for (
unsigned int line = 0;
1603 i + line * (this->degree - 1) +
1606 for (
unsigned int j = 0;
j < this->degree - 1;)
1609 i + (
j + 4) * this->degree -
j) =
1619 for (
unsigned int vertex = 0;
1626 for (
unsigned int i = 2; i <= this->degree; ++i)
1631 i + 3 *
this->degree - 1) = tmp;
1636 i +
this->degree + 1) = tmp;
1638 i + 2 *
this->degree) = tmp;
1640 i + 3 *
this->degree - 1) = tmp;
1643 for (
unsigned int j = 0;
j < this->degree - 1;)
1646 j + (i + 2) *
this->degree + 2 -
1649 i + (
j + 4) *
this->degree -
j -
1656 for (
unsigned int j = 2;
j <= this->degree; ++
j)
1660 i + (
j + 2) *
this->degree -
j) =
1666 for (
unsigned int k =
j + 1;
k <= this->degree; ++
k)
1672 i + (
k + 2) *
this->degree -
k) =
1682 for (
unsigned int j = i + 1;
j <= this->degree; ++
j)
1691 for (
unsigned int k = 0;
k < this->degree - 1;)
1694 k + (
j + 2) *
this->degree +
1703 j +
this->degree + 1) = tmp;
1705 if (((i +
j) & 1) != 0
u)
1709 j + 2 *
this->degree) = tmp;
1711 j + 3 *
this->degree - 1) = tmp;
1714 j + 3 *
this->degree - 1) = tmp;
1717 for (
unsigned int k = 2;
k <= this->degree; ++
k)
1721 j + (
k + 2) *
this->degree -
k) =
1727 for (
unsigned int l =
k + 1; l <= this->degree;
1734 j + (l + 2) *
this->degree - l) =
1742 for (
unsigned int k = 0;
k < this->degree - 1;)
1745 j + (
k + 4) *
this->degree -
1757 for (
unsigned int vertex = 0;
1762 for (
unsigned int i = 0; i < this->degree - 1;)
1764 for (
unsigned int line = 0;
1768 i + line * (this->degree - 1) +
1771 for (
unsigned int j = 0;
j < this->degree - 1;)
1774 i + (
j + 4) * this->degree -
j) =
1792 for (
unsigned int i = 2; i <= this->degree; ++i)
1798 i + 2 *
this->degree) = tmp;
1800 i + 3 *
this->degree - 1) = tmp;
1803 for (
unsigned int j = 0;
j < this->degree - 1;)
1806 j + (i + 2) * this->degree + 2 -
1809 i + (
j + 4) *
this->degree - 2) =
1816 i +
this->degree + 1) = tmp;
1818 i + 3 *
this->degree - 1) = tmp;
1821 for (
unsigned int j = 2;
j <= this->degree; ++
j)
1825 i + (
j + 2) *
this->degree -
j) =
1831 for (
unsigned int k =
j + 1;
k <= this->degree; ++
k)
1837 i + (
k + 2) *
this->degree -
k) =
1847 for (
unsigned int j = i + 1;
j <= this->degree; ++
j)
1857 j + 2 *
this->degree) = tmp;
1859 j + 3 *
this->degree - 1) = tmp;
1862 j +
this->degree + 1) = tmp;
1864 j + 3 *
this->degree - 1) = tmp;
1867 for (
unsigned int k = 2;
k <= this->degree; ++
k)
1871 j + (
k + 2) *
this->degree -
k) =
1877 for (
unsigned int l =
k + 1; l <= this->degree;
1884 j + (l + 2) *
this->degree - l) =
1892 for (
unsigned int k = 0;
k < this->degree - 1;)
1895 k + (
j + 2) * this->degree +
1898 j + (
k + 4) *
this->degree -
1986 std::vector<unsigned int>
h2l(fe.n_dofs_per_cell());
1989 const unsigned int degree = fe.n_dofs_per_line() + 1;
1992 const unsigned int n = degree + 1;
2010 for (
unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
2031 unsigned int next_index = 0;
2033 h2l[next_index++] = 0;
2034 h2l[next_index++] = 1;
2035 h2l[next_index++] = n;
2036 h2l[next_index++] = n + 1;
2038 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2039 h2l[next_index++] = (2 + i) * n;
2041 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2042 h2l[next_index++] = (2 + i) * n + 1;
2044 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2045 h2l[next_index++] = 2 + i;
2047 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2048 h2l[next_index++] = n + 2 + i;
2050 Assert(fe.n_dofs_per_quad(0 ) ==
2051 fe.n_dofs_per_line() * fe.n_dofs_per_line(),
2053 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2054 for (
unsigned int j = 0;
j < fe.n_dofs_per_line(); ++
j)
2055 h2l[next_index++] = (2 + i) * n + 2 +
j;
2064 unsigned int next_index = 0;
2065 const unsigned int n2 = n * n;
2068 h2l[next_index++] = 0;
2069 h2l[next_index++] = 1;
2070 h2l[next_index++] = n;
2071 h2l[next_index++] = n + 1;
2073 h2l[next_index++] =
n2;
2074 h2l[next_index++] =
n2 + 1;
2075 h2l[next_index++] =
n2 + n;
2076 h2l[next_index++] =
n2 + n + 1;
2080 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2081 h2l[next_index++] = (2 + i) * n;
2082 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2083 h2l[next_index++] = (2 + i) * n + 1;
2084 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2085 h2l[next_index++] = 2 + i;
2086 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2087 h2l[next_index++] = n + 2 + i;
2089 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2090 h2l[next_index++] =
n2 + (2 + i) * n;
2091 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2092 h2l[next_index++] =
n2 + (2 + i) * n + 1;
2093 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2094 h2l[next_index++] =
n2 + 2 + i;
2095 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2096 h2l[next_index++] =
n2 + n + 2 + i;
2098 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2099 h2l[next_index++] = (2 + i) *
n2;
2100 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2101 h2l[next_index++] = (2 + i) *
n2 + 1;
2102 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2103 h2l[next_index++] = (2 + i) *
n2 + n;
2104 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2105 h2l[next_index++] = (2 + i) *
n2 + n + 1;
2110 const unsigned int face_no = 0;
2115 fe.n_dofs_per_line() * fe.n_dofs_per_line(),
2118 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2119 for (
unsigned int j = 0;
j < fe.n_dofs_per_line(); ++
j)
2120 h2l[next_index++] = (2 + i) *
n2 + (2 +
j) * n;
2122 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2123 for (
unsigned int j = 0;
j < fe.n_dofs_per_line(); ++
j)
2124 h2l[next_index++] = (2 + i) *
n2 + (2 +
j) * n + 1;
2126 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2127 for (
unsigned int j = 0;
j < fe.n_dofs_per_line(); ++
j)
2128 h2l[next_index++] = (2 + i) *
n2 + 2 +
j;
2130 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2131 for (
unsigned int j = 0;
j < fe.n_dofs_per_line(); ++
j)
2132 h2l[next_index++] = (2 + i) *
n2 + n + 2 +
j;
2134 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2135 for (
unsigned int j = 0;
j < fe.n_dofs_per_line(); ++
j)
2136 h2l[next_index++] = (2 + i) * n + 2 +
j;
2138 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2139 for (
unsigned int j = 0;
j < fe.n_dofs_per_line(); ++
j)
2140 h2l[next_index++] =
n2 + (2 + i) * n + 2 +
j;
2143 Assert(fe.n_dofs_per_hex() ==
2144 fe.n_dofs_per_quad(
face_no) * fe.n_dofs_per_line(),
2146 for (
unsigned int i = 0; i < fe.n_dofs_per_line(); ++i)
2147 for (
unsigned int j = 0;
j < fe.n_dofs_per_line(); ++
j)
2148 for (
unsigned int k = 0;
k < fe.n_dofs_per_line(); ++
k)
2149 h2l[next_index++] = (2 + i) *
n2 + (2 +
j) * n + 2 +
k;