Nektar++
Public Member Functions | Protected Member Functions | List of all members
Nektar::SpatialDomains::RefRegionCylinder Class Reference

Derived class for the refinement surface region. More...

#include <RefRegionCylinder.h>

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

Public Member Functions

 RefRegionCylinder (const unsigned int coordim, NekDouble radius, std::vector< NekDouble > coord1, std::vector< NekDouble > coord2, std::vector< unsigned int > numModes, std::vector< unsigned int > numPoints)
 Constructor. More...
 
 ~RefRegionCylinder () override=default
 Destructor. More...
 
- Public Member Functions inherited from Nektar::SpatialDomains::RefRegion
 RefRegion (const unsigned int coordim, NekDouble m_radius, std::vector< NekDouble > coord1, std::vector< NekDouble > coord2, std::vector< unsigned int > numModes, std::vector< unsigned int > numPoints)
 Constructor. More...
 
virtual ~RefRegion ()=default
 Destructor. More...
 
virtual bool v_Contains (const Array< OneD, NekDouble > &coords)=0
 Pure virtual fuction. More...
 
std::vector< unsigned int > GetNumModes ()
 Get the number of modes to update expansion. More...
 
std::vector< unsigned int > GetNumPoints ()
 Get the number of quadrature points to update expansion. More...
 

Protected Member Functions

bool v_Contains (const Array< OneD, NekDouble > &coords) override
 Check if vertex is inside the surface region. More...
 

Additional Inherited Members

- Protected Attributes inherited from Nektar::SpatialDomains::RefRegion
std::vector< NekDoublem_coord1
 Coordinate 1. More...
 
std::vector< NekDoublem_coord2
 Coordinate 2. More...
 
std::vector< unsigned int > m_numModes
 Number of modes. More...
 
std::vector< unsigned int > m_numPoints
 Number of quadrature points. More...
 
NekDouble m_radius
 Radius of the surface region. More...
 
unsigned int m_coordim
 Dimension of the coordinate (space dimension) More...
 

Detailed Description

Derived class for the refinement surface region.

Definition at line 48 of file RefRegionCylinder.h.

Constructor & Destructor Documentation

◆ RefRegionCylinder()

Nektar::SpatialDomains::RefRegionCylinder::RefRegionCylinder ( const unsigned int  coordim,
NekDouble  radius,
std::vector< NekDouble coord1,
std::vector< NekDouble coord2,
std::vector< unsigned int >  numModes,
std::vector< unsigned int >  numPoints 
)

Constructor.

Definition at line 40 of file RefRegionCylinder.cpp.

46 : RefRegion(coordim, radius, coord1, coord2, numModes, numPoints)
47{
48}
RefRegion(const unsigned int coordim, NekDouble m_radius, std::vector< NekDouble > coord1, std::vector< NekDouble > coord2, std::vector< unsigned int > numModes, std::vector< unsigned int > numPoints)
Constructor.
Definition: RefRegion.cpp:45

◆ ~RefRegionCylinder()

Nektar::SpatialDomains::RefRegionCylinder::~RefRegionCylinder ( )
overridedefault

Destructor.

Member Function Documentation

◆ v_Contains()

bool Nektar::SpatialDomains::RefRegionCylinder::v_Contains ( const Array< OneD, NekDouble > &  coords)
overrideprotectedvirtual

Check if vertex is inside the surface region.

Check if vertex is inside the cylinder.

Parameters
coordscoordinates of the vertex
Returns
true or false depending on if the vertex is inside or not of the surface defined by the user.

Implements Nektar::SpatialDomains::RefRegion.

Definition at line 57 of file RefRegionCylinder.cpp.

