Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Private Attributes | List of all members
Nektar::MultiRegions::LocTraceToTraceMap Class Reference

A helper class to deal with trace operations in the discontinuous Galerkin code. More...

#include <LocTraceToTraceMap.h>

Collaboration diagram for Nektar::MultiRegions::LocTraceToTraceMap:
Collaboration graph
[legend]

Public Member Functions

 LocTraceToTraceMap (const ExpList &locExp, const ExpListSharedPtr &trace, const Array< OneD, Array< OneD, LocalRegions::ExpansionSharedPtr > > &elmtToTrace, const std::vector< bool > &LeftAdjacents)
 Set up trace to trace mapping components. More...
 
virtual ~LocTraceToTraceMap ()
 
void Setup2D (const ExpList &locExp, const ExpListSharedPtr &trace, const Array< OneD, Array< OneD, LocalRegions::ExpansionSharedPtr > > &elmtToTrace, const std::vector< bool > &LeftAdjacents)
 Set up member variables for a two-dimensional problem. More...
 
void Setup3D (const ExpList &locExp, const ExpListSharedPtr &trace, const Array< OneD, Array< OneD, LocalRegions::ExpansionSharedPtr > > &elmtToTrace, const std::vector< bool > &LeftAdjacents)
 Set up member variables for a three-dimensional problem. More...
 
void LocTracesFromField (const Array< OneD, const NekDouble > &field, Array< OneD, NekDouble > faces)
 Gather the local traces in physical space from field using m_fieldToLocTraceMap. More...
 
void FwdLocTracesFromField (const Array< OneD, const NekDouble > &field, Array< OneD, NekDouble > faces)
 Gather the forwards-oriented local traces in physical space from field using m_fieldToLocTraceMap. More...
 
void InterpLocEdgesToTrace (const int dir, const Array< OneD, const NekDouble > &locfaces, Array< OneD, NekDouble > edges)
 Interpolate local trace edges to global trace edge point distributions where required. More...
 
void InterpLocFacesToTrace (const int dir, const Array< OneD, const NekDouble > &locfaces, Array< OneD, NekDouble > faces)
 Interpolate local faces to trace face point distributions where required. More...
 
void AddTraceCoeffsToFieldCoeffs (const Array< OneD, const NekDouble > &trace, Array< OneD, NekDouble > &field)
 Add contributions from trace coefficients to the elemental field storage. More...
 
void AddTraceCoeffsToFieldCoeffs (const int dir, const Array< OneD, const NekDouble > &race, Array< OneD, NekDouble > &field)
 Add contributions from backwards or forwards oriented trace coefficients to the elemental field storage. More...
 
int GetNFwdLocTracePts ()
 Return the number of `forward' local trace points. More...
 
int GetNLocTracePts ()
 Return the number of local trace points. More...
 

Private Attributes

int m_nFwdLocTracePts
 The number of forward trace points. A local trace element is `forward' if it is the side selected for the global trace. More...
 
int m_nLocTracePts
 The number of local trace points. More...
 
int m_nTracePts
 The number of global trace points. More...
 
Array< OneD, int > m_fieldToLocTraceMap
 A mapping from the local trace points, arranged as all forwards traces followed by backwards traces, to elemental storage. More...
 
Array< OneD, Array< OneD, int > > m_LocTraceToTraceMap
 A mapping from local trace points to the global trace. Dimension 0 holds forward traces, dimension 1 backward. More...
 
Array< OneD, Array< OneD,
InterpLocTraceToTrace > > 
m_interpTrace
 A mapping holding the type of interpolation needed for each local trace. Dimension 0 holds forward traces, dimension 1 backward. More...
 
Array< OneD, Array< OneD,
DNekMatSharedPtr > > 
m_interpTraceI0
 Interpolation matrices for either 2D edges or first coordinate of 3D face. More...
 
Array< OneD, Array< OneD,
DNekMatSharedPtr > > 
m_interpTraceI1
 Interpolation matrices for the second coordinate of 3D face, not used in 2D. More...
 
Array< OneD, Array< OneD,
TraceInterpPoints > > 
m_interpPoints
 Interpolation points key distributions for each of the local to global mappings. More...
 
Array< OneD, Array< OneD,
Array< OneD, NekDouble > > > 
m_interpEndPtI0
 Mapping to hold first coordinate direction endpoint interpolation, which can be more optimal if using Gauss-Radau distribution for triangles. More...
 
Array< OneD, Array< OneD,
Array< OneD, NekDouble > > > 
m_interpEndPtI1
 Mapping to hold second coordinate direction endpoint interpolation, which can be more optimal if using Gauss-Radau distribution for triangles. More...
 
Array< OneD, Array< OneD, int > > m_interpNfaces
 Number of edges/faces on a 2D/3D element that require interpolation. More...
 
int m_nTraceCoeffs [2]
 Number of forwards/backwards trace coefficients. More...
 
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtMap
 Mapping from forwards/backwards trace coefficients to elemental coefficient storage. More...
 
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtTrace
 Mapping from forwards/backwards trace coefficients to the position of the trace element in global storage. More...
 
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtSign
 Sign array for mapping from forwards/backwards trace coefficients to local trace storage. More...
 

Detailed Description

A helper class to deal with trace operations in the discontinuous Galerkin code.

This class sets up a number of mappings to deal with operations that take the "local trace" of an expansion list – i.e. the concatenation of all elemental facets – to the "global trace" – where the duplicate facets between connected elements have been removed.

Elements:      Local trace:              Global trace:
+----+----+    + +---+ +   + +---+ +     + +---+ + +---+ +
|    |    |    |       |   |       |     |       |       |
|    |    |    |       |   |       |     |       |       |
+----+----+    + +---+ +   + +---+ +     + +---+ + +---+ +

There are a number of mappings that are required that this class provides maps for:

These are documented in the member variables and class functions.

Definition at line 154 of file LocTraceToTraceMap.h.

Constructor & Destructor Documentation

Nektar::MultiRegions::LocTraceToTraceMap::LocTraceToTraceMap ( const ExpList locExp,
const ExpListSharedPtr trace,
const Array< OneD, Array< OneD, LocalRegions::ExpansionSharedPtr > > &  elmtToTrace,
const std::vector< bool > &  LeftAdjacents 
)

Set up trace to trace mapping components.

Parameters
locExpExpansion list of full dimension problem.
traceExpansion list of one dimension lower trace.
elmtToTraceMapping from elemental facets to trace.
leftAdjacentsVector of bools denoting forwards-oriented traces.
Todo:
Add 1D support

Definition at line 64 of file LocTraceToTraceMap.cpp.

References ASSERTL0, and Nektar::MultiRegions::ExpList::GetExp().

70 {
71  const LocalRegions::ExpansionVector &locExpVector = *(locExp.GetExp());
72 
73  // Assume that all the elements have same dimension
74  int m_expdim = locExpVector[0]->GetShapeDimension();
75 
76  // Switch between 1D, 2D and 3D
77  switch (m_expdim)
78  {
79  case 1:
80  break;
81  case 2:
82  Setup2D(locExp, trace, elmtToTrace, LeftAdjacents);
83  break;
84  case 3:
85  Setup3D(locExp, trace, elmtToTrace, LeftAdjacents);
86  break;
87  default:
88  ASSERTL0(false, "Number of dimensions greater than 3")
89  break;
90  }
91 }
void Setup2D(const ExpList &locExp, const ExpListSharedPtr &trace, const Array< OneD, Array< OneD, LocalRegions::ExpansionSharedPtr > > &elmtToTrace, const std::vector< bool > &LeftAdjacents)
Set up member variables for a two-dimensional problem.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
std::vector< ExpansionSharedPtr > ExpansionVector
Definition: Expansion.h:70
void Setup3D(const ExpList &locExp, const ExpListSharedPtr &trace, const Array< OneD, Array< OneD, LocalRegions::ExpansionSharedPtr > > &elmtToTrace, const std::vector< bool > &LeftAdjacents)
Set up member variables for a three-dimensional problem.
Nektar::MultiRegions::LocTraceToTraceMap::~LocTraceToTraceMap ( )
virtual

Definition at line 93 of file LocTraceToTraceMap.cpp.

94 {
95 }

Member Function Documentation

void Nektar::MultiRegions::LocTraceToTraceMap::AddTraceCoeffsToFieldCoeffs ( const Array< OneD, const NekDouble > &  trace,
Array< OneD, NekDouble > &  field 
)

Add contributions from trace coefficients to the elemental field storage.

Parameters
traceArray of global trace coefficients.
fieldArray containing field coefficients storage.

