From bbe37dced510c22bb238563099eab7dae9220813 Mon Sep 17 00:00:00 2001
From: Wenqing Wang <wenqing.wang@ufz.de>
Date: Thu, 3 Aug 2017 12:24:43 +0200
Subject: [PATCH] [utest] Added the analytic element volumes of curved edge
triangle and quadrilateral
elements, respectively, in the tests for gradient of shape functions.
---
Tests/NumLib/FeTestData/TestFeQUAD8.h | 25 ++++++++++++++----------
Tests/NumLib/FeTestData/TestFeQUAD9.h | 27 +++++++++++++++-----------
Tests/NumLib/FeTestData/TestFeTRI6.h | 28 +++++++++++++++------------
3 files changed, 47 insertions(+), 33 deletions(-)
diff --git a/Tests/NumLib/FeTestData/TestFeQUAD8.h b/Tests/NumLib/FeTestData/TestFeQUAD8.h
index 1de4e8bba9b..16231ab0c55 100644
--- a/Tests/NumLib/FeTestData/TestFeQUAD8.h
+++ b/Tests/NumLib/FeTestData/TestFeQUAD8.h
@@ -40,26 +40,31 @@ public:
static const unsigned n_sample_pt_order3 = 3 * 3;
static const unsigned global_dim = MeshElementType::dimension;
- const double perturbation = -0.1;
+ const double a = 1.0;
+ const double perturbation = 0.1;
- /// create a mesh element
+ /// create a 8 node quadrilateral element
MeshElementType* createMeshElement()
{
+ // Convex 8 node quadrilateral element with curved edges.
auto** nodes = new MeshLib::Node*[e_nnodes];
nodes[0] = new MeshLib::Node(0.0, 0.0, 0.0);
- nodes[1] = new MeshLib::Node(1.0, 0.0, 0.0);
- nodes[2] = new MeshLib::Node(1.0, 1.0, 0.0);
- nodes[3] = new MeshLib::Node(0.0, 1.0, 0.0);
+ nodes[1] = new MeshLib::Node(a, 0.0, 0.0);
+ nodes[2] = new MeshLib::Node(a, a, 0.0);
+ nodes[3] = new MeshLib::Node(0.0, a, 0.0);
- nodes[4] = new MeshLib::Node(0.5, perturbation, 0.0);
- nodes[5] = new MeshLib::Node(1.0 - perturbation, 0.5, 0.0);
- nodes[6] = new MeshLib::Node(0.5, 1.0 - perturbation, 0.0);
- nodes[7] = new MeshLib::Node(perturbation, 0.5, 0.0);
+ nodes[4] = new MeshLib::Node(0.5 * a, -perturbation, 0.0);
+ nodes[5] = new MeshLib::Node(a + perturbation, 0.5 * a, 0.0);
+ nodes[6] = new MeshLib::Node(0.5 * a, a + perturbation, 0.0);
+ nodes[7] = new MeshLib::Node(-perturbation, 0.5 * a, 0.0);
return new MeshElementType(nodes);
}
- double getVolume() const { return 1.26666666667; }
+ double getVolume() const
+ {
+ return a * a + 4 * (2.0 * a * perturbation / 3.0);
+ }
};
} // namespace
diff --git a/Tests/NumLib/FeTestData/TestFeQUAD9.h b/Tests/NumLib/FeTestData/TestFeQUAD9.h
index bde49072424..1da551dfabd 100644
--- a/Tests/NumLib/FeTestData/TestFeQUAD9.h
+++ b/Tests/NumLib/FeTestData/TestFeQUAD9.h
@@ -40,27 +40,32 @@ public:
static const unsigned n_sample_pt_order3 = 3 * 3;
static const unsigned global_dim = MeshElementType::dimension;
- const double perturbation = -0.1;
+ const double a = 1.0;
+ const double perturbation = 0.1;
- /// create a mesh element
+ /// create a 9 node quadrilateral element
MeshElementType* createMeshElement()
{
+ // Convex 9 node quadrilateral element with curved edges.
auto** nodes = new MeshLib::Node*[e_nnodes];
nodes[0] = new MeshLib::Node(0.0, 0.0, 0.0);
- nodes[1] = new MeshLib::Node(1.0, 0.0, 0.0);
- nodes[2] = new MeshLib::Node(1.0, 1.0, 0.0);
- nodes[3] = new MeshLib::Node(0.0, 1.0, 0.0);
+ nodes[1] = new MeshLib::Node(a, 0.0, 0.0);
+ nodes[2] = new MeshLib::Node(a, a, 0.0);
+ nodes[3] = new MeshLib::Node(0.0, a, 0.0);
- nodes[4] = new MeshLib::Node(0.5, perturbation, 0.0);
- nodes[5] = new MeshLib::Node(1.0 - perturbation, 0.5, 0.0);
- nodes[6] = new MeshLib::Node(0.5, 1.0 - perturbation, 0.0);
- nodes[7] = new MeshLib::Node(perturbation, 0.5, 0.0);
+ nodes[4] = new MeshLib::Node(0.5 * a, -perturbation, 0.0);
+ nodes[5] = new MeshLib::Node(a + perturbation, 0.5 * a, 0.0);
+ nodes[6] = new MeshLib::Node(0.5 * a, a + perturbation, 0.0);
+ nodes[7] = new MeshLib::Node(-perturbation, 0.5 * a, 0.0);
- nodes[8] = new MeshLib::Node(0.45, 0.25, 0.0);
+ nodes[8] = new MeshLib::Node(0.45 * a, 0.25 * a, 0.0);
return new MeshElementType(nodes);
}
- double getVolume() const { return 1.2666666666667; }
+ double getVolume() const
+ {
+ return a * a + 4 * (2.0 * a * perturbation / 3.0);
+ }
};
} // namespace
diff --git a/Tests/NumLib/FeTestData/TestFeTRI6.h b/Tests/NumLib/FeTestData/TestFeTRI6.h
index 70a4a0ab8e9..121ce5f3b88 100644
--- a/Tests/NumLib/FeTestData/TestFeTRI6.h
+++ b/Tests/NumLib/FeTestData/TestFeTRI6.h
@@ -41,32 +41,36 @@ public:
static const unsigned n_sample_pt_order2 = 3;
static const unsigned n_sample_pt_order3 = 4;
static const unsigned global_dim = MeshElementType::dimension;
+
+ const double a = 1.0;
const double perturbation = 0.1;
- /// create a mesh element
+ /// create a 6 node triangle element
MeshElementType* createMeshElement()
{
+ // Concave 6 node triangle element with curved edges.
auto** nodes = new MeshLib::Node*[e_nnodes];
nodes[0] = new MeshLib::Node(0.0, 0.0, 0.0);
- nodes[1] = new MeshLib::Node(1.0, 0.0, 0.0);
- nodes[2] = new MeshLib::Node(0.0, 1.0, 0.0);
+ nodes[1] = new MeshLib::Node(a, 0.0, 0.0);
+ nodes[2] = new MeshLib::Node(0.0, a, 0.0);
- nodes[3] = new MeshLib::Node(0.5, perturbation, 0.0);
+ nodes[3] = new MeshLib::Node(0.5 * a, perturbation, 0.0);
const double perturbation_hypotenuse =
- 0.5 * perturbation * std::sqrt(2.0);
- nodes[4] = new MeshLib::Node(0.5 - perturbation_hypotenuse,
- 0.5 - perturbation_hypotenuse, 0.0);
- nodes[5] = new MeshLib::Node(perturbation, 0.5, 0.0);
+ 0.5 * perturbation * std::sqrt(2 * a * a);
+ nodes[4] = new MeshLib::Node(0.5 * a - perturbation_hypotenuse,
+ 0.5 * a - perturbation_hypotenuse, 0.0);
+ nodes[5] = new MeshLib::Node(perturbation, 0.5 * a, 0.0);
return new MeshElementType(nodes);
}
double getVolume() const
{
- // The length of hypotenuse is a_h=sqrt(2)
- // Area = 0.5 - 2 * 2 * a * b /3 - 2 * a_h * b /3
- // where a=1.0, b=perturbation
- return 0.27238576250846036;
+ // The length of hypotenuse is a_h sqrt(2 * a * a)
+ // Area = a^2 /2 - 2 * 2 * a * b /3 - 2 * a_h * b /3
+ // where b=perturbation
+ return 0.5 * a * a - 4.0 * a * perturbation / 3.0 -
+ 2.0 * std::sqrt(2 * a * a) * perturbation / 3.0;
}
};
--
GitLab