58{
59 const size_t dim = coords.size();
60 Array<OneD, NekDouble> e(dim, 0.0); // direction: rb - ra
61 Array<OneD, NekDouble> m(dim, 0.0); // momentum: ra x rb
62 NekDouble d = 0.0; // distance
63
64 // Direction
65 e[0] = m_coord2[0] - m_coord1[0];
66 e[1] = m_coord2[1] - m_coord1[1];
67 e[2] = m_coord2[2] - m_coord1[2];
68
69 // Cross product of vectors 'coord1' and 'coord2'
70 m[0] = m_coord1[1] * m_coord2[2] - m_coord1[2] * m_coord2[1];
71 m[1] = m_coord1[2] * m_coord2[0] - m_coord1[0] * m_coord2[2];
72 m[2] = m_coord1[0] * m_coord2[1] - m_coord1[1] * m_coord2[0];
73
74 // 1. Distance of P (coords) to line AB (coord1coord2) is equal or less than
75 // R
76 // d = || e x (rp - ra) || / || e ||
77
78 // rA - rP
79 Array<OneD, NekDouble> rpa(dim, 0.0);
80 rpa[0] = coords[0] - m_coord1[0];
81 rpa[1] = coords[1] - m_coord1[1];
82 rpa[2] = coords[2] - m_coord1[2];
83
84 // || e ||
85 NekDouble e_mod = sqrt(e[0] * e[0] + e[1] * e[1] + e[2] * e[2]);
86
87 // || e x (rp - ra) ||
88 Array<OneD, NekDouble> exrpa(dim, 0.0);
89 exrpa[0] = e[1] * rpa[2] - e[2] * rpa[1];
90 exrpa[1] = e[2] * rpa[0] - e[0] * rpa[2];
91 exrpa[2] = e[0] * rpa[1] - e[1] * rpa[0];
92
93 NekDouble exrpa_mod =
94 sqrt(exrpa[0] * exrpa[0] + exrpa[1] * exrpa[1] + exrpa[2] * exrpa[2]);
95
96 d = exrpa_mod / e_mod;
97 if (d >= m_radius)
98 {
99 return false;
100 }
101
102 // 2. Closest point Q on line AB to P
103 // rq = rp + (e x (m + e x rp)) / ||e||^{2}
104
105 // (m + e x rp)
106 Array<OneD, NekDouble> mpexrp(dim, 0.0);
107 mpexrp[0] = m[0] + e[1] * coords[2] - e[2] * coords[1];
108 mpexrp[1] = m[1] + e[2] * coords[0] - e[0] * coords[2];
109 mpexrp[2] = m[2] + e[0] * coords[1] - e[1] * coords[0];
110
111 // e x (m + e x rp) = numerator
112 Array<OneD, NekDouble> numerator(dim, 0.0);
113 numerator[0] = e[1] * mpexrp[2] - e[2] * mpexrp[1];
114 numerator[1] = e[2] * mpexrp[0] - e[0] * mpexrp[2];
115 numerator[2] = e[0] * mpexrp[1] - e[1] * mpexrp[0];
116
117 // rq
118 Array<OneD, NekDouble> rq(dim, 0.0);
119 rq[0] = coords[0] + (numerator[0] / pow(e_mod, 2));
120 rq[1] = coords[1] + (numerator[1] / pow(e_mod, 2));
121 rq[2] = coords[2] + (numerator[2] / pow(e_mod, 2));
122
123 // 3. The baricentric coordinates of Q(wa,wb) such that rq = wa*ra+wb*rb
124 // wa = ||rq x rb||/||m||, wb = ||rq x ra||/||m||
125 NekDouble wa = 0.0;
126 NekDouble wb = 0.0;
127
128 // ||m||
129 NekDouble m_mod = sqrt(m[0] * m[0] + m[1] * m[1] + m[2] * m[2]);
130 // m_mod > 0: condition
131 ASSERTL0(m_mod, "The cylinder axis must not go through the origin");
132
133 // ||rq x rb||
134 Array<OneD, NekDouble> rqxrb(dim, 0.0);
135 rqxrb[0] = rq[1] * m_coord2[2] - rq[2] * m_coord2[1];
136 rqxrb[1] = rq[2] * m_coord2[0] - rq[0] * m_coord2[2];
137 rqxrb[2] = rq[0] * m_coord2[1] - rq[1] * m_coord2[0];
138 NekDouble rqxrb_mod =
139 sqrt(rqxrb[0] * rqxrb[0] + rqxrb[1] * rqxrb[1] + rqxrb[2] * rqxrb[2]);
140
141 // ||rq x ra||
142 Array<OneD, NekDouble> rqxra(dim, 0.0);
143 rqxra[0] = rq[1] * m_coord1[2] - rq[2] * m_coord1[1];
144 rqxra[1] = rq[2] * m_coord1[0] - rq[0] * m_coord1[2];
145 rqxra[2] = rq[0] * m_coord1[1] - rq[1] * m_coord1[0];
146 NekDouble rqxra_mod =
147 sqrt(rqxra[0] * rqxra[0] + rqxra[1] * rqxra[1] + rqxra[2] * rqxra[2]);
148
149 // wa
150 wa = rqxrb_mod / m_mod;
151 // wb
152 wb = rqxra_mod / m_mod;
153
154 // 4. Check that point Q between A and B by making sure the baricentric
155 // coordinates are between 0 and 1.
156 if ((wa >= 0) && (wa <= 1) && (wb >= 0) && (wb <= 1))
157 {
158 return true;
159 }
160 else
161 {
162 return false;
163 }
164}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
std::vector< NekDouble > m_coord2
Coordinate 2.
Definition: RefRegion.h:83
std::vector< NekDouble > m_coord1
Coordinate 1.
Definition: RefRegion.h:81
NekDouble m_radius
Radius of the surface region.
Definition: RefRegion.h:89
std::vector< double > d(NPUPPER *NPUPPER)
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References ASSERTL0, Nektar::UnitTests::d(), Nektar::SpatialDomains::RefRegion::m_coord1, Nektar::SpatialDomains::RefRegion::m_coord2, Nektar::SpatialDomains::RefRegion::m_radius, and tinysimd::sqrt().