Definition at line 1213 of file LocTraceToTraceMap.cpp.

1215 {
1216  int nvals = m_nTraceCoeffs[0] + m_nTraceCoeffs[1];
1217  for (int i = 0; i < nvals; ++i)
1218  {
1219  field[m_traceCoeffsToElmtMap[0][i]] +=
1220  m_traceCoeffsToElmtSign[0][i] *
1221  trace[m_traceCoeffsToElmtTrace[0][i]];
1222  }
1223 }
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtTrace
Mapping from forwards/backwards trace coefficients to the position of the trace element in global sto...
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtSign
Sign array for mapping from forwards/backwards trace coefficients to local trace storage.
int m_nTraceCoeffs[2]
Number of forwards/backwards trace coefficients.
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtMap
Mapping from forwards/backwards trace coefficients to elemental coefficient storage.
void Nektar::MultiRegions::LocTraceToTraceMap::AddTraceCoeffsToFieldCoeffs ( const int  dir,
const Array< OneD, const NekDouble > &  trace,
Array< OneD, NekDouble > &  field 
)

Add contributions from backwards or forwards oriented trace coefficients to the elemental field storage.

Parameters
dirSelects forwards (0) or backwards (1) direction
traceArray of global trace coefficients.
fieldArray containing field coefficients storage.

Definition at line 1233 of file LocTraceToTraceMap.cpp.

1237 {
1238  int nvals = m_nTraceCoeffs[dir];
1239  for (int i = 0; i < nvals; ++i)
1240  {
1241  field[m_traceCoeffsToElmtMap[dir][i]] +=
1242  m_traceCoeffsToElmtSign[dir][i] *
1243  trace[m_traceCoeffsToElmtTrace[dir][i]];
1244  }
1245 }
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtTrace
Mapping from forwards/backwards trace coefficients to the position of the trace element in global sto...
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtSign
Sign array for mapping from forwards/backwards trace coefficients to local trace storage.
int m_nTraceCoeffs[2]
Number of forwards/backwards trace coefficients.
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtMap
Mapping from forwards/backwards trace coefficients to elemental coefficient storage.
void Nektar::MultiRegions::LocTraceToTraceMap::FwdLocTracesFromField ( const Array< OneD, const NekDouble > &  field,
Array< OneD, NekDouble faces 
)

Gather the forwards-oriented local traces in physical space from field using m_fieldToLocTraceMap.

Parameters
fieldSolution field in physical space
facesResulting local forwards-oriented traces.

Definition at line 866 of file LocTraceToTraceMap.cpp.

References Vmath::Gathr().

868 {
870 }
void Gathr(int n, const T *x, const int *y, T *z)
Gather vector z[i] = x[y[i]].
Definition: Vmath.cpp:630
int m_nFwdLocTracePts
The number of forward trace points. A local trace element is `forward' if it is the side selected for...
Array< OneD, int > m_fieldToLocTraceMap
A mapping from the local trace points, arranged as all forwards traces followed by backwards traces...
int Nektar::MultiRegions::LocTraceToTraceMap::GetNFwdLocTracePts ( )
inline

Return the number of `forward' local trace points.

Definition at line 212 of file LocTraceToTraceMap.h.

References m_nFwdLocTracePts.

