Nektar++
Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | List of all members
Nektar::SpatialDomains::Geometry3D Class Reference

3D geometry information More...

#include <Geometry3D.h>

Inheritance diagram for Nektar::SpatialDomains::Geometry3D:
[legend]

Public Member Functions

 Geometry3D ()
 
 Geometry3D (const int coordim)
 
 ~Geometry3D () override=default
 
- Public Member Functions inherited from Nektar::SpatialDomains::Geometry
 Geometry ()
 Default constructor. More...
 
 Geometry (int coordim)
 Constructor when supplied a coordinate dimension. More...
 
virtual ~Geometry ()=default
 
int GetCoordim () const
 Return the coordinate dimension of this object (i.e. the dimension of the space in which this object is embedded). More...
 
void SetCoordim (int coordim)
 Sets the coordinate dimension of this object (i.e. the dimension of the space in which this object is embedded). More...
 
GeomFactorsSharedPtr GetGeomFactors ()
 Get the geometric factors for this object, generating them if required. More...
 
GeomFactorsSharedPtr GetRefGeomFactors (const Array< OneD, const LibUtilities::BasisSharedPtr > &tbasis)
 
GeomFactorsSharedPtr GetMetricInfo ()
 Get the geometric factors for this object. More...
 
LibUtilities::ShapeType GetShapeType (void)
 Get the geometric shape type of this object. More...
 
int GetGlobalID (void) const
 Get the ID of this object. More...
 
void SetGlobalID (int globalid)
 Set the ID of this object. More...
 
int GetVid (int i) const
 Get the ID of vertex i of this object. More...
 
int GetEid (int i) const
 Get the ID of edge i of this object. More...
 
int GetFid (int i) const
 Get the ID of face i of this object. More...
 
int GetTid (int i) const
 Get the ID of trace i of this object. More...
 
PointGeomSharedPtr GetVertex (int i) const
 Returns vertex i of this object. More...
 
Geometry1DSharedPtr GetEdge (int i) const
 Returns edge i of this object. More...
 
Geometry2DSharedPtr GetFace (int i) const
 Returns face i of this object. More...
 
StdRegions::Orientation GetEorient (const int i) const
 Returns the orientation of edge i with respect to the ordering of edges in the standard element. More...
 
StdRegions::Orientation GetForient (const int i) const
 Returns the orientation of face i with respect to the ordering of faces in the standard element. More...
 
int GetNumVerts () const
 Get the number of vertices of this object. More...
 
int GetNumEdges () const
 Get the number of edges of this object. More...
 
int GetNumFaces () const
 Get the number of faces of this object. More...
 
int GetShapeDim () const
 Get the object's shape dimension. More...
 
StdRegions::StdExpansionSharedPtr GetXmap () const
 Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standard element to physical element. More...
 
const Array< OneD, const NekDouble > & GetCoeffs (const int i) const
 Return the coefficients of the transformation Geometry::m_xmap in coordinate direction i. More...
 
void FillGeom ()
 Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements. More...
 
std::array< NekDouble, 6 > GetBoundingBox ()
 Generates the bounding box for the element. More...
 
void ClearBoundingBox ()
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
 Determine whether an element contains a particular Cartesian coordinate \((x,y,z)\). More...
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol)
 Determine whether an element contains a particular Cartesian coordinate \((x,y,z)\). More...
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &dist)
 Determine whether an element contains a particular Cartesian coordinate \(\vec{x} = (x,y,z)\). More...
 
