Skip to content

Commit

Permalink
Use new tangential vectors
Browse files Browse the repository at this point in the history
  • Loading branch information
@kronbichler
kronbichler committed Feb 29, 2020
1 parent 9d62f2b commit 72c67b6
Show file tree
Hide file tree
Showing 3 changed files with 37 additions and 162 deletions.
52 changes: 9 additions & 43 deletions source/fe/mapping_fe_field.cc
View file
Expand Up @@ -581,56 +581,22 @@ MappingFEField<dim, spacedim, VectorType, DoFHandlerType>::compute_face_data(
dim - 1, std::vector<Tensor<1, spacedim>>(n_original_q_points));

// Compute tangentials to the unit cell.
for (unsigned int i = 0; i < data.unit_tangentials.size(); ++i)
data.unit_tangentials[i].resize(n_original_q_points);

if (dim == 2)
{
// ensure a counterclockwise
// orientation of tangentials
static const int tangential_orientation[4] = {-1, 1, 1, -1};
for (const unsigned int i : GeometryInfo<dim>::face_indices())
{
Tensor<1, dim> tang;
tang[1 - i / 2] = tangential_orientation[i];
std::fill(data.unit_tangentials[i].begin(),
data.unit_tangentials[i].end(),
tang);
}
}
else if (dim == 3)
for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
{
for (const unsigned int i : GeometryInfo<dim>::face_indices())
data.unit_tangentials[i].resize(n_original_q_points);
std::fill(data.unit_tangentials[i].begin(),
data.unit_tangentials[i].end(),
GeometryInfo<dim>::unit_tangential_vectors[i][0]);
if (dim > 2)
{
Tensor<1, dim> tang1, tang2;

const unsigned int nd =
GeometryInfo<dim>::unit_normal_direction[i];

// first tangential
// vector in direction
// of the (nd+1)%3 axis
// and inverted in case
// of unit inward normal
tang1[(nd + 1) % dim] =
GeometryInfo<dim>::unit_normal_orientation[i];
// second tangential
// vector in direction
// of the (nd+2)%3 axis
tang2[(nd + 2) % dim] = 1.;

// same unit tangents
// for all quadrature
// points on this face
std::fill(data.unit_tangentials[i].begin(),
data.unit_tangentials[i].end(),
tang1);
data.unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.resize(n_original_q_points);
std::fill(
data.unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.begin(),
data.unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.end(),
tang2);
GeometryInfo<dim>::unit_tangential_vectors[i][1]);
}
}
}
Expand Down
72 changes: 14 additions & 58 deletions source/fe/mapping_manifold.cc
View file
Expand Up @@ -117,67 +117,23 @@ MappingManifold<dim, spacedim>::InternalData::initialize_face(
std::vector<Tensor<1, spacedim>>(n_original_q_points));

// Compute tangentials to the unit cell.
for (unsigned int i = 0; i < unit_tangentials.size(); ++i)
unit_tangentials[i].resize(n_original_q_points);
switch (dim)
for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
{
case 2:
unit_tangentials[i].resize(n_original_q_points);
std::fill(unit_tangentials[i].begin(),
unit_tangentials[i].end(),
GeometryInfo<dim>::unit_tangential_vectors[i][0]);
if (dim > 2)
{
// ensure a counterclockwise
// orientation of tangentials
static const int tangential_orientation[4] = {-1, 1, 1, -1};
for (unsigned int i = 0;
i < GeometryInfo<dim>::faces_per_cell;
++i)
{
Tensor<1, dim> tang;
tang[1 - i / 2] = tangential_orientation[i];
std::fill(unit_tangentials[i].begin(),
unit_tangentials[i].end(),
tang);
}
break;
}
case 3:
{
for (unsigned int i = 0;
i < GeometryInfo<dim>::faces_per_cell;
++i)
{
Tensor<1, dim> tang1, tang2;

const unsigned int nd =
GeometryInfo<dim>::unit_normal_direction[i];

// first tangential
// vector in direction
// of the (nd+1)%3 axis
// and inverted in case
// of unit inward normal
tang1[(nd + 1) % dim] =
GeometryInfo<dim>::unit_normal_orientation[i];
// second tangential
// vector in direction
// of the (nd+2)%3 axis
tang2[(nd + 2) % dim] = 1.;

// same unit tangents
// for all quadrature
// points on this face
std::fill(unit_tangentials[i].begin(),
unit_tangentials[i].end(),
tang1);
std::fill(
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.begin(),
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.end(),
tang2);
}
break;
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.resize(n_original_q_points);
std::fill(
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.begin(),
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.end(),
GeometryInfo<dim>::unit_tangential_vectors[i][1]);
}
default:
Assert(false, ExcNotImplemented());
}
}
}
Expand Down
75 changes: 14 additions & 61 deletions source/fe/mapping_q_generic.cc
View file
Expand Up @@ -841,70 +841,23 @@ MappingQGeneric<dim, spacedim>::InternalData::initialize_face(
std::vector<Tensor<1, spacedim>>(n_original_q_points));

// Compute tangentials to the unit cell.
for (unsigned int i = 0; i < unit_tangentials.size(); ++i)
unit_tangentials[i].resize(n_original_q_points);
switch (dim)
for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
{
case 2:
unit_tangentials[i].resize(n_original_q_points);
std::fill(unit_tangentials[i].begin(),
unit_tangentials[i].end(),
GeometryInfo<dim>::unit_tangential_vectors[i][0]);
if (dim > 2)
{
// ensure a counterclockwise orientation of tangentials
static const int tangential_orientation[4] = {-1, 1, 1, -1};
for (unsigned int i = 0;
i < GeometryInfo<dim>::faces_per_cell;
++i)
{
Tensor<1, dim> tang;
tang[1 - i / 2] = tangential_orientation[i];
std::fill(unit_tangentials[i].begin(),
unit_tangentials[i].end(),
tang);
}

break;
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.resize(n_original_q_points);
std::fill(
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.begin(),
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.end(),
GeometryInfo<dim>::unit_tangential_vectors[i][1]);
}

case 3:
{
for (unsigned int i = 0;
i < GeometryInfo<dim>::faces_per_cell;
++i)
{
Tensor<1, dim> tang1, tang2;

const unsigned int nd =
GeometryInfo<dim>::unit_normal_direction[i];

// first tangential
// vector in direction
// of the (nd+1)%3 axis
// and inverted in case
// of unit inward normal
tang1[(nd + 1) % dim] =
GeometryInfo<dim>::unit_normal_orientation[i];
// second tangential
// vector in direction
// of the (nd+2)%3 axis
tang2[(nd + 2) % dim] = 1.;

// same unit tangents
// for all quadrature
// points on this face
std::fill(unit_tangentials[i].begin(),
unit_tangentials[i].end(),
tang1);
std::fill(
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.begin(),
unit_tangentials[GeometryInfo<dim>::faces_per_cell + i]
.end(),
tang2);
}

break;
}

default:
Assert(false, ExcNotImplemented());
}
}
}
Expand Down

0 comments on commit 72c67b6

Please sign in to comment.