213  {
214  return m_nFwdLocTracePts;
215  }
int m_nFwdLocTracePts
The number of forward trace points. A local trace element is `forward' if it is the side selected for...
int Nektar::MultiRegions::LocTraceToTraceMap::GetNLocTracePts ( )
inline

Return the number of local trace points.

Definition at line 220 of file LocTraceToTraceMap.h.

References m_nLocTracePts.

221  {
222  return m_nLocTracePts;
223  }
int m_nLocTracePts
The number of local trace points.
void Nektar::MultiRegions::LocTraceToTraceMap::InterpLocEdgesToTrace ( const int  dir,
const Array< OneD, const NekDouble > &  locedges,
Array< OneD, NekDouble edges 
)

Interpolate local trace edges to global trace edge point distributions where required.

Parameters
dirSelects forwards (0) or backwards (1) direction.
locfacesLocal trace edge storage.
facesGlobal trace edge storage

Definition at line 880 of file LocTraceToTraceMap.cpp.

References ASSERTL0, ASSERTL1, Vmath::Ddot(), Nektar::MultiRegions::eInterpDir0, Nektar::MultiRegions::eInterpEndPtDir0, Nektar::MultiRegions::eNoInterp, Nektar::LibUtilities::PointsKey::GetNumPoints(), Vmath::Scatr(), and Vmath::Vcopy().

884 {
885  ASSERTL1(dir < 2,
886  "option dir out of range, "
887  " dir=0 is fwd, dir=1 is bwd");
888 
889  int cnt = 0;
890  int cnt1 = 0;
891 
892  // tmp space assuming forward map is of size of trace
893  Array<OneD, NekDouble> tmp(m_nTracePts);
894 
895  for (int i = 0; i < m_interpTrace[dir].num_elements(); ++i)
896  {
897  // Check if there are edges to interpolate
898  if (m_interpNfaces[dir][i])
899  {
900  // Get to/from points
901  LibUtilities::PointsKey fromPointsKey0 =
902  m_interpPoints[dir][i].get<0>();
903  LibUtilities::PointsKey toPointsKey0 =
904  m_interpPoints[dir][i].get<2>();
905 
906  int fnp = fromPointsKey0.GetNumPoints();
907  int tnp = toPointsKey0.GetNumPoints();
908  int nedges = m_interpNfaces[dir][i];
909 
910  // Do interpolation here if required
911  switch (m_interpTrace[dir][i])
912  {
913  case eNoInterp: // Just copy
914  {
915  Vmath::Vcopy(nedges * fnp,
916  locedges.get() + cnt,
917  1,
918  tmp.get() + cnt1,
919  1);
920  }
921  break;
922  case eInterpDir0:
923  {
924  DNekMatSharedPtr I0 = m_interpTraceI0[dir][i];
925  Blas::Dgemm('N',
926  'N',
927  tnp,
928  nedges,
929  fnp,
930  1.0,
931  I0->GetPtr().get(),
932  tnp,
933  locedges.get() + cnt,
934  fnp,
935  0.0,
936  tmp.get() + cnt1,
937  tnp);
938  }
939  break;
940  case eInterpEndPtDir0:
941  {
942  Array<OneD, NekDouble> I0 = m_interpEndPtI0[dir][i];
943 
944  for (int k = 0; k < nedges; ++k)
945  {
946  Vmath::Vcopy(fnp,
947  &locedges[cnt + k * fnp],
948  1,
949  &tmp[cnt1 + k * tnp],
950  1);
951 
952  tmp[cnt1 + k * tnp + tnp - 1] = Blas::Ddot(
953  fnp, locedges.get() + cnt + k * fnp, 1, &I0[0], 1);
954  }
955  }
956  break;
957  default:
958  ASSERTL0(false,
959  "Invalid interpolation type for 2D elements");
960  break;
961  }
962 
963  cnt += nedges * fnp;
964  cnt1 += nedges * tnp;
965  }
966  }
967 
968  Vmath::Scatr(m_LocTraceToTraceMap[dir].num_elements(),
969  tmp.get(),
970  m_LocTraceToTraceMap[dir].get(),
971  edges.get());
972 }
int m_nTracePts
The number of global trace points.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
Array< OneD, Array< OneD, DNekMatSharedPtr > > m_interpTraceI0
Interpolation matrices for either 2D edges or first coordinate of 3D face.
boost::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:70
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_interpEndPtI0
Mapping to hold first coordinate direction endpoint interpolation, which can be more optimal if using...
void Scatr(int n, const T *x, const int *y, T *z)
Scatter vector z[y[i]] = x[i].
Definition: Vmath.cpp:659
T Ddot(int n, const Array< OneD, const T > &w, const int incw, const Array< OneD, const T > &x, const int incx, const Array< OneD, const int > &y, const int incy)
Definition: VmathArray.hpp:434
Array< OneD, Array< OneD, int > > m_interpNfaces
Number of edges/faces on a 2D/3D element that require interpolation.
Array< OneD, Array< OneD, int > > m_LocTraceToTraceMap
A mapping from local trace points to the global trace. Dimension 0 holds forward traces, dimension 1 backward.
Array< OneD, Array< OneD, InterpLocTraceToTrace > > m_interpTrace
A mapping holding the type of interpolation needed for each local trace. Dimension 0 holds forward tr...
Array< OneD, Array< OneD, TraceInterpPoints > > m_interpPoints
Interpolation points key distributions for each of the local to global mappings.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:218
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Nektar::MultiRegions::LocTraceToTraceMap::InterpLocFacesToTrace ( const int  dir,
const Array< OneD, const NekDouble > &  locfaces,
Array< OneD, NekDouble faces 
)

Interpolate local faces to trace face point distributions where required.

Parameters
dirSelects forwards (0) or backwards (1) direction.
locfacesLocal trace face storage.
facesGlobal trace face storage

Definition at line 982 of file LocTraceToTraceMap.cpp.

References ASSERTL1, Vmath::Ddot(), Nektar::MultiRegions::eInterpBothDirs, Nektar::MultiRegions::eInterpDir0, Nektar::MultiRegions::eInterpDir1, Nektar::MultiRegions::eInterpEndPtDir0, Nektar::MultiRegions::eInterpEndPtDir0InterpDir1, Nektar::MultiRegions::eInterpEndPtDir1, Nektar::MultiRegions::eNoInterp, Nektar::LibUtilities::PointsKey::GetNumPoints(), Vmath::Scatr(), and Vmath::Vcopy().

986 {
987  ASSERTL1(dir < 2,
988  "option dir out of range, "
989  " dir=0 is fwd, dir=1 is bwd");
990 
991  int cnt1 = 0;
992  int cnt = 0;
993 
994  // tmp space assuming forward map is of size of trace
995  Array<OneD, NekDouble> tmp(m_nTracePts);
996 
997  for (int i = 0; i < m_interpTrace[dir].num_elements(); ++i)
998  {
999  // Check if there are faces to interpolate
1000  if (m_interpNfaces[dir][i])
1001  {
1002  // Get to/from points
1003  LibUtilities::PointsKey fromPointsKey0 =
1004  m_interpPoints[dir][i].get<0>();
1005  LibUtilities::PointsKey fromPointsKey1 =
1006  m_interpPoints[dir][i].get<1>();
1007  LibUtilities::PointsKey toPointsKey0 =
1008  m_interpPoints[dir][i].get<2>();
1009  LibUtilities::PointsKey toPointsKey1 =
1010  m_interpPoints[dir][i].get<3>();
1011 
1012  int fnp0 = fromPointsKey0.GetNumPoints();
1013  int fnp1 = fromPointsKey1.GetNumPoints();
1014  int tnp0 = toPointsKey0.GetNumPoints();
1015  int tnp1 = toPointsKey1.GetNumPoints();
1016  int nfromfacepts = m_interpNfaces[dir][i] * fnp0 * fnp1;
1017  ;
1018 
1019  // Do interpolation here if required
1020  switch (m_interpTrace[dir][i])
1021  {
1022  case eNoInterp: // Just copy
1023  {
1024  Vmath::Vcopy(nfromfacepts,
1025  locfaces.get() + cnt,
1026  1,
1027  tmp.get() + cnt1,
1028  1);
1029  }
1030  break;
1031  case eInterpDir0:
1032  {
1033  DNekMatSharedPtr I0 = m_interpTraceI0[dir][i];
1034  Blas::Dgemm('N',
1035  'N',
1036  tnp0,
1037  tnp1,
1038  fnp0,
1039  1.0,
1040  I0->GetPtr().get(),
1041  tnp0,
1042  locfaces.get() + cnt,
1043  fnp0,
1044  0.0,
1045  tmp.get() + cnt1,
1046  tnp0);
1047  }
1048  break;
1049  case eInterpEndPtDir0:
1050  {
1051  int nfaces = m_interpNfaces[dir][i];
1052  for (int k = 0; k < fnp0; ++k)
1053  {
1054  Vmath::Vcopy(nfaces * fnp1,
1055  locfaces.get() + cnt + k,
1056  fnp0,
1057  tmp.get() + cnt1 + k,
1058  tnp0);
1059  }
1060  Array<OneD, NekDouble> I0 = m_interpEndPtI0[dir][i];
1061  Blas::Dgemv('T',
1062  fnp0,
1063  tnp1 * m_interpNfaces[dir][i],
1064  1.0,
1065  tmp.get() + cnt1,
1066  tnp0,
1067  I0.get(),
1068  1,
1069  0.0,
1070  tmp.get() + cnt1 + tnp0 - 1,
1071  tnp0);
1072  }
1073  break;
1074  case eInterpDir1:
1075  {
1076  DNekMatSharedPtr I1 = m_interpTraceI1[dir][i];
1077  for (int j = 0; j < m_interpNfaces[dir][i]; ++j)
1078  {
1079  Blas::Dgemm('N',
1080  'T',
1081  tnp0,
1082  tnp1,
1083  fnp1,
1084  1.0,
1085  locfaces.get() + cnt + j * fnp0 * fnp1,
1086  tnp0,
1087  I1->GetPtr().get(),
1088  tnp1,
1089  0.0,
1090  tmp.get() + cnt1 + j * tnp0 * tnp1,
1091  tnp0);
1092  }
1093  }
1094  break;
1095  case eInterpEndPtDir1:
1096  {
1097  Array<OneD, NekDouble> I1 = m_interpEndPtI1[dir][i];
1098  for (int j = 0; j < m_interpNfaces[dir][i]; ++j)
1099  {
1100  // copy all points
1101  Vmath::Vcopy(fnp0 * fnp1,
1102  locfaces.get() + cnt + j * fnp0 * fnp1,
1103  1,
1104  tmp.get() + cnt1 + j * tnp0 * tnp1,
1105  1);
1106 
1107  // interpolate end points
1108  for (int k = 0; k < tnp0; ++k)
1109  {
1110  tmp[cnt1 + k + (j + 1) * tnp0 * tnp1 - tnp0] =
1111  Blas::Ddot(fnp1,
1112  locfaces.get() + cnt +
1113  j * fnp0 * fnp1 + k,
1114  fnp0,
1115  &I1[0],
1116  1);
1117  }
1118  }
1119  }
1120  break;
1121  case eInterpBothDirs:
1122  {
1123  DNekMatSharedPtr I0 = m_interpTraceI0[dir][i];
1124  DNekMatSharedPtr I1 = m_interpTraceI1[dir][i];
1125  Array<OneD, NekDouble> wsp(m_interpNfaces[dir][i] * fnp0 *
1126  tnp1 * fnp0);
1127 
1128  for (int j = 0; j < m_interpNfaces[dir][i]; ++j)
1129  {
1130  Blas::Dgemm('N',
1131  'T',
1132  fnp0,
1133  tnp1,
1134  fnp1,
1135  1.0,
1136  locfaces.get() + cnt + j * fnp0 * fnp1,
1137  fnp0,
1138  I1->GetPtr().get(),
1139  tnp1,
1140  0.0,
1141  wsp.get() + j * fnp0 * tnp1,
1142  fnp0);
1143  }
1144  Blas::Dgemm('N',
1145  'N',
1146  tnp0,
1147  tnp1 * m_interpNfaces[dir][i],
1148  fnp0,
1149  1.0,
1150  I0->GetPtr().get(),
1151  tnp0,
1152  wsp.get(),
1153  fnp0,
1154  0.0,
1155  tmp.get() + cnt1,
1156  tnp0);
1157  }
1158  break;
1160  {
1161  DNekMatSharedPtr I1 = m_interpTraceI1[dir][i];
1162 
1163  for (int j = 0; j < m_interpNfaces[dir][i]; ++j)
1164  {
1165  Blas::Dgemm('N',
1166  'T',
1167  fnp0,
1168  tnp1,
1169  fnp1,
1170  1.0,
1171  locfaces.get() + cnt + j * fnp0 * fnp1,
1172  fnp0,
1173  I1->GetPtr().get(),
1174  tnp1,
1175  0.0,
1176  tmp.get() + cnt1 + j * tnp0 * tnp1,
1177  tnp0);
1178  }
1179 
1180  Array<OneD, NekDouble> I0 = m_interpEndPtI0[dir][i];
1181  Blas::Dgemv('T',
1182  fnp0,
1183  tnp1 * m_interpNfaces[dir][i],
1184  1.0,
1185  tmp.get() + cnt1,
1186  tnp0,
1187  I0.get(),
1188  1,
1189  0.0,
1190  tmp.get() + cnt1 + tnp0 - 1,
1191  tnp0);
1192  }
1193  break;
1194  }
1195  cnt += nfromfacepts;
1196  cnt1 += m_interpNfaces[dir][i] * tnp0 * tnp1;
1197  }
1198  }
1199 
1200  Vmath::Scatr(m_LocTraceToTraceMap[dir].num_elements(),
1201  tmp.get(),
1202  m_LocTraceToTraceMap[dir].get(),
1203  faces.get());
1204 }
int m_nTracePts
The number of global trace points.
Array< OneD, Array< OneD, DNekMatSharedPtr > > m_interpTraceI0
Interpolation matrices for either 2D edges or first coordinate of 3D face.
boost::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:70
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_interpEndPtI0
Mapping to hold first coordinate direction endpoint interpolation, which can be more optimal if using...
void Scatr(int n, const T *x, const int *y, T *z)
Scatter vector z[y[i]] = x[i].
Definition: Vmath.cpp:659
Array< OneD, Array< OneD, DNekMatSharedPtr > > m_interpTraceI1
Interpolation matrices for the second coordinate of 3D face, not used in 2D.
T Ddot(int n, const Array< OneD, const T > &w, const int incw, const Array< OneD, const T > &x, const int incx, const Array< OneD, const int > &y, const int incy)
Definition: VmathArray.hpp:434
Array< OneD, Array< OneD, int > > m_interpNfaces
Number of edges/faces on a 2D/3D element that require interpolation.
Array< OneD, Array< OneD, int > > m_LocTraceToTraceMap
A mapping from local trace points to the global trace. Dimension 0 holds forward traces, dimension 1 backward.
Array< OneD, Array< OneD, InterpLocTraceToTrace > > m_interpTrace
A mapping holding the type of interpolation needed for each local trace. Dimension 0 holds forward tr...
Array< OneD, Array< OneD, TraceInterpPoints > > m_interpPoints
Interpolation points key distributions for each of the local to global mappings.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:218
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_interpEndPtI1
Mapping to hold second coordinate direction endpoint interpolation, which can be more optimal if usin...
void Nektar::MultiRegions::LocTraceToTraceMap::LocTracesFromField ( const Array< OneD, const NekDouble > &  field,
Array< OneD, NekDouble faces 
)

Gather the local traces in physical space from field using m_fieldToLocTraceMap.

Parameters
fieldSolution field in physical space
facesResulting local traces.

Definition at line 850 of file LocTraceToTraceMap.cpp.

References Vmath::Gathr().

852 {
853  Vmath::Gathr(m_fieldToLocTraceMap.num_elements(),
854  field,
856  faces);
857 }
void Gathr(int n, const T *x, const int *y, T *z)
Gather vector z[i] = x[y[i]].
Definition: Vmath.cpp:630
Array< OneD, int > m_fieldToLocTraceMap
A mapping from the local trace points, arranged as all forwards traces followed by backwards traces...
void Nektar::MultiRegions::LocTraceToTraceMap::Setup2D ( const ExpList locExp,
const ExpListSharedPtr trace,
const Array< OneD, Array< OneD, LocalRegions::ExpansionSharedPtr > > &  elmtToTrace,
const std::vector< bool > &  LeftAdjacents 
)

Set up member variables for a two-dimensional problem.

Parameters
locExpExpansion list of 2D elements
traceExpansion list of the one-dimensional trace.
elmtToTraceMapping from elemental edges to trace.
leftAdjacentsVector of bools denoting forwards-oriented traces.

Definition at line 105 of file LocTraceToTraceMap.cpp.

References Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eGaussRadauMAlpha1Beta0, Nektar::MultiRegions::eInterpDir0, Nektar::MultiRegions::eInterpEndPtDir0, Nektar::MultiRegions::eNoInterp, Nektar::LibUtilities::eNoPointsType, Nektar::LibUtilities::eTriangle, Nektar::MultiRegions::ExpList::GetCoeff_Offset(), Nektar::MultiRegions::ExpList::GetExp(), Nektar::LibUtilities::PointsKey::GetNumPoints(), Nektar::MultiRegions::ExpList::GetPhys_Offset(), Nektar::LibUtilities::PointsKey::GetPointsType(), Nektar::iterator, Nektar::LibUtilities::PointsManager(), sign, and Vmath::Vcopy().

111 {
112  m_LocTraceToTraceMap = Array<OneD, Array<OneD, int> >(2);
113  m_interpTrace = Array<OneD, Array<OneD, InterpLocTraceToTrace> >(2);
114  m_interpTraceI0 = Array<OneD, Array<OneD, DNekMatSharedPtr> >(2);
115  m_interpEndPtI0 = Array<OneD, Array<OneD, Array<OneD, NekDouble> > >(2);
116  m_interpPoints = Array<OneD, Array<OneD, TraceInterpPoints> >(2);
117  m_interpNfaces = Array<OneD, Array<OneD, int> >(2);
118 
119  m_traceCoeffsToElmtMap = Array<OneD, Array<OneD, int> >(2);
120  m_traceCoeffsToElmtTrace = Array<OneD, Array<OneD, int> >(2);
121  m_traceCoeffsToElmtSign = Array<OneD, Array<OneD, int> >(2);
122 
124  const boost::shared_ptr<LocalRegions::ExpansionVector> exp =
125  locExp.GetExp();
126 
127  int cnt, n, e, phys_offset;
128 
129  int nexp = exp->size();
130  m_nTracePts = trace->GetTotPoints();
131 
132  // Count number of edges and points required for maps
133  int nFwdPts = 0;
134  int nBwdPts = 0;
135  int nFwdCoeffs = 0;
136  int nBwdCoeffs = 0;
137  m_nFwdLocTracePts = 0;
138  m_nLocTracePts = 0;
139 
140  for (cnt = n = 0; n < nexp; ++n)
141  {
142  exp2d = (*exp)[n]->as<LocalRegions::Expansion2D>();
143 
144  for (int i = 0; i < exp2d->GetNedges(); ++i, ++cnt)
145  {
146  int nLocPts = exp2d->GetEdgeNumPoints(i);
147  m_nLocTracePts += nLocPts;
148 
149  if (LeftAdjacents[cnt])
150  {
151  nFwdPts += elmtToTrace[n][i]->GetTotPoints();
152  nFwdCoeffs += elmtToTrace[n][i]->GetNcoeffs();
153  m_nFwdLocTracePts += nLocPts;
154  }
155  else
156  {
157  nBwdPts += elmtToTrace[n][i]->GetTotPoints();
158  nBwdCoeffs += elmtToTrace[n][i]->GetNcoeffs();
159  }
160  }
161  }
162 
163  m_fieldToLocTraceMap = Array<OneD, int>(m_nLocTracePts);
164 
165  m_LocTraceToTraceMap[0] = Array<OneD, int>(nFwdPts);
166  m_LocTraceToTraceMap[1] = Array<OneD, int>(nBwdPts);
167 
168  m_nTraceCoeffs[0] = nFwdCoeffs;
169  m_nTraceCoeffs[1] = nBwdCoeffs;
170 
171  m_traceCoeffsToElmtMap[0] = Array<OneD, int>(nFwdCoeffs + nBwdCoeffs);
172  m_traceCoeffsToElmtMap[1] = m_traceCoeffsToElmtMap[0] + nFwdCoeffs;
173  m_traceCoeffsToElmtTrace[0] = Array<OneD, int>(nFwdCoeffs + nBwdCoeffs);
174  m_traceCoeffsToElmtTrace[1] = m_traceCoeffsToElmtTrace[0] + nFwdCoeffs;
175  m_traceCoeffsToElmtSign[0] = Array<OneD, int>(nFwdCoeffs + nBwdCoeffs);
176  m_traceCoeffsToElmtSign[1] = m_traceCoeffsToElmtSign[0] + nFwdCoeffs;
177 
178  // Gather information about trace interpolations
179  map<TraceInterpPoints, vector<pair<int, int> >, cmpop> TraceInterpMap;
180  map<TraceInterpPoints, vector<pair<int, int> >, cmpop>::iterator it;
181 
182  vector<vector<int> > TraceOrder;
183  TraceOrder.resize(nexp);
184  int nedge;
185  int fwdcnt = 0;
186  int bwdcnt = 0;
187 
188  // Generate a map of similar traces with the same
189  // interpolation requirements
190  for (cnt = n = 0; n < nexp; ++n)
191  {
192  exp2d = (*exp)[n]->as<LocalRegions::Expansion2D>();
193  nedge = exp2d->GetNedges();
194  TraceOrder[n].resize(nedge);
195 
196  int coeffoffset = locExp.GetCoeff_Offset(n);
197  for (e = 0; e < nedge; ++e, ++cnt)
198  {
199  StdRegions::StdExpansionSharedPtr edge = elmtToTrace[n][e];
200  StdRegions::Orientation orient = exp2d->GetEorient(e);
201 
202  LibUtilities::PointsKey fromPointsKey0, fromPointsKey1;
203  LibUtilities::PointsKey toPointsKey0, toPointsKey1;
204 
205  // For Spencer's data structure
206  int basisDir = 0;
207  if (exp2d->DetShapeType() == LibUtilities::eTriangle)
208  {
209  basisDir = e == 0 ? 0 : 1;
210  }
211  else
212  {
213  basisDir = e % 2;
214  }
215 
216  fromPointsKey0 = exp2d->GetBasis(basisDir)->GetPointsKey();
217  fromPointsKey1 =
218  LibUtilities::PointsKey(0, LibUtilities::eNoPointsType);
219  toPointsKey0 = edge->GetBasis(0)->GetPointsKey();
220  toPointsKey1 =
221  LibUtilities::PointsKey(0, LibUtilities::eNoPointsType);
222 
223  // Spencer's data structure
224  TraceInterpPoints fpoint(
225  fromPointsKey0, fromPointsKey1, toPointsKey0, toPointsKey1);
226 
227  pair<int, int> epf(n, e);
228  TraceInterpMap[fpoint].push_back(epf);
229  TraceOrder[n][e] = cnt;
230 
231  // Setup for coefficient mapping from trace normal flux
232  // to elements
233  Array<OneD, unsigned int> map;
234  Array<OneD, int> sign;
235  // Not sure about the call GetBasisNumModes passing (0)!
236  exp2d->GetEdgeToElementMap(
237  e, orient, map, sign, edge->GetBasisNumModes(0));
238 
239  int order_f = edge->GetNcoeffs();
240  int foffset = trace->GetCoeff_Offset(edge->GetElmtId());
241 
242  double fac = (*exp)[n]->EdgeNormalNegated(e) ? -1.0 : 1.0;
243 
245  elmtToTrace[n][e]->as<LocalRegions::Expansion1D>();
246 
247  if (exp2d->GetEdgeExp(e)->GetRightAdjacentElementExp())
248  {
249  if (locExp1d->GetRightAdjacentElementExp()
250  ->GetGeom()
251  ->GetGlobalID() == exp2d->GetGeom()->GetGlobalID())
252  {
253  fac = -1.0;
254  }
255  }
256 
257  if (LeftAdjacents[cnt])
258  {
259  for (int i = 0; i < order_f; ++i)
260  {
261  m_traceCoeffsToElmtMap[0][fwdcnt] = coeffoffset + map[i];
262  m_traceCoeffsToElmtTrace[0][fwdcnt] = foffset + i;
263  m_traceCoeffsToElmtSign[0][fwdcnt++] = fac * sign[i];
264  }
265  }
266  else
267  {
268  for (int i = 0; i < order_f; ++i)
269  {
270  m_traceCoeffsToElmtMap[1][bwdcnt] = coeffoffset + map[i];
271  m_traceCoeffsToElmtTrace[1][bwdcnt] = foffset + i;
272  m_traceCoeffsToElmtSign[1][bwdcnt++] = fac * sign[i];
273  }
274  }
275  }
276  }
277 
278  int nInterpType = TraceInterpMap.size();
279 
280  for (int i = 0; i < 2; ++i)
281  {
282  m_interpTrace[i] = Array<OneD, InterpLocTraceToTrace>(nInterpType);
283  m_interpTraceI0[i] = Array<OneD, DNekMatSharedPtr>(nInterpType);
284  m_interpEndPtI0[i] = Array<OneD, Array<OneD, NekDouble> >(nInterpType);
285  m_interpPoints[i] = Array<OneD, TraceInterpPoints>(nInterpType);
286  m_interpNfaces[i] = Array<OneD, int>(nInterpType, 0);
287  }
288 
289  int nedgepts, nedgepts1;
290  int cnt1 = 0;
291  int cnt2 = 0;
292  int cntFwd = 0;
293  int cntBwd = 0;
294  int cntFwd1 = 0;
295  int cntBwd1 = 0;
296  int set;
297  Array<OneD, int> edgeids;
298  Array<OneD, int> locTraceToTraceMap;
299  cnt = 0;
300 
301  for (it = TraceInterpMap.begin(); it != TraceInterpMap.end(); ++it, ++cnt1)
302  {
303  LibUtilities::PointsKey fromPointsKey0 = it->first.get<0>();
304  LibUtilities::PointsKey toPointsKey0 = it->first.get<2>();
305 
306  bool fwdSet = false;
307  bool bwdSet = false;
308 
309  for (int f = 0; f < it->second.size(); ++f, ++cnt2)
310  {
311  n = it->second[f].first;
312  e = it->second[f].second;
313 
314  StdRegions::StdExpansionSharedPtr edge = elmtToTrace[n][e];
315 
316  exp2d = (*exp)[n]->as<LocalRegions::Expansion2D>();
317  phys_offset = locExp.GetPhys_Offset(n);
318 
319  // Mapping of new edge order to one that loops over elmts
320  // then set up mapping of faces in standard cartesian order
321  exp2d->GetEdgePhysMap(e, edgeids);
322 
323  nedgepts = exp2d->GetEdgeNumPoints(e);
324  nedgepts1 = edge->GetTotPoints();
325 
326  StdRegions::Orientation orient = exp2d->GetCartesianEorient(e);
327 
328  // Account for eBackwards orientation
329  exp2d->ReOrientEdgePhysMap(elmtToTrace[n][e]->GetNverts(),
330  orient,
331  toPointsKey0.GetNumPoints(),
332  locTraceToTraceMap);
333 
334  int offset = trace->GetPhys_Offset(elmtToTrace[n][e]->GetElmtId());
335 
336  if (LeftAdjacents[TraceOrder[n][e]])
337  {
338  for (int i = 0; i < nedgepts; ++i)
339  {
340  m_fieldToLocTraceMap[cntFwd + i] = phys_offset + edgeids[i];
341  }
342 
343  for (int i = 0; i < nedgepts1; ++i)
344  {
345  m_LocTraceToTraceMap[0][cntFwd1 + i] =
346  offset + locTraceToTraceMap[i];
347  }
348 
349  cntFwd += nedgepts;
350  cntFwd1 += nedgepts1;
351  set = 0;
352  }
353  else
354  {
355  for (int i = 0; i < nedgepts; ++i)
356  {
358  phys_offset + edgeids[i];
359  }
360 
361  for (int i = 0; i < nedgepts1; ++i)
362  {
363  m_LocTraceToTraceMap[1][cntBwd1 + i] =
364  offset + locTraceToTraceMap[i];
365  }
366 
367  cntBwd += nedgepts;
368  cntBwd1 += nedgepts1;
369  set = 1;
370  }
371 
372  m_interpNfaces[set][cnt1] += 1;
373 
374  if ((fwdSet == false && set == 0) || (bwdSet == false && set == 1))
375  {
376  m_interpPoints[set][cnt1] = it->first;
377 
378  if (fromPointsKey0 == toPointsKey0)
379  {
380  m_interpTrace[set][cnt1] = eNoInterp;
381  }
382  else
383  {
384  m_interpTrace[set][cnt1] = eInterpDir0;
385  m_interpTraceI0[set][cnt1] =
386  LibUtilities::PointsManager()[fromPointsKey0]->GetI(
387  toPointsKey0);
388 
389  // Check to see if we can
390  // just interpolate endpoint
391  if ((fromPointsKey0.GetPointsType() ==
393  (toPointsKey0.GetPointsType() ==
395  {
396  if (fromPointsKey0.GetNumPoints() + 1 ==
397  toPointsKey0.GetNumPoints())
398  {
399  m_interpTrace[set][cnt1] = eInterpEndPtDir0;
400 
401  int fnp0 = fromPointsKey0.GetNumPoints();
402 
403  int tnp0 = toPointsKey0.GetNumPoints();
404 
405  m_interpEndPtI0[set][cnt1] =
406  Array<OneD, NekDouble>(fnp0);
407 
408  Vmath::Vcopy(
409  fnp0,
410  m_interpTraceI0[set][cnt1]->GetPtr().get() +
411  tnp0 - 1,
412  tnp0,
413  &m_interpEndPtI0[set][cnt1][0],
414  1);
415  }
416  }
417  }
418 
419  if (set == 0)
420  {
421  fwdSet = true;
422  }
423  else
424  {
425  bwdSet = true;
426  }
427  }
428  }
429  }
430 }
int m_nTracePts
The number of global trace points.
boost::tuple< LibUtilities::PointsKey, LibUtilities::PointsKey, LibUtilities::PointsKey, LibUtilities::PointsKey > TraceInterpPoints
Map holding points distributions required for interpolation of local traces onto global trace in two ...
Array< OneD, Array< OneD, DNekMatSharedPtr > > m_interpTraceI0
Interpolation matrices for either 2D edges or first coordinate of 3D face.
#define sign(a, b)
return the sign(b)*a
Definition: Polylib.cpp:22
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtTrace
Mapping from forwards/backwards trace coefficients to the position of the trace element in global sto...
Gauss Radau pinned at x=-1, .
Definition: PointsType.h:57
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_interpEndPtI0
Mapping to hold first coordinate direction endpoint interpolation, which can be more optimal if using...
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtSign
Sign array for mapping from forwards/backwards trace coefficients to local trace storage.
int m_nLocTracePts
The number of local trace points.
int m_nTraceCoeffs[2]
Number of forwards/backwards trace coefficients.
int m_nFwdLocTracePts
The number of forward trace points. A local trace element is `forward' if it is the side selected for...
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtMap
Mapping from forwards/backwards trace coefficients to elemental coefficient storage.
PointsManagerT & PointsManager(void)
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
Array< OneD, Array< OneD, int > > m_interpNfaces
Number of edges/faces on a 2D/3D element that require interpolation.
Array< OneD, Array< OneD, int > > m_LocTraceToTraceMap
A mapping from local trace points to the global trace. Dimension 0 holds forward traces, dimension 1 backward.
Array< OneD, Array< OneD, InterpLocTraceToTrace > > m_interpTrace
A mapping holding the type of interpolation needed for each local trace. Dimension 0 holds forward tr...
Array< OneD, Array< OneD, TraceInterpPoints > > m_interpPoints
Interpolation points key distributions for each of the local to global mappings.
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
boost::shared_ptr< Expansion1D > Expansion1DSharedPtr
Definition: Expansion1D.h:53
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
Array< OneD, int > m_fieldToLocTraceMap
A mapping from the local trace points, arranged as all forwards traces followed by backwards traces...
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:50
boost::shared_ptr< Expansion2D > Expansion2DSharedPtr
Definition: Expansion1D.h:49
void Nektar::MultiRegions::LocTraceToTraceMap::Setup3D ( const ExpList locExp,
const ExpListSharedPtr trace,
const Array< OneD, Array< OneD, LocalRegions::ExpansionSharedPtr > > &  elmtToTrace,
const std::vector< bool > &  LeftAdjacents 
)

Set up member variables for a three-dimensional problem.

Parameters
locExpExpansion list of 3D elements
traceExpansion list of the two-dimensional trace.
elmtToTraceMapping from elemental faces to trace.
leftAdjacentsVector of bools denoting forwards-oriented traces.

Definition at line 440 of file LocTraceToTraceMap.cpp.

References Nektar::StdRegions::eDir1FwdDir2_Dir2FwdDir1, Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eGaussRadauMAlpha1Beta0, Nektar::MultiRegions::eInterpBothDirs, Nektar::MultiRegions::eInterpDir0, Nektar::MultiRegions::eInterpDir1, Nektar::MultiRegions::eInterpEndPtDir0, Nektar::MultiRegions::eInterpEndPtDir0InterpDir1, Nektar::MultiRegions::eInterpEndPtDir1, Nektar::MultiRegions::eNoInterp, Nektar::MultiRegions::ExpList::GetCoeff_Offset(), Nektar::MultiRegions::ExpList::GetExp(), Nektar::LibUtilities::PointsKey::GetNumPoints(), Nektar::MultiRegions::ExpList::GetPhys_Offset(), Nektar::LibUtilities::PointsKey::GetPointsType(), Nektar::iterator, Nektar::LibUtilities::PointsManager(), sign, and Vmath::Vcopy().

446 {
447  m_LocTraceToTraceMap = Array<OneD, Array<OneD, int> >(2);
448 
449  m_interpTrace = Array<OneD, Array<OneD, InterpLocTraceToTrace> >(2);
450  m_interpTraceI0 = Array<OneD, Array<OneD, DNekMatSharedPtr> >(2);
451  m_interpTraceI1 = Array<OneD, Array<OneD, DNekMatSharedPtr> >(2);
452  m_interpEndPtI0 = Array<OneD, Array<OneD, Array<OneD, NekDouble> > >(2);
453  m_interpEndPtI1 = Array<OneD, Array<OneD, Array<OneD, NekDouble> > >(2);
454  m_interpPoints = Array<OneD, Array<OneD, TraceInterpPoints> >(2);
455  m_interpNfaces = Array<OneD, Array<OneD, int> >(2);
456 
457  m_traceCoeffsToElmtMap = Array<OneD, Array<OneD, int> >(2);
458  m_traceCoeffsToElmtTrace = Array<OneD, Array<OneD, int> >(2);
459  m_traceCoeffsToElmtSign = Array<OneD, Array<OneD, int> >(2);
460 
462  const boost::shared_ptr<LocalRegions::ExpansionVector> exp =
463  locExp.GetExp();
464 
465  int cnt, n, e, phys_offset;
466 
467  int nexp = exp->size();
468  m_nTracePts = trace->GetTotPoints();
469 
470  // Count number of faces and points required for maps
471  int nFwdPts = 0;
472  int nBwdPts = 0;
473  int nFwdCoeffs = 0;
474  int nBwdCoeffs = 0;
475  m_nFwdLocTracePts = 0;
476  m_nLocTracePts = 0;
477 
478  for (cnt = n = 0; n < nexp; ++n)
479  {
480  exp3d = (*exp)[n]->as<LocalRegions::Expansion3D>();
481 
482  for (int i = 0; i < exp3d->GetNfaces(); ++i, ++cnt)
483  {
484  int nLocPts = exp3d->GetFaceNumPoints(i);
485  m_nLocTracePts += nLocPts;
486 
487  if (LeftAdjacents[cnt])
488  {
489  nFwdPts += elmtToTrace[n][i]->GetTotPoints();
490  nFwdCoeffs += elmtToTrace[n][i]->GetNcoeffs();
491  m_nFwdLocTracePts += nLocPts;
492  }
493  else
494  {
495  nBwdPts += elmtToTrace[n][i]->GetTotPoints();
496  nBwdCoeffs += elmtToTrace[n][i]->GetNcoeffs();
497  }
498  }
499  }
500 
501  m_fieldToLocTraceMap = Array<OneD, int>(m_nLocTracePts);
502 
503  m_LocTraceToTraceMap[0] = Array<OneD, int>(nFwdPts);
504  m_LocTraceToTraceMap[1] = Array<OneD, int>(nBwdPts);
505 
506  m_nTraceCoeffs[0] = nFwdCoeffs;
507  m_nTraceCoeffs[1] = nBwdCoeffs;
508 
509  m_traceCoeffsToElmtMap[0] = Array<OneD, int>(nFwdCoeffs + nBwdCoeffs);
510  m_traceCoeffsToElmtMap[1] = m_traceCoeffsToElmtMap[0] + nFwdCoeffs;
511  m_traceCoeffsToElmtTrace[0] = Array<OneD, int>(nFwdCoeffs + nBwdCoeffs);
512  m_traceCoeffsToElmtTrace[1] = m_traceCoeffsToElmtTrace[0] + nFwdCoeffs;
513  m_traceCoeffsToElmtSign[0] = Array<OneD, int>(nFwdCoeffs + nBwdCoeffs);
514  m_traceCoeffsToElmtSign[1] = m_traceCoeffsToElmtSign[0] + nFwdCoeffs;
515 
516  // Gather information about trace interpolations
517  map<TraceInterpPoints, vector<pair<int, int> >, cmpop> TraceInterpMap;
518  map<TraceInterpPoints, vector<pair<int, int> >, cmpop>::iterator it;
519 
520  vector<vector<int> > TraceOrder;
521  TraceOrder.resize(nexp);
522  int nface;
523  int fwdcnt = 0;
524  int bwdcnt = 0;
525 
526  // Generate a map of similar traces with the same
527  // interpolation or projection requirements
528  for (cnt = n = 0; n < nexp; ++n)
529  {
530  exp3d = (*exp)[n]->as<LocalRegions::Expansion3D>();
531  nface = exp3d->GetNfaces();
532  TraceOrder[n].resize(nface);
533 
534  int coeffoffset = locExp.GetCoeff_Offset(n);
535  for (e = 0; e < nface; ++e, ++cnt)
536  {
537  StdRegions::StdExpansionSharedPtr face = elmtToTrace[n][e];
538  StdRegions::Orientation orient = exp3d->GetForient(e);
539 
540  LibUtilities::PointsKey fromPointsKey0, fromPointsKey1;
541  LibUtilities::PointsKey toPointsKey0, toPointsKey1;
542 
543  // 3D specific
544  int dir0 = exp3d->GetGeom3D()->GetDir(e, 0);
545  int dir1 = exp3d->GetGeom3D()->GetDir(e, 1);
546 
547  fromPointsKey0 = exp3d->GetBasis(dir0)->GetPointsKey();
548  fromPointsKey1 = exp3d->GetBasis(dir1)->GetPointsKey();
549 
551  {
552  toPointsKey0 = face->GetBasis(0)->GetPointsKey();
553  toPointsKey1 = face->GetBasis(1)->GetPointsKey();
554  }
555  else // transpose points key evaluation
556  {
557  toPointsKey0 = face->GetBasis(1)->GetPointsKey();
558  toPointsKey1 = face->GetBasis(0)->GetPointsKey();
559  }
560 
561  TraceInterpPoints fpoint(
562  fromPointsKey0, fromPointsKey1, toPointsKey0, toPointsKey1);
563 
564  pair<int, int> epf(n, e);
565  TraceInterpMap[fpoint].push_back(epf);
566  TraceOrder[n][e] = cnt;
567 
568  // Setup for coefficient mapping from trace normal
569  // flux to elements
570  Array<OneD, unsigned int> map;
571  Array<OneD, int> sign;
572  exp3d->GetFaceToElementMap(e,
573  orient,
574  map,
575  sign,
576  face->GetBasisNumModes(0),
577  face->GetBasisNumModes(1));
578 
579  int order_f = face->GetNcoeffs();
580  int foffset = trace->GetCoeff_Offset(face->GetElmtId());
581 
582  int fac = (*exp)[n]->FaceNormalNegated(e) ? -1.0 : 1.0;
583 
584  if (exp3d->GetFaceExp(e)->GetRightAdjacentElementExp())
585  {
586  if (exp3d->GetFaceExp(e)
587  ->GetRightAdjacentElementExp()
588  ->GetGeom3D()
589  ->GetGlobalID() == exp3d->GetGeom3D()->GetGlobalID())
590  {
591  fac = -1.0;
592  }
593  }
594 
595  if (LeftAdjacents[cnt])
596  {
597  for (int i = 0; i < order_f; ++i)
598  {
599  m_traceCoeffsToElmtMap[0][fwdcnt] = coeffoffset + map[i];
600  m_traceCoeffsToElmtTrace[0][fwdcnt] = foffset + i;
601  m_traceCoeffsToElmtSign[0][fwdcnt++] = fac * sign[i];
602  }
603  }
604  else
605  {
606  for (int i = 0; i < order_f; ++i)
607  {
608  m_traceCoeffsToElmtMap[1][bwdcnt] = coeffoffset + map[i];
609  m_traceCoeffsToElmtTrace[1][bwdcnt] = foffset + i;
610  m_traceCoeffsToElmtSign[1][bwdcnt++] = fac * sign[i];
611  }
612  }
613  }
614  }
615 
616  int nInterpType = TraceInterpMap.size();
617  for (int i = 0; i < 2; ++i)
618  {
619  m_interpTrace[i] = Array<OneD, InterpLocTraceToTrace>(nInterpType);
620  m_interpTraceI0[i] = Array<OneD, DNekMatSharedPtr>(nInterpType);
621  m_interpTraceI1[i] = Array<OneD, DNekMatSharedPtr>(nInterpType);
622  m_interpEndPtI0[i] = Array<OneD, Array<OneD, NekDouble> >(nInterpType);
623  m_interpEndPtI1[i] = Array<OneD, Array<OneD, NekDouble> >(nInterpType);
624  m_interpPoints[i] = Array<OneD, TraceInterpPoints>(nInterpType);
625  m_interpNfaces[i] = Array<OneD, int>(nInterpType, 0);
626  }
627 
628  int nfacepts, nfacepts1;
629  int cnt1 = 0;
630  int cnt2 = 0;
631  int cntFwd = 0;
632  int cntBwd = 0;
633  int cntFwd1 = 0;
634  int cntBwd1 = 0;
635  int set;
636  Array<OneD, int> faceids;
637  Array<OneD, int> locTraceToTraceMap;
638  cnt = 0;
639  for (it = TraceInterpMap.begin(); it != TraceInterpMap.end(); ++it, ++cnt1)
640  {
641  LibUtilities::PointsKey fromPointsKey0 = it->first.get<0>();
642  LibUtilities::PointsKey fromPointsKey1 = it->first.get<1>();
643  LibUtilities::PointsKey toPointsKey0 = it->first.get<2>();
644  LibUtilities::PointsKey toPointsKey1 = it->first.get<3>();
645 
646  bool fwdSet = false;
647  bool bwdSet = false;
648 
649  for (int f = 0; f < it->second.size(); ++f, ++cnt2)
650  {
651  n = it->second[f].first;
652  e = it->second[f].second;
653 
654  StdRegions::StdExpansionSharedPtr face = elmtToTrace[n][e];
655 
656  exp3d = (*exp)[n]->as<LocalRegions::Expansion3D>();
657  phys_offset = locExp.GetPhys_Offset(n);
658 
659  // mapping of new face order to one that loops over elmts
660  // then faces set up mapping of faces in standard cartesian
661  // order
662  exp3d->GetFacePhysMap(e, faceids);
663  nfacepts = exp3d->GetFaceNumPoints(e);
664  nfacepts1 = face->GetTotPoints();
665 
666  StdRegions::Orientation orient = exp3d->GetForient(e);
667 
668  exp3d->ReOrientFacePhysMap(elmtToTrace[n][e]->GetNverts(),
669  orient,
670  toPointsKey0.GetNumPoints(),
671  toPointsKey1.GetNumPoints(),
672  locTraceToTraceMap);
673 
674  int offset = trace->GetPhys_Offset(elmtToTrace[n][e]->GetElmtId());
675 
676  if (LeftAdjacents[TraceOrder[n][e]])
677  {
678  for (int i = 0; i < nfacepts; ++i)
679  {
680  m_fieldToLocTraceMap[cntFwd + i] = phys_offset + faceids[i];
681  }
682 
683  for (int i = 0; i < nfacepts1; ++i)
684  {
685  m_LocTraceToTraceMap[0][cntFwd1 + i] =
686  offset + locTraceToTraceMap[i];
687  }
688 
689  cntFwd += nfacepts;
690  cntFwd1 += nfacepts1;
691  set = 0;
692  }
693  else
694  {
695  for (int i = 0; i < nfacepts; ++i)
696  {
698  phys_offset + faceids[i];
699  }
700 
701  for (int i = 0; i < nfacepts1; ++i)
702  {
703  m_LocTraceToTraceMap[1][cntBwd1 + i] =
704  offset + locTraceToTraceMap[i];
705  }
706 
707  cntBwd += nfacepts;
708  cntBwd1 += nfacepts1;
709  set = 1;
710  }
711 
712  m_interpNfaces[set][cnt1] += 1;
713 
714  if (((fwdSet == false) && (set == 0)) ||
715  ((bwdSet == false) && (set == 1)))
716  {
717  m_interpPoints[set][cnt1] = it->first;
718 
719  if (fromPointsKey0 == toPointsKey0)
720  {
721  if (fromPointsKey1 == toPointsKey1)
722  {
723  m_interpTrace[set][cnt1] = eNoInterp;
724  }
725  else
726  {
727  m_interpTrace[set][cnt1] = eInterpDir1;
728  m_interpTraceI1[set][cnt1] =
729  LibUtilities::PointsManager()[fromPointsKey1]->GetI(
730  toPointsKey1);
731 
732  // Check to see if we can just
733  // interpolate endpoint
734  if ((fromPointsKey1.GetPointsType() ==
736  (toPointsKey1.GetPointsType() ==
738  {
739  if (fromPointsKey1.GetNumPoints() + 1 ==
740  toPointsKey1.GetNumPoints())
741  {
742  m_interpTrace[set][cnt1] = eInterpEndPtDir1;
743  int fnp1 = fromPointsKey1.GetNumPoints();
744  int tnp1 = toPointsKey1.GetNumPoints();
745  m_interpEndPtI1[set][cnt1] =
746  Array<OneD, NekDouble>(fnp1);
747  Vmath::Vcopy(
748  fnp1,
749  m_interpTraceI1[set][cnt1]->GetPtr().get() +
750  tnp1 - 1,
751  tnp1,
752  &m_interpEndPtI1[set][cnt1][0],
753  1);
754  }
755  }
756  }
757  }
758  else
759  {
760  if (fromPointsKey1 == toPointsKey1)
761  {
762  m_interpTrace[set][cnt1] = eInterpDir0;
763  m_interpTraceI0[set][cnt1] =
764  LibUtilities::PointsManager()[fromPointsKey0]->GetI(
765  toPointsKey0);
766 
767  // Check to see if we can just
768  // interpolate endpoint
769  if ((fromPointsKey0.GetPointsType() ==
771  (toPointsKey0.GetPointsType() ==
773  {
774  if (fromPointsKey0.GetNumPoints() + 1 ==
775  toPointsKey0.GetNumPoints())
776  {
777  m_interpTrace[set][cnt1] = eInterpEndPtDir0;
778  int fnp0 = fromPointsKey0.GetNumPoints();
779  int tnp0 = toPointsKey0.GetNumPoints();
780  m_interpEndPtI0[set][cnt1] =
781  Array<OneD, NekDouble>(fnp0);
782  Vmath::Vcopy(
783  fnp0,
784  m_interpTraceI0[set][cnt1]->GetPtr().get() +
785  tnp0 - 1,
786  tnp0,
787  &m_interpEndPtI0[set][cnt1][0],
788  1);
789  }
790  }
791  }
792  else
793  {
794  m_interpTrace[set][cnt1] = eInterpBothDirs;
795  m_interpTraceI0[set][cnt1] =
796  LibUtilities::PointsManager()[fromPointsKey0]->GetI(
797  toPointsKey0);
798  m_interpTraceI1[set][cnt1] =
799  LibUtilities::PointsManager()[fromPointsKey1]->GetI(
800  toPointsKey1);
801 
802  // check to see if we can just
803  // interpolate endpoint
804  if ((fromPointsKey0.GetPointsType() ==
806  (toPointsKey0.GetPointsType() ==
808  {
809  if (fromPointsKey0.GetNumPoints() + 1 ==
810  toPointsKey0.GetNumPoints())
811  {
812  m_interpTrace[set][cnt1] =
814  int fnp0 = fromPointsKey0.GetNumPoints();
815  int tnp0 = toPointsKey0.GetNumPoints();
816  m_interpEndPtI0[set][cnt1] =
817  Array<OneD, NekDouble>(fnp0);
818  Vmath::Vcopy(
819  fnp0,
820  m_interpTraceI0[set][cnt1]->GetPtr().get() +
821  tnp0 - 1,
822  tnp0,
823  &m_interpEndPtI0[set][cnt1][0],
824  1);
825  }
826  }
827  }
828  }
829 
830  if (set == 0)
831  {
832  fwdSet = true;
833  }
834  else
835  {
836  bwdSet = true;
837  }
838  }
839  }
840  }
841 }
int m_nTracePts
The number of global trace points.
boost::tuple< LibUtilities::PointsKey, LibUtilities::PointsKey, LibUtilities::PointsKey, LibUtilities::PointsKey > TraceInterpPoints
Map holding points distributions required for interpolation of local traces onto global trace in two ...
Array< OneD, Array< OneD, DNekMatSharedPtr > > m_interpTraceI0
Interpolation matrices for either 2D edges or first coordinate of 3D face.
#define sign(a, b)
return the sign(b)*a
Definition: Polylib.cpp:22
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtTrace
Mapping from forwards/backwards trace coefficients to the position of the trace element in global sto...
Gauss Radau pinned at x=-1, .
Definition: PointsType.h:57
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_interpEndPtI0
Mapping to hold first coordinate direction endpoint interpolation, which can be more optimal if using...
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtSign
Sign array for mapping from forwards/backwards trace coefficients to local trace storage.
int m_nLocTracePts
The number of local trace points.
int m_nTraceCoeffs[2]
Number of forwards/backwards trace coefficients.
boost::shared_ptr< Expansion3D > Expansion3DSharedPtr
Definition: Expansion2D.h:48
int m_nFwdLocTracePts
The number of forward trace points. A local trace element is `forward' if it is the side selected for...
Array< OneD, Array< OneD, int > > m_traceCoeffsToElmtMap
Mapping from forwards/backwards trace coefficients to elemental coefficient storage.
PointsManagerT & PointsManager(void)
Array< OneD, Array< OneD, DNekMatSharedPtr > > m_interpTraceI1
Interpolation matrices for the second coordinate of 3D face, not used in 2D.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
Array< OneD, Array< OneD, int > > m_interpNfaces
Number of edges/faces on a 2D/3D element that require interpolation.
Array< OneD, Array< OneD, int > > m_LocTraceToTraceMap
A mapping from local trace points to the global trace. Dimension 0 holds forward traces, dimension 1 backward.
Array< OneD, Array< OneD, InterpLocTraceToTrace > > m_interpTrace
A mapping holding the type of interpolation needed for each local trace. Dimension 0 holds forward tr...
Array< OneD, Array< OneD, TraceInterpPoints > > m_interpPoints
Interpolation points key distributions for each of the local to global mappings.
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
Array< OneD, int > m_fieldToLocTraceMap
A mapping from the local trace points, arranged as all forwards traces followed by backwards traces...
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:50
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_interpEndPtI1
Mapping to hold second coordinate direction endpoint interpolation, which can be more optimal if usin...

Member Data Documentation

Array<OneD, int> Nektar::MultiRegions::LocTraceToTraceMap::m_fieldToLocTraceMap
private

A mapping from the local trace points, arranged as all forwards traces followed by backwards traces, to elemental storage.

Definition at line 235 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::MultiRegions::LocTraceToTraceMap::m_interpEndPtI0
private

Mapping to hold first coordinate direction endpoint interpolation, which can be more optimal if using Gauss-Radau distribution for triangles.

Definition at line 253 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::MultiRegions::LocTraceToTraceMap::m_interpEndPtI1
private

Mapping to hold second coordinate direction endpoint interpolation, which can be more optimal if using Gauss-Radau distribution for triangles.

Definition at line 257 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, int> > Nektar::MultiRegions::LocTraceToTraceMap::m_interpNfaces
private

Number of edges/faces on a 2D/3D element that require interpolation.

Definition at line 259 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, TraceInterpPoints> > Nektar::MultiRegions::LocTraceToTraceMap::m_interpPoints
private

Interpolation points key distributions for each of the local to global mappings.

Definition at line 250 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, InterpLocTraceToTrace> > Nektar::MultiRegions::LocTraceToTraceMap::m_interpTrace
private

A mapping holding the type of interpolation needed for each local trace. Dimension 0 holds forward traces, dimension 1 backward.

Definition at line 241 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, DNekMatSharedPtr> > Nektar::MultiRegions::LocTraceToTraceMap::m_interpTraceI0
private

Interpolation matrices for either 2D edges or first coordinate of 3D face.

Definition at line 244 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, DNekMatSharedPtr> > Nektar::MultiRegions::LocTraceToTraceMap::m_interpTraceI1
private

Interpolation matrices for the second coordinate of 3D face, not used in 2D.

Definition at line 247 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, int> > Nektar::MultiRegions::LocTraceToTraceMap::m_LocTraceToTraceMap
private

A mapping from local trace points to the global trace. Dimension 0 holds forward traces, dimension 1 backward.

Definition at line 238 of file LocTraceToTraceMap.h.

int Nektar::MultiRegions::LocTraceToTraceMap::m_nFwdLocTracePts
private

The number of forward trace points. A local trace element is `forward' if it is the side selected for the global trace.

Definition at line 228 of file LocTraceToTraceMap.h.

Referenced by GetNFwdLocTracePts().

int Nektar::MultiRegions::LocTraceToTraceMap::m_nLocTracePts
private

The number of local trace points.

Definition at line 230 of file LocTraceToTraceMap.h.

Referenced by GetNLocTracePts().

int Nektar::MultiRegions::LocTraceToTraceMap::m_nTraceCoeffs[2]
private

Number of forwards/backwards trace coefficients.

Definition at line 261 of file LocTraceToTraceMap.h.

int Nektar::MultiRegions::LocTraceToTraceMap::m_nTracePts
private

The number of global trace points.

Definition at line 232 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, int> > Nektar::MultiRegions::LocTraceToTraceMap::m_traceCoeffsToElmtMap
private

Mapping from forwards/backwards trace coefficients to elemental coefficient storage.

Definition at line 264 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, int> > Nektar::MultiRegions::LocTraceToTraceMap::m_traceCoeffsToElmtSign
private

Sign array for mapping from forwards/backwards trace coefficients to local trace storage.

Definition at line 270 of file LocTraceToTraceMap.h.

Array<OneD, Array<OneD, int> > Nektar::MultiRegions::LocTraceToTraceMap::m_traceCoeffsToElmtTrace
private

Mapping from forwards/backwards trace coefficients to the position of the trace element in global storage.

Definition at line 267 of file LocTraceToTraceMap.h.