NekDouble GetLocCoords (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
 Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this geometry object. More...
 
NekDouble GetCoord (const int i, const Array< OneD, const NekDouble > &Lcoord)
 Given local collapsed coordinate Lcoord, return the value of physical coordinate in direction i. More...
 
int PreliminaryCheck (const Array< OneD, const NekDouble > &gloCoord)
 A fast and robust check if a given global coord is outside of a deformed element. For regular elements, this check is unnecessary. More...
 
bool MinMaxCheck (const Array< OneD, const NekDouble > &gloCoord)
 Check if given global coord is within the BoundingBox of the element. More...
 
bool ClampLocCoords (Array< OneD, NekDouble > &locCoord, NekDouble tol=std::numeric_limits< NekDouble >::epsilon())
 Clamp local coords to be within standard regions [-1, 1]^dim. More...
 
NekDouble FindDistance (const Array< OneD, const NekDouble > &xs, Array< OneD, NekDouble > &xi)
 
int GetVertexEdgeMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given vertex. More...
 
int GetVertexFaceMap (int i, int j) const
 Returns the standard element face IDs that are connected to a given vertex. More...
 
int GetEdgeFaceMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given face. More...
 
int GetEdgeNormalToFaceVert (int i, int j) const
 Returns the standard lement edge IDs that are normal to a given face vertex. More...
 
int GetDir (const int i, const int j=0) const
 Returns the element coordinate direction corresponding to a given face coordinate direction. More...
 
void Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors. More...
 
void ResetNonRecursive (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object non-recursively: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors. More...
 
void Setup ()
 
void GenGeomFactors ()
 Handles generation of geometry factors. More...
 

Static Public Attributes

static const int kDim = 3
 

Protected Member Functions

NekDouble v_GetLocCoords (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords) override
 Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this geometry object. More...
 
void NewtonIterationForLocCoord (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &ptsx, const Array< OneD, const NekDouble > &ptsy, const Array< OneD, const NekDouble > &ptsz, Array< OneD, NekDouble > &Lcoords, NekDouble &dist)
 
void NewtonIterationForLocCoord (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
 
void v_FillGeom () override
 Put all quadrature information into face/edge structure and backward transform. More...
 
NekDouble v_GetCoord (const int i, const Array< OneD, const NekDouble > &Lcoord) override
 Given local collapsed coordinate Lcoord return the value of physical coordinate in direction i. More...
 
void v_CalculateInverseIsoParam () override
 
int v_AllLeftCheck (const Array< OneD, const NekDouble > &gloCoord) override
 
int v_GetShapeDim () const override
 Get the object's shape dimension. More...
 
int v_GetNumVerts () const override
 Get the number of vertices of this object. More...
 
int v_GetNumEdges () const override
 Get the number of edges of this object. More...
 
int v_GetNumFaces () const override
 Get the number of faces of this object. More...
 
PointGeomSharedPtr v_GetVertex (int i) const override
 
Geometry1DSharedPtr v_GetEdge (int i) const override
 Returns edge i of this object. More...
 
Geometry2DSharedPtr v_GetFace (int i) const override
 Returns face i of this object. More...
 
StdRegions::Orientation v_GetEorient (const int i) const override
 Returns the orientation of edge i with respect to the ordering of edges in the standard element. More...
 
StdRegions::Orientation v_GetForient (const int i) const override
 Returns the orientation of face i with respect to the ordering of faces in the standard element. More...
 
- Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
virtual PointGeomSharedPtr v_GetVertex (int i) const =0
 
virtual Geometry1DSharedPtr v_GetEdge (int i) const
 Returns edge i of this object. More...
 
virtual Geometry2DSharedPtr v_GetFace (int i) const
 Returns face i of this object. More...
 
virtual StdRegions::Orientation v_GetEorient (const int i) const
 Returns the orientation of edge i with respect to the ordering of edges in the standard element. More...
 
virtual StdRegions::Orientation v_GetForient (const int i) const
 Returns the orientation of face i with respect to the ordering of faces in the standard element. More...
 
virtual int v_GetNumVerts () const
 Get the number of vertices of this object. More...
 
virtual int v_GetNumEdges () const
 Get the number of edges of this object. More...
 
virtual int v_GetNumFaces () const
 Get the number of faces of this object. More...
 
virtual int v_GetShapeDim () const
 Get the object's shape dimension. More...
 
virtual StdRegions::StdExpansionSharedPtr v_GetXmap () const
 Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standard element to physical element. More...
 
virtual void v_FillGeom ()
 Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements. More...
 
virtual bool v_ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &dist)
 Determine whether an element contains a particular Cartesian coordinate \(\vec{x} = (x,y,z)\). More...
 
virtual int v_AllLeftCheck (const Array< OneD, const NekDouble > &gloCoord)
 
virtual NekDouble v_GetCoord (const int i, const Array< OneD, const NekDouble > &Lcoord)
 Given local collapsed coordinate Lcoord, return the value of physical coordinate in direction i. More...
 
virtual NekDouble v_GetLocCoords (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
 Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this geometry object. More...
 
virtual NekDouble v_FindDistance (const Array< OneD, const NekDouble > &xs, Array< OneD, NekDouble > &xi)
 
virtual int v_GetVertexEdgeMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given vertex. More...
 
virtual int v_GetVertexFaceMap (int i, int j) const
 Returns the standard element face IDs that are connected to a given vertex. More...
 
virtual int v_GetEdgeFaceMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given face. More...
 
virtual int v_GetEdgeNormalToFaceVert (const int i, const int j) const
 Returns the standard lement edge IDs that are normal to a given face vertex. More...
 
virtual int v_GetDir (const int faceidx, const int facedir) const
 Returns the element coordinate direction corresponding to a given face coordinate direction. More...
 
virtual void v_Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors. More...
 
virtual void v_Setup ()
 
virtual void v_GenGeomFactors ()=0
 
void SetUpCoeffs (const int nCoeffs)
 Initialise the Geometry::m_coeffs array. More...
 
virtual void v_CalculateInverseIsoParam ()
 

Protected Attributes

PointGeomVector m_verts
 
SegGeomVector m_edges
 
Geometry2DVector m_faces
 
std::vector< StdRegions::Orientationm_eorient
 
std::vector< StdRegions::Orientationm_forient
 
int m_eid
 
bool m_ownverts
 
- Protected Attributes inherited from Nektar::SpatialDomains::Geometry
int m_coordim
 Coordinate dimension of this geometry object. More...
 
GeomFactorsSharedPtr m_geomFactors
 Geometric factors. More...
 
GeomState m_geomFactorsState
 State of the geometric factors. More...
 
StdRegions::StdExpansionSharedPtr m_xmap
 \(\chi\) mapping containing isoparametric transformation. More...
 
GeomState m_state
 Enumeration to dictate whether coefficients are filled. More...
 
bool m_setupState
 Wether or not the setup routines have been run. More...
 
GeomType m_geomType
 Type of geometry. More...
 
LibUtilities::ShapeType m_shapeType
 Type of shape. More...
 
int m_globalID
 Global ID. More...
 
Array< OneD, Array< OneD, NekDouble > > m_coeffs
 Array containing expansion coefficients of m_xmap. More...
 
Array< OneD, NekDoublem_boundingBox
 Array containing bounding box. More...
 
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
 
Array< OneD, Array< OneD, NekDouble > > m_invIsoParam
 
int m_straightEdge
 

Additional Inherited Members

- Static Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsSharedPtr ValidateRegGeomFactor (GeomFactorsSharedPtr geomFactor)
 Check to see if a geometric factor has already been created that contains the same regular information. More...
 
- Static Protected Attributes inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsVector m_regGeomFactorsManager
 

Detailed Description

3D geometry information

Definition at line 64 of file Geometry3D.h.

Constructor & Destructor Documentation

◆ Geometry3D() [1/2]

Nektar::SpatialDomains::Geometry3D::Geometry3D ( )

Definition at line 47 of file Geometry3D.cpp.

48{
49}

◆ Geometry3D() [2/2]

Nektar::SpatialDomains::Geometry3D::Geometry3D ( const int  coordim)

Definition at line 51 of file Geometry3D.cpp.

51 : Geometry(coordim)
52{
54 "Coordinate dimension should be at least 3 for a 3D geometry.");
55}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
Geometry()
Default constructor.
Definition: Geometry.cpp:50
int m_coordim
Coordinate dimension of this geometry object.
Definition: Geometry.h:188

References ASSERTL0, and Nektar::SpatialDomains::Geometry::m_coordim.

◆ ~Geometry3D()

Nektar::SpatialDomains::Geometry3D::~Geometry3D ( )
overridedefault

Member Function Documentation

◆ NewtonIterationForLocCoord() [1/2]

void Nektar::SpatialDomains::Geometry3D::NewtonIterationForLocCoord ( const Array< OneD, const NekDouble > &  coords,
Array< OneD, NekDouble > &  Lcoords 
)
protected

Definition at line 64 of file Geometry3D.cpp.

66{
67 // maximum iterations for convergence
68 const int MaxIterations = 51;
69 // |x-xp|^2 < EPSILON error tolerance
70 const NekDouble Tol = 1.e-8;
71 // |r,s| > LcoordDIV stop the search
72 const NekDouble LcoordDiv = 15.0;
73
74 Array<OneD, const NekDouble> Jac =
75 m_geomFactors->GetJac(m_xmap->GetPointsKeys());
76
77 NekDouble ScaledTol =
78 Vmath::Vsum(Jac.size(), Jac, 1) / ((NekDouble)Jac.size());
79 ScaledTol *= Tol * Tol;
80
81 NekDouble xmap, ymap, zmap, res;
82 int cnt = 0;
83 Array<OneD, NekDouble> var(8, 1.);
84 Array<OneD, NekDouble> deriv(9);
85 Array<OneD, NekDouble> tmp(3);
86 while (cnt++ < MaxIterations)
87 {
88 var[1] = Lcoords[0];
89 var[2] = Lcoords[1];
90 var[3] = Lcoords[2];
91 var[4] = Lcoords[0] * Lcoords[1];
92 var[5] = Lcoords[1] * Lcoords[2];
93 var[6] = Lcoords[0] * Lcoords[2];
94 var[7] = var[4] * Lcoords[2];
95 // calculate the global point corresponding to Lcoords
96 xmap =
97 Vmath::Dot(m_isoParameter[0].size(), var, 1, m_isoParameter[0], 1);
98 ymap =
99 Vmath::Dot(m_isoParameter[0].size(), var, 1, m_isoParameter[1], 1);
100 zmap =
101 Vmath::Dot(m_isoParameter[0].size(), var, 1, m_isoParameter[2], 1);
102
103 tmp[0] = coords[0] - xmap;
104 tmp[1] = coords[1] - ymap;
105 tmp[2] = coords[2] - zmap;
106
107 res = tmp[0] * tmp[0] + tmp[1] * tmp[1] + tmp[2] * tmp[2];
108 if (res < ScaledTol)
109 {
110 break;
111 }
112
113 // Interpolate derivative metric at Lcoords (ddx1, ddx2, ddx3)
114 deriv[0] = m_isoParameter[0][1] + m_isoParameter[0][4] * Lcoords[1];
115 deriv[1] = m_isoParameter[0][2] + m_isoParameter[0][4] * Lcoords[0];
116 deriv[2] = m_isoParameter[0][3];
117 deriv[3] = m_isoParameter[1][1] + m_isoParameter[1][4] * Lcoords[1];
118 deriv[4] = m_isoParameter[1][2] + m_isoParameter[1][4] * Lcoords[0];
119 deriv[5] = m_isoParameter[1][3];
120 deriv[6] = m_isoParameter[2][1] + m_isoParameter[2][4] * Lcoords[1];
121 deriv[7] = m_isoParameter[2][2] + m_isoParameter[2][4] * Lcoords[0];
122 deriv[8] = m_isoParameter[2][3];
123 if (m_isoParameter[0].size() >= 6)
124 {
125 deriv[1] += m_isoParameter[0][5] * Lcoords[2];
126 deriv[2] += m_isoParameter[0][5] * Lcoords[1];
127 deriv[4] += m_isoParameter[1][5] * Lcoords[2];
128 deriv[5] += m_isoParameter[1][5] * Lcoords[1];
129 deriv[7] += m_isoParameter[2][5] * Lcoords[2];
130 deriv[8] += m_isoParameter[2][5] * Lcoords[1];
131 }
132 if (m_isoParameter[0].size() >= 8)
133 {
134 deriv[0] += m_isoParameter[0][6] * Lcoords[2] +
135 m_isoParameter[0][7] * var[5];
136 deriv[1] += m_isoParameter[0][7] * var[6];
137 deriv[2] += m_isoParameter[0][6] * Lcoords[0] +
138 m_isoParameter[0][7] * var[4];
139 deriv[3] += m_isoParameter[1][6] * Lcoords[2] +
140 m_isoParameter[1][7] * var[5];
141 deriv[4] += m_isoParameter[1][7] * var[6];
142 deriv[5] += m_isoParameter[1][6] * Lcoords[0] +
143 m_isoParameter[1][7] * var[4];
144 deriv[6] += m_isoParameter[2][6] * Lcoords[2] +
145 m_isoParameter[2][7] * var[5];
146 deriv[7] += m_isoParameter[2][7] * var[6];
147 deriv[8] += m_isoParameter[2][6] * Lcoords[0] +
148 m_isoParameter[2][7] * var[4];
149 }
150 DNekMatSharedPtr mat =
152 Vmath::Vcopy(9, deriv, 1, mat->GetPtr(), 1);
153 mat->Invert();
154 Lcoords[0] += Vmath::Dot(3, mat->GetPtr(), 1, tmp, 1);
155 Lcoords[1] += Vmath::Dot(3, mat->GetPtr() + 3, 1, tmp, 1);
156 Lcoords[2] += Vmath::Dot(3, mat->GetPtr() + 6, 1, tmp, 1);
157
158 if (!(std::isfinite(Lcoords[0]) && std::isfinite(Lcoords[1]) &&
159 std::isfinite(Lcoords[2])) ||
160 fabs(Lcoords[0]) > LcoordDiv || fabs(Lcoords[1]) > LcoordDiv ||
161 fabs(Lcoords[2]) > LcoordDiv)
162 {
163 std::ostringstream ss;
164 ss << "Iteration has diverged in NewtonIterationForLocCoord in "
165 "element "
166 << GetGlobalID();
167 WARNINGL1(false, ss.str());
168 return;
169 }
170 }
171
172 if (cnt >= MaxIterations)
173 {
174 std::ostringstream ss;
175
176 ss << "Reached MaxIterations (" << MaxIterations
177 << ") in Newton iteration ";
178
179 WARNINGL1(cnt < MaxIterations, ss.str());
180 }
181}
#define WARNINGL1(condition, msg)
Definition: ErrorUtil.hpp:243
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
Definition: Geometry.h:209
int GetGlobalID(void) const
Get the ID of this object.
Definition: Geometry.h:326
StdRegions::StdExpansionSharedPtr m_xmap
mapping containing isoparametric transformation.
Definition: Geometry.h:194
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
Definition: Geometry.h:190
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75
double NekDouble
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.hpp:608
T Dot(int n, const T *w, const T *x)
dot product
Definition: Vmath.hpp:761
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Vmath::Dot(), Nektar::eFULL, Nektar::SpatialDomains::Geometry::GetGlobalID(), Nektar::SpatialDomains::Geometry::m_geomFactors, Nektar::SpatialDomains::Geometry::m_isoParameter, Nektar::SpatialDomains::Geometry::m_xmap, Vmath::Vcopy(), Vmath::Vsum(), and WARNINGL1.

◆ NewtonIterationForLocCoord() [2/2]

void Nektar::SpatialDomains::Geometry3D::NewtonIterationForLocCoord ( const Array< OneD, const NekDouble > &  coords,
const Array< OneD, const NekDouble > &  ptsx,
const Array< OneD, const NekDouble > &  ptsy,
const Array< OneD, const NekDouble > &  ptsz,
Array< OneD, NekDouble > &  Lcoords,
NekDouble dist 
)
protected

Definition at line 190 of file Geometry3D.cpp.

196{
197 // maximum iterations for convergence
198 const int MaxIterations = NekConstants::kNewtonIterations;
199 // |x-xp|^2 < EPSILON error tolerance
200 const NekDouble Tol = 1.e-8;
201 // |r,s| > LcoordDIV stop the search
202 const NekDouble LcoordDiv = 15.0;
203
204 Array<OneD, const NekDouble> Jac =
205 m_geomFactors->GetJac(m_xmap->GetPointsKeys());
206
207 NekDouble ScaledTol =
208 Vmath::Vsum(Jac.size(), Jac, 1) / ((NekDouble)Jac.size());
209 ScaledTol *= Tol;
210
211 NekDouble xmap, ymap, zmap, F1, F2, F3;
212
213 NekDouble derx_1, derx_2, derx_3, dery_1, dery_2, dery_3, derz_1, derz_2,
214 derz_3, jac;
215
216 // save intiial guess for later reference if required.
217 NekDouble init0 = Lcoords[0], init1 = Lcoords[1], init2 = Lcoords[2];
218
219 Array<OneD, NekDouble> DxD1(ptsx.size());
220 Array<OneD, NekDouble> DxD2(ptsx.size());
221 Array<OneD, NekDouble> DxD3(ptsx.size());
222 Array<OneD, NekDouble> DyD1(ptsx.size());
223 Array<OneD, NekDouble> DyD2(ptsx.size());
224 Array<OneD, NekDouble> DyD3(ptsx.size());
225 Array<OneD, NekDouble> DzD1(ptsx.size());
226 Array<OneD, NekDouble> DzD2(ptsx.size());
227 Array<OneD, NekDouble> DzD3(ptsx.size());
228
229 // Ideally this will be stored in m_geomfactors
230 m_xmap->PhysDeriv(ptsx, DxD1, DxD2, DxD3);
231 m_xmap->PhysDeriv(ptsy, DyD1, DyD2, DyD3);
232 m_xmap->PhysDeriv(ptsz, DzD1, DzD2, DzD3);
233
234 int cnt = 0;
235 Array<OneD, DNekMatSharedPtr> I(3);
236 Array<OneD, NekDouble> eta(3);
237
238 F1 = F2 = F3 = 2000; // Starting value of Function
239 NekDouble resid = sqrt(F1 * F1 + F2 * F2 + F3 * F3);
240 while (cnt++ < MaxIterations)
241 {
242 // evaluate lagrange interpolant at Lcoords
243 m_xmap->LocCoordToLocCollapsed(Lcoords, eta);
244 I[0] = m_xmap->GetBasis(0)->GetI(eta);
245 I[1] = m_xmap->GetBasis(1)->GetI(eta + 1);
246 I[2] = m_xmap->GetBasis(2)->GetI(eta + 2);
247
248 // calculate the global point `corresponding to Lcoords
249 xmap = m_xmap->PhysEvaluate(I, ptsx);
250 ymap = m_xmap->PhysEvaluate(I, ptsy);
251 zmap = m_xmap->PhysEvaluate(I, ptsz);
252
253 F1 = coords[0] - xmap;
254 F2 = coords[1] - ymap;
255 F3 = coords[2] - zmap;
256
257 if (F1 * F1 + F2 * F2 + F3 * F3 < ScaledTol)
258 {
259 resid = sqrt(F1 * F1 + F2 * F2 + F3 * F3);
260 break;
261 }
262
263 // Interpolate derivative metric at Lcoords
264 derx_1 = m_xmap->PhysEvaluate(I, DxD1);
265 derx_2 = m_xmap->PhysEvaluate(I, DxD2);
266 derx_3 = m_xmap->PhysEvaluate(I, DxD3);
267 dery_1 = m_xmap->PhysEvaluate(I, DyD1);
268 dery_2 = m_xmap->PhysEvaluate(I, DyD2);
269 dery_3 = m_xmap->PhysEvaluate(I, DyD3);
270 derz_1 = m_xmap->PhysEvaluate(I, DzD1);
271 derz_2 = m_xmap->PhysEvaluate(I, DzD2);
272 derz_3 = m_xmap->PhysEvaluate(I, DzD3);
273
274 jac = derx_1 * (dery_2 * derz_3 - dery_3 * derz_2) -
275 derx_2 * (dery_1 * derz_3 - dery_3 * derz_1) +
276 derx_3 * (dery_1 * derz_2 - dery_2 * derz_1);
277
278 // use analytical inverse of derivitives which are also similar to
279 // those of metric factors.
280 Lcoords[0] =
281 Lcoords[0] +
282 ((dery_2 * derz_3 - dery_3 * derz_2) * (coords[0] - xmap) -
283 (derx_2 * derz_3 - derx_3 * derz_2) * (coords[1] - ymap) +
284 (derx_2 * dery_3 - derx_3 * dery_2) * (coords[2] - zmap)) /
285 jac;
286
287 Lcoords[1] =
288 Lcoords[1] -
289 ((dery_1 * derz_3 - dery_3 * derz_1) * (coords[0] - xmap) -
290 (derx_1 * derz_3 - derx_3 * derz_1) * (coords[1] - ymap) +
291 (derx_1 * dery_3 - derx_3 * dery_1) * (coords[2] - zmap)) /
292 jac;
293
294 Lcoords[2] =
295 Lcoords[2] +
296 ((dery_1 * derz_2 - dery_2 * derz_1) * (coords[0] - xmap) -
297 (derx_1 * derz_2 - derx_2 * derz_1) * (coords[1] - ymap) +
298 (derx_1 * dery_2 - derx_2 * dery_1) * (coords[2] - zmap)) /
299 jac;
300
301 if (!(std::isfinite(Lcoords[0]) && std::isfinite(Lcoords[1]) &&
302 std::isfinite(Lcoords[2])))
303 {
304 dist = 1e16;
305 std::ostringstream ss;
306 ss << "nan or inf found in NewtonIterationForLocCoord in element "
307 << GetGlobalID();
308 WARNINGL1(false, ss.str());
309 return;
310 }
311 if (fabs(Lcoords[0]) > LcoordDiv || fabs(Lcoords[1]) > LcoordDiv ||
312 fabs(Lcoords[2]) > LcoordDiv)
313 {
314 break; // lcoords have diverged so stop iteration
315 }
316 }
317
318 m_xmap->LocCoordToLocCollapsed(Lcoords, eta);
319 if (ClampLocCoords(eta, 0.))
320 {
321 I[0] = m_xmap->GetBasis(0)->GetI(eta);
322 I[1] = m_xmap->GetBasis(1)->GetI(eta + 1);
323 I[2] = m_xmap->GetBasis(2)->GetI(eta + 2);
324 // calculate the global point corresponding to Lcoords
325 xmap = m_xmap->PhysEvaluate(I, ptsx);
326 ymap = m_xmap->PhysEvaluate(I, ptsy);
327 zmap = m_xmap->PhysEvaluate(I, ptsz);
328 F1 = coords[0] - xmap;
329 F2 = coords[1] - ymap;
330 F3 = coords[2] - zmap;
331 dist = sqrt(F1 * F1 + F2 * F2 + F3 * F3);
332 }
333 else
334 {
335 dist = 0.;
336 }
337
338 if (cnt >= MaxIterations)
339 {
340 Array<OneD, NekDouble> collCoords(3);
341 m_xmap->LocCoordToLocCollapsed(Lcoords, collCoords);
342
343 // if coordinate is inside element dump error!
344 if ((collCoords[0] >= -1.0 && collCoords[0] <= 1.0) &&
345 (collCoords[1] >= -1.0 && collCoords[1] <= 1.0) &&
346 (collCoords[2] >= -1.0 && collCoords[2] <= 1.0))
347 {
348 std::ostringstream ss;
349
350 ss << "Reached MaxIterations (" << MaxIterations
351 << ") in Newton iteration ";
352 ss << "Init value (" << std::setprecision(4) << init0 << ","
353 << init1 << "," << init2 << ") ";
354 ss << "Fin value (" << Lcoords[0] << "," << Lcoords[1] << ","
355 << Lcoords[2] << ") ";
356 ss << "Resid = " << resid
357 << " Tolerance = " << std::sqrt(ScaledTol);
358
359 WARNINGL1(cnt < MaxIterations, ss.str());
360 }
361 }
362}
bool ClampLocCoords(Array< OneD, NekDouble > &locCoord, NekDouble tol=std::numeric_limits< NekDouble >::epsilon())
Clamp local coords to be within standard regions [-1, 1]^dim.
Definition: Geometry.cpp:530
static const unsigned int kNewtonIterations
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References Nektar::SpatialDomains::Geometry::ClampLocCoords(), Nektar::SpatialDomains::Geometry::GetGlobalID(), Nektar::NekConstants::kNewtonIterations, Nektar::SpatialDomains::Geometry::m_geomFactors, Nektar::SpatialDomains::Geometry::m_xmap, tinysimd::sqrt(), Vmath::Vsum(), and WARNINGL1.

Referenced by v_GetLocCoords().

◆ v_AllLeftCheck()

int Nektar::SpatialDomains::Geometry3D::v_AllLeftCheck ( const Array< OneD, const NekDouble > &  gloCoord)
overrideprotectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 183 of file Geometry3D.cpp.

185{
186 // todo: check only plane surfaces
187 return 0;
188}

◆ v_CalculateInverseIsoParam()

void Nektar::SpatialDomains::Geometry3D::v_CalculateInverseIsoParam ( )
overrideprotectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 565 of file Geometry3D.cpp.

566{
567 DNekMatSharedPtr mat =
569 for (int i = 0; i < 3; ++i)
570 {
571 for (int j = 0; j < 3; ++j)
572 {
573
574 mat->SetValue(i, j, m_isoParameter[i][j + 1]);
575 }
576 }
577 mat->Invert();
578 m_invIsoParam = Array<OneD, Array<OneD, NekDouble>>(3);
579 for (int i = 0; i < 3; ++i)
580 {
581 m_invIsoParam[i] = Array<OneD, NekDouble>(3);
582 for (int j = 0; j < 3; ++j)
583 {
584 m_invIsoParam[i][j] = mat->GetValue(i, j);
585 }
586 }
587}
Array< OneD, Array< OneD, NekDouble > > m_invIsoParam
Definition: Geometry.h:210

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::eFULL, Nektar::SpatialDomains::Geometry::m_invIsoParam, and Nektar::SpatialDomains::Geometry::m_isoParameter.

Referenced by Nektar::SpatialDomains::HexGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PrismGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PyrGeom::v_GenGeomFactors(), and Nektar::SpatialDomains::TetGeom::v_GenGeomFactors().

◆ v_FillGeom()

void Nektar::SpatialDomains::Geometry3D::v_FillGeom ( )
overrideprotectedvirtual

Put all quadrature information into face/edge structure and backward transform.

Note verts, edges, and faces are listed according to anticlockwise convention but points in _coeffs have to be in array format from left to right.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 442 of file Geometry3D.cpp.

443{
444 if (m_state == ePtsFilled)
445 {
446 return;
447 }
448
449 int i, j, k;
450
451 for (i = 0; i < m_forient.size(); i++)
452 {
453 m_faces[i]->FillGeom();
454
455 int nFaceCoeffs = m_faces[i]->GetXmap()->GetNcoeffs();
456
457 Array<OneD, unsigned int> mapArray(nFaceCoeffs);
458 Array<OneD, int> signArray(nFaceCoeffs);
459
460 if (m_forient[i] < 9)
461 {
462 m_xmap->GetTraceToElementMap(
463 i, mapArray, signArray, m_forient[i],
464 m_faces[i]->GetXmap()->GetTraceNcoeffs(0),
465 m_faces[i]->GetXmap()->GetTraceNcoeffs(1));
466 }
467 else
468 {
469 m_xmap->GetTraceToElementMap(
470 i, mapArray, signArray, m_forient[i],
471 m_faces[i]->GetXmap()->GetTraceNcoeffs(1),
472 m_faces[i]->GetXmap()->GetTraceNcoeffs(0));
473 }
474
475 for (j = 0; j < m_coordim; j++)
476 {
477 const Array<OneD, const NekDouble> &coeffs =
478 m_faces[i]->GetCoeffs(j);
479
480 for (k = 0; k < nFaceCoeffs; k++)
481 {
482 NekDouble v = signArray[k] * coeffs[k];
483 m_coeffs[j][mapArray[k]] = v;
484 }
485 }
486 }
487
489}
std::vector< StdRegions::Orientation > m_forient
Definition: Geometry3D.h:81
GeomState m_state
Enumeration to dictate whether coefficients are filled.
Definition: Geometry.h:196
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
Definition: Geometry.h:206
StdRegions::StdExpansionSharedPtr GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition: Geometry.h:435
@ ePtsFilled
Geometric information has been generated.

References Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::Geometry::GetXmap(), Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, m_faces, m_forient, Nektar::SpatialDomains::Geometry::m_state, and Nektar::SpatialDomains::Geometry::m_xmap.

Referenced by Nektar::SpatialDomains::HexGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PrismGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PyrGeom::v_GenGeomFactors(), Nektar::SpatialDomains::TetGeom::v_GenGeomFactors(), and v_GetLocCoords().

◆ v_GetCoord()

NekDouble Nektar::SpatialDomains::Geometry3D::v_GetCoord ( const int  i,
const Array< OneD, const NekDouble > &  Lcoord 
)
overrideprotectedvirtual

Given local collapsed coordinate Lcoord return the value of physical coordinate in direction i.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 495 of file Geometry3D.cpp.

497{
498 ASSERTL1(m_state == ePtsFilled, "Geometry is not in physical space");
499
500 Array<OneD, NekDouble> tmp(m_xmap->GetTotPoints());
501 m_xmap->BwdTrans(m_coeffs[i], tmp);
502
503 return m_xmap->PhysEvaluate(Lcoord, tmp);
504}
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242

References ASSERTL1, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_state, and Nektar::SpatialDomains::Geometry::m_xmap.

◆ v_GetEdge()

Geometry1DSharedPtr Nektar::SpatialDomains::Geometry3D::v_GetEdge ( int  i) const
overrideprotectedvirtual

Returns edge i of this object.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 535 of file Geometry3D.cpp.

536{
537 ASSERTL2(i >= 0 && i <= m_edges.size() - 1,
538 "Edge ID must be between 0 and " +
539 std::to_string(m_edges.size() - 1));
540 return m_edges[i];
541}
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed to...
Definition: ErrorUtil.hpp:265

References ASSERTL2, and m_edges.

◆ v_GetEorient()

StdRegions::Orientation Nektar::SpatialDomains::Geometry3D::v_GetEorient ( const int  i) const
inlineoverrideprotectedvirtual

Returns the orientation of edge i with respect to the ordering of edges in the standard element.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 549 of file Geometry3D.cpp.

550{
551 ASSERTL2(i >= 0 && i <= m_edges.size() - 1,
552 "Edge ID must be between 0 and " +
553 std::to_string(m_edges.size() - 1));
554 return m_eorient[i];
555}
std::vector< StdRegions::Orientation > m_eorient
Definition: Geometry3D.h:80

References ASSERTL2, m_edges, and m_eorient.

◆ v_GetFace()

Geometry2DSharedPtr Nektar::SpatialDomains::Geometry3D::v_GetFace ( int  i) const
overrideprotectedvirtual

Returns face i of this object.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 543 of file Geometry3D.cpp.

544{
545 ASSERTL2((i >= 0) && (i <= 5), "Edge id must be between 0 and 4");
546 return m_faces[i];
547}

References ASSERTL2, and m_faces.

◆ v_GetForient()

StdRegions::Orientation Nektar::SpatialDomains::Geometry3D::v_GetForient ( const int  i) const
overrideprotectedvirtual

Returns the orientation of face i with respect to the ordering of faces in the standard element.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 557 of file Geometry3D.cpp.

558{
559 ASSERTL2(i >= 0 && i <= m_faces.size() - 1,
560 "Face ID must be between 0 and " +
561 std::to_string(m_faces.size() - 1));
562 return m_forient[i];
563}

References ASSERTL2, m_faces, and m_forient.

◆ v_GetLocCoords()

NekDouble Nektar::SpatialDomains::Geometry3D::v_GetLocCoords ( const Array< OneD, const NekDouble > &  coords,
Array< OneD, NekDouble > &  Lcoords 
)
overrideprotectedvirtual

Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this geometry object.

For curvilinear and non-affine elements (i.e. where the Jacobian varies as a function of the standard element coordinates), this is a non-linear optimisation problem that requires the use of a Newton iteration. Note therefore that this can be an expensive operation.

Note that, clearly, the provided Cartesian coordinate lie outside the element. The function therefore returns the minimum distance from some position in the element to . Lcoords will also be constrained to fit within the range \([-1,1]^d\) where \( d \) is the dimension of the element.

Parameters
coordsInput Cartesian global coordinates
LcoordsCorresponding local coordinates
Returns
Distance between obtained coordinates and provided ones.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 364 of file Geometry3D.cpp.

366{
367 NekDouble dist = std::numeric_limits<double>::max();
368 Array<OneD, NekDouble> tmpcoords(3);
369 tmpcoords[0] = coords[0];
370 tmpcoords[1] = coords[1];
371 tmpcoords[2] = coords[2];
372 if (GetMetricInfo()->GetGtype() == eRegular)
373 {
374 tmpcoords[0] -= m_isoParameter[0][0];
375 tmpcoords[1] -= m_isoParameter[1][0];
376 tmpcoords[2] -= m_isoParameter[2][0];
377 Lcoords[0] = m_invIsoParam[0][0] * tmpcoords[0] +
378 m_invIsoParam[0][1] * tmpcoords[1] +
379 m_invIsoParam[0][2] * tmpcoords[2];
380 Lcoords[1] = m_invIsoParam[1][0] * tmpcoords[0] +
381 m_invIsoParam[1][1] * tmpcoords[1] +
382 m_invIsoParam[1][2] * tmpcoords[2];
383 Lcoords[2] = m_invIsoParam[2][0] * tmpcoords[0] +
384 m_invIsoParam[2][1] * tmpcoords[1] +
385 m_invIsoParam[2][2] * tmpcoords[2];
386 }
387 else if (m_straightEdge)
388 {
389 ClampLocCoords(Lcoords, 0.);
390 NewtonIterationForLocCoord(tmpcoords, Lcoords);
391 }
392 else
393 {
394 v_FillGeom();
395 // Determine nearest point of coords to values in m_xmap
396 int npts = m_xmap->GetTotPoints();
397 Array<OneD, NekDouble> ptsx(npts), ptsy(npts), ptsz(npts);
398 Array<OneD, NekDouble> tmp1(npts), tmp2(npts);
399
400 m_xmap->BwdTrans(m_coeffs[0], ptsx);
401 m_xmap->BwdTrans(m_coeffs[1], ptsy);
402 m_xmap->BwdTrans(m_coeffs[2], ptsz);
403
404 const Array<OneD, const NekDouble> za = m_xmap->GetPoints(0);
405 const Array<OneD, const NekDouble> zb = m_xmap->GetPoints(1);
406 const Array<OneD, const NekDouble> zc = m_xmap->GetPoints(2);
407
408 // guess the first local coords based on nearest point
409 Vmath::Sadd(npts, -coords[0], ptsx, 1, tmp1, 1);
410 Vmath::Vmul(npts, tmp1, 1, tmp1, 1, tmp1, 1);
411 Vmath::Sadd(npts, -coords[1], ptsy, 1, tmp2, 1);
412 Vmath::Vvtvp(npts, tmp2, 1, tmp2, 1, tmp1, 1, tmp1, 1);
413 Vmath::Sadd(npts, -coords[2], ptsz, 1, tmp2, 1);
414 Vmath::Vvtvp(npts, tmp2, 1, tmp2, 1, tmp1, 1, tmp1, 1);
415
416 int min_i = Vmath::Imin(npts, tmp1, 1);
417
418 // Get Local coordinates
419 int qa = za.size(), qb = zb.size();
420
421 Array<OneD, NekDouble> eta(3, 0.);
422 eta[2] = zc[min_i / (qa * qb)];
423 min_i = min_i % (qa * qb);
424 eta[1] = zb[min_i / qa];
425 eta[0] = za[min_i % qa];
426 m_xmap->LocCollapsedToLocCoord(eta, Lcoords);
427
428 // Perform newton iteration to find local coordinates
429 NewtonIterationForLocCoord(coords, ptsx, ptsy, ptsz, Lcoords, dist);
430 }
431 return dist;
432}
void v_FillGeom() override
Put all quadrature information into face/edge structure and backward transform.
Definition: Geometry3D.cpp:442
void NewtonIterationForLocCoord(const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &ptsx, const Array< OneD, const NekDouble > &ptsy, const Array< OneD, const NekDouble > &ptsz, Array< OneD, NekDouble > &Lcoords, NekDouble &dist)
Definition: Geometry3D.cpp:190
GeomFactorsSharedPtr GetMetricInfo()
Get the geometric factors for this object.
Definition: Geometry.h:310
@ eRegular
Geometry is straight-sided with constant geometric factors.
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.hpp:72
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.hpp:366
int Imin(int n, const T *x, const int incx)
Return the index of the minimum element in x.
Definition: Vmath.hpp:704
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
Definition: Vmath.hpp:194

References Nektar::SpatialDomains::Geometry::ClampLocCoords(), Nektar::SpatialDomains::eRegular, Nektar::SpatialDomains::Geometry::GetMetricInfo(), Vmath::Imin(), Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_invIsoParam, Nektar::SpatialDomains::Geometry::m_isoParameter, Nektar::SpatialDomains::Geometry::m_straightEdge, Nektar::SpatialDomains::Geometry::m_xmap, NewtonIterationForLocCoord(), Vmath::Sadd(), v_FillGeom(), Vmath::Vmul(), and Vmath::Vvtvp().

◆ v_GetNumEdges()

int Nektar::SpatialDomains::Geometry3D::v_GetNumEdges ( ) const
overrideprotectedvirtual

Get the number of edges of this object.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 520 of file Geometry3D.cpp.

521{
522 return m_edges.size();
523}

References m_edges.

◆ v_GetNumFaces()

int Nektar::SpatialDomains::Geometry3D::v_GetNumFaces ( ) const
overrideprotectedvirtual

Get the number of faces of this object.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 525 of file Geometry3D.cpp.

526{
527 return m_faces.size();
528}

References m_faces.

◆ v_GetNumVerts()

int Nektar::SpatialDomains::Geometry3D::v_GetNumVerts ( ) const
overrideprotectedvirtual

Get the number of vertices of this object.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 515 of file Geometry3D.cpp.

516{
517 return m_verts.size();
518}

References m_verts.

◆ v_GetShapeDim()

int Nektar::SpatialDomains::Geometry3D::v_GetShapeDim ( ) const
overrideprotectedvirtual

Get the object's shape dimension.

For example, a segment is one dimensional and quadrilateral is two dimensional.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 510 of file Geometry3D.cpp.

511{
512 return 3;
513}

◆ v_GetVertex()

PointGeomSharedPtr Nektar::SpatialDomains::Geometry3D::v_GetVertex ( int  i) const
overrideprotectedvirtual

Implements Nektar::SpatialDomains::Geometry.

Definition at line 530 of file Geometry3D.cpp.

531{
532 return m_verts[i];
533}

References m_verts.

Member Data Documentation

◆ kDim

const int Nektar::SpatialDomains::Geometry3D::kDim = 3
static

Definition at line 74 of file Geometry3D.h.

◆ m_edges

SegGeomVector Nektar::SpatialDomains::Geometry3D::m_edges
protected

◆ m_eid

int Nektar::SpatialDomains::Geometry3D::m_eid
protected

Definition at line 82 of file Geometry3D.h.

◆ m_eorient

std::vector<StdRegions::Orientation> Nektar::SpatialDomains::Geometry3D::m_eorient
protected

◆ m_faces

Geometry2DVector Nektar::SpatialDomains::Geometry3D::m_faces
protected

◆ m_forient

std::vector<StdRegions::Orientation> Nektar::SpatialDomains::Geometry3D::m_forient
protected

◆ m_ownverts

bool Nektar::SpatialDomains::Geometry3D::m_ownverts
protected

Definition at line 83 of file Geometry3D.h.

◆ m_verts

PointGeomVector Nektar::SpatialDomains::Geometry3D::m_verts
protected