Skip to content
GitLab
Menu
Projects
Groups
Snippets
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
iwr
dune-vtk
Commits
a9022011
Commit
a9022011
authored
Jun 23, 2021
by
Praetorius, Simon
Browse files
Update LagrangePoints due to changes in dune-geometry and dune-localfunctions
parent
cb00bc20
Changes
5
Hide whitespace changes
Inline
Side-by-side
dune/vtk/gridcreators/lagrangegridcreator.hh
View file @
a9022011
...
...
@@ -226,7 +226,11 @@ namespace Dune
int
order
(
GeometryType
type
,
std
::
size_t
nNodes
)
const
{
for
(
int
o
=
1
;
o
<=
int
(
nNodes
);
++
o
)
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
if
(
numLagrangePoints
(
type
.
id
(),
type
.
dim
(),
o
)
==
std
::
size_t
(
nNodes
))
#else
if
(
numLagrangePoints
(
type
,
o
)
==
std
::
size_t
(
nNodes
))
#endif
return
o
;
return
1
;
...
...
dune/vtk/gridfunctions/lagrangegridfunction.hh
View file @
a9022011
...
...
@@ -4,6 +4,7 @@
#include <vector>
#include <dune/common/dynvector.hh>
#include <dune/common/version.hh>
#include <dune/localfunctions/lagrange.hh>
#include <dune/vtk/gridfunctions/common.hh>
#include <dune/vtk/gridfunctions/continuousgridfunction.hh>
...
...
@@ -83,7 +84,11 @@ namespace Dune
std
::
int64_t
shift
=
(
insertionIndex
==
0
?
0
:
(
*
offsets_
)[
insertionIndex
-
1
]);
std
::
int64_t
numNodes
=
(
*
offsets_
)[
insertionIndex
]
-
shift
;
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
[[
maybe_unused
]]
std
::
int64_t
maxNumNodes
=
numLagrangePoints
(
element
.
type
().
id
(),
element
.
type
().
dim
(),
20
);
#else
[[
maybe_unused
]]
std
::
int64_t
maxNumNodes
=
numLagrangePoints
(
element
.
type
(),
20
);
#endif
VTK_ASSERT
(
numNodes
>
0
&&
numNodes
<
maxNumNodes
);
int
order
=
creator_
->
order
(
element
.
type
(),
numNodes
);
...
...
dune/vtk/utility/lagrangepoints.hh
View file @
a9022011
...
...
@@ -4,6 +4,7 @@
#include <array>
#include <dune/common/exceptions.hh>
#include <dune/common/version.hh>
#include <dune/geometry/type.hh>
#include <dune/localfunctions/lagrange/equidistantpoints.hh>
...
...
@@ -47,16 +48,29 @@ namespace Dune
}
/// Fill the lagrange points for the given topology type `Topology`
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
template
<
class
Topology
>
bool
build
()
{
build
(
GeometryType
(
Topology
{}));
return
true
;
}
#else
template
<
GeometryType
::
Id
geometryId
>
bool
build
()
{
build
(
GeometryType
(
geometryId
));
return
true
;
}
#endif
/// Returns whether the point set support the given topology type `Topology` and can
/// generate point for the given order.
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
template
<
class
Topology
>
#else
template
<
GeometryType
::
Id
geometryId
>
#endif
static
bool
supports
(
std
::
size_t
order
)
{
return
true
;
...
...
dune/vtk/utility/lagrangepoints.impl.hh
View file @
a9022011
...
...
@@ -4,6 +4,7 @@
#include <array>
#include <dune/common/exceptions.hh>
#include <dune/common/version.hh>
#include <dune/geometry/type.hh>
#include <dune/localfunctions/lagrange/equidistantpoints.hh>
...
...
@@ -89,7 +90,11 @@ template <class K>
template
<
class
Points
>
void
LagrangePointSetBuilder
<
K
,
2
>::
operator
()
(
GeometryType
gt
,
int
order
,
Points
&
points
)
const
{
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
std
::
size_t
nPoints
=
numLagrangePoints
(
gt
.
id
(),
dim
,
order
);
#else
std
::
size_t
nPoints
=
numLagrangePoints
(
gt
,
order
);
#endif
if
(
gt
.
isTriangle
())
buildTriangle
(
nPoints
,
order
,
points
);
...
...
@@ -200,8 +205,8 @@ void LagrangePointSetBuilder<K,2>::buildQuad(std::size_t nPoints, int order, Poi
const
K
r
=
K
(
m
)
/
orders
[
0
];
const
K
s
=
K
(
n
)
/
orders
[
1
];
Vec
p
=
(
1.0
-
r
)
*
(
nodes
[
3
]
*
s
+
nodes
[
0
]
*
(
1.0
-
s
))
+
r
*
(
nodes
[
2
]
*
s
+
nodes
[
1
]
*
(
1.0
-
s
));
Vec
p
=
K
(
1.0
-
r
)
*
(
nodes
[
3
]
*
s
+
nodes
[
0
]
*
K
(
1.0
-
s
))
+
r
*
(
nodes
[
2
]
*
s
+
nodes
[
1
]
*
K
(
1.0
-
s
));
auto
[
idx
,
key
]
=
calcQuadKey
(
m
,
n
,
orders
);
points
[
idx
]
=
LP
{
p
,
key
};
...
...
@@ -262,7 +267,11 @@ template <class K>
template
<
class
Points
>
void
LagrangePointSetBuilder
<
K
,
3
>::
operator
()
(
GeometryType
gt
,
unsigned
int
order
,
Points
&
points
)
const
{
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
std
::
size_t
nPoints
=
numLagrangePoints
(
gt
.
id
(),
dim
,
order
);
#else
std
::
size_t
nPoints
=
numLagrangePoints
(
gt
,
order
);
#endif
if
(
gt
.
isTetrahedron
())
buildTetra
(
nPoints
,
order
,
points
);
...
...
@@ -410,8 +419,8 @@ void LagrangePointSetBuilder<K,3>::buildHex (std::size_t nPoints, int order, Poi
const
K
r
=
K
(
m
)
/
orders
[
0
];
const
K
s
=
K
(
n
)
/
orders
[
1
];
const
K
t
=
K
(
p
)
/
orders
[
2
];
Vec
point
=
(
1.0
-
r
)
*
((
nodes
[
3
]
*
(
1.0
-
t
)
+
nodes
[
7
]
*
t
)
*
s
+
(
nodes
[
0
]
*
(
1.0
-
t
)
+
nodes
[
4
]
*
t
)
*
(
1.0
-
s
))
+
r
*
((
nodes
[
2
]
*
(
1.0
-
t
)
+
nodes
[
6
]
*
t
)
*
s
+
(
nodes
[
1
]
*
(
1.0
-
t
)
+
nodes
[
5
]
*
t
)
*
(
1.0
-
s
));
Vec
point
=
K
(
1.0
-
r
)
*
((
nodes
[
3
]
*
K
(
1.0
-
t
)
+
nodes
[
7
]
*
t
)
*
s
+
(
nodes
[
0
]
*
K
(
1.0
-
t
)
+
nodes
[
4
]
*
t
)
*
K
(
1.0
-
s
))
+
r
*
((
nodes
[
2
]
*
K
(
1.0
-
t
)
+
nodes
[
6
]
*
t
)
*
s
+
(
nodes
[
1
]
*
K
(
1.0
-
t
)
+
nodes
[
5
]
*
t
)
*
K
(
1.0
-
s
));
auto
[
idx
,
key
]
=
calcHexKey
(
m
,
n
,
p
,
orders
);
points
[
idx
]
=
LP
{
point
,
key
};
...
...
dune/vtk/utility/test/test-lagrange.cc
View file @
a9022011
...
...
@@ -2,6 +2,7 @@
// vi: set et ts=4 sw=2 sts=2:
#include <config.h>
#include <dune/common/version.hh>
#include <dune/localfunctions/utility/field.hh>
#include <dune/localfunctions/utility/basisprint.hh>
...
...
@@ -72,24 +73,43 @@ bool test(const Basis &basis, const Points &points, bool verbose)
return
ret
;
}
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
template
<
class
Topology
>
#else
template
<
Dune
::
GeometryType
::
Id
geometryId
>
#endif
bool
test
(
unsigned
int
order
,
bool
verbose
=
false
)
{
typedef
Dune
::
LagrangeBasisFactory
<
Dune
::
Vtk
::
LagrangePointSet
,
Topology
::
dimension
,
StorageField
,
ComputeField
>
BasisFactory
;
typedef
Dune
::
LagrangeCoefficientsFactory
<
Dune
::
Vtk
::
LagrangePointSet
,
Topology
::
dimension
,
double
>
LagrangeCoefficientsFactory
;
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
typedef
Dune
::
LagrangeBasisFactory
<
Dune
::
Vtk
::
LagrangePointSet
,
Topology
::
dimension
,
StorageField
,
ComputeField
>
BasisFactory
;
typedef
Dune
::
LagrangeCoefficientsFactory
<
Dune
::
Vtk
::
LagrangePointSet
,
Topology
::
dimension
,
double
>
LagrangeCoefficientsFactory
;
#else
constexpr
Dune
::
GeometryType
geometry
=
geometryId
;
typedef
Dune
::
LagrangeBasisFactory
<
Dune
::
Vtk
::
LagrangePointSet
,
geometry
.
dim
(),
StorageField
,
ComputeField
>
BasisFactory
;
typedef
Dune
::
LagrangeCoefficientsFactory
<
Dune
::
Vtk
::
LagrangePointSet
,
geometry
.
dim
(),
double
>
LagrangeCoefficientsFactory
;
#endif
bool
ret
=
true
;
for
(
unsigned
int
o
=
0
;
o
<=
order
;
++
o
)
{
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
const
typename
LagrangeCoefficientsFactory
::
Object
*
pointsPtr
=
LagrangeCoefficientsFactory
::
template
create
<
Topology
>(
o
);
#else
const
typename
LagrangeCoefficientsFactory
::
Object
*
pointsPtr
=
LagrangeCoefficientsFactory
::
template
create
<
geometry
>(
o
);
#endif
if
(
pointsPtr
==
0
)
continue
;
std
::
cout
<<
"# Testing "
<<
Topology
::
name
()
<<
" in dimension "
<<
Topology
::
dimension
<<
" with order "
<<
o
<<
std
::
endl
;
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
std
::
cout
<<
"Testing "
<<
Topology
::
name
()
<<
" in dimension "
<<
Topology
::
dimension
<<
" with order "
<<
o
<<
std
::
endl
;
typename
BasisFactory
::
Object
&
basis
=
*
BasisFactory
::
template
create
<
Topology
>(
o
);
#else
std
::
cout
<<
"Testing "
<<
geometry
<<
" with order "
<<
o
<<
std
::
endl
;
typename
BasisFactory
::
Object
&
basis
=
*
BasisFactory
::
template
create
<
geometry
>(
o
);
#endif
ret
|=
test
(
basis
,
*
pointsPtr
,
verbose
);
...
...
@@ -97,10 +117,17 @@ bool test(unsigned int order, bool verbose = false)
// derivatives in a human readabible form (aka LaTeX source)
#ifdef TEST_OUTPUT_FUNCTIONS
std
::
stringstream
name
;
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
name
<<
"lagrange_"
<<
Topology
::
name
()
<<
"_p"
<<
o
<<
".basis"
;
std
::
ofstream
out
(
name
.
str
().
c_str
());
Dune
::
basisPrint
<
0
,
BasisFactory
,
typename
BasisFactory
::
StorageField
>
(
out
,
basis
);
Dune
::
basisPrint
<
1
,
BasisFactory
,
typename
BasisFactory
::
StorageField
>
(
out
,
basis
);
#else
name
<<
"lagrange_"
<<
geometry
<<
"_p"
<<
o
<<
".basis"
;
std
::
ofstream
out
(
name
.
str
().
c_str
());
Dune
::
basisPrint
<
0
,
BasisFactory
,
typename
BasisFactory
::
StorageField
,
geometry
>
(
out
,
basis
);
Dune
::
basisPrint
<
1
,
BasisFactory
,
typename
BasisFactory
::
StorageField
,
geometry
>
(
out
,
basis
);
#endif
#endif // TEST_OUTPUT_FUNCTIONS
LagrangeCoefficientsFactory
::
release
(
pointsPtr
);
...
...
@@ -149,19 +176,36 @@ int main ( int argc, char **argv )
bool
tests
=
true
;
#ifdef CHECKDIM1
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
tests
&=
test
<
Prism
<
Point
>
>
(
order
);
// line
tests
&=
test
<
Pyramid
<
Point
>
>
(
order
);
// line
#else
tests
&=
test
<
GeometryTypes
::
cube
(
1
)
>
(
order
);
// line
tests
&=
test
<
GeometryTypes
::
simplex
(
1
)
>
(
order
);
// line
#endif
#endif
#ifdef CHECKDIM2
//tests &= test<Prism<Prism<Point> > > (order); // quad
//tests &= test<Pyramid<Pyramid<Point> > >(order); // triangle
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
tests
&=
test
<
Prism
<
Prism
<
Point
>
>
>
(
order
);
// quad
tests
&=
test
<
Pyramid
<
Pyramid
<
Point
>
>
>
(
order
);
// triangle
#else
tests
&=
test
<
GeometryTypes
::
cube
(
2
)
>
(
order
);
// quad
tests
&=
test
<
GeometryTypes
::
simplex
(
2
)
>
(
order
);
// triangle
#endif
#endif
#ifdef CHECKDIM3
// tests &= test<Prism<Prism<Prism<Point> > > >(order); // hexahedron
#if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
tests
&=
test
<
Prism
<
Prism
<
Prism
<
Point
>
>
>
>
(
order
);
// hexahedron
// tests &= test<Prism<Pyramid<Pyramid<Point> > > >(order);
// tests &= test<Pyramid<Pyramid<Pyramid<Point> > > >(order); // tetrahedron
tests
&=
test
<
Pyramid
<
Pyramid
<
Pyramid
<
Point
>
>
>
>
(
order
);
// tetrahedron
#else
tests
&=
test
<
GeometryTypes
::
cube
(
3
)
>
(
order
);
// hexahedron
// tests &= test<GeometryTypes::prism>(order);
// tests &= test<GeometryTypes::pyramid>(order);
tests
&=
test
<
GeometryTypes
::
simplex
(
3
)
>
(
order
);
#endif
#endif
return
(
tests
?
0
:
1
);
...
...
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
.
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment