Commit f4c76885 authored by Dave Moxey's avatar Dave Moxey
Browse files

Tidying 3D geometry objects, removing deprecated constructors (ticket #7).

Added helmholtz3D-prism.xml to MultiRegions demos.


git-svn-id: https://gforge.sci.utah.edu/svn/nektar/trunk@3613 305cdda6-5ce1-45b3-a98d-dfc68c8b3305
parent 800484f6
<?xml version="1.0" encoding="utf-8" ?>
<NEKTAR>
<GEOMETRY DIM="3" SPACE="3">
<VERTEX>
<V ID="0">0.00000000e+00 0.00000000e+00 5.00000000e-01</V>
<V ID="1">5.00000000e-01 0.00000000e+00 5.00000000e-01</V>
<V ID="2">5.00000000e-01 0.00000000e+00 0.00000000e+00</V>
<V ID="3">0.00000000e+00 0.00000000e+00 0.00000000e+00</V>
<V ID="4">5.00000000e-01 5.00000000e-01 5.00000000e-01</V>
<V ID="5">5.00000000e-01 5.00000000e-01 0.00000000e+00</V>
<V ID="6">0.00000000e+00 0.00000000e+00 1.00000000e+00</V>
<V ID="7">5.00000000e-01 0.00000000e+00 1.00000000e+00</V>
<V ID="8">5.00000000e-01 5.00000000e-01 1.00000000e+00</V>
<V ID="9">0.00000000e+00 5.00000000e-01 5.00000000e-01</V>
<V ID="10">0.00000000e+00 5.00000000e-01 0.00000000e+00</V>
<V ID="11">0.00000000e+00 5.00000000e-01 1.00000000e+00</V>
<V ID="12">5.00000000e-01 1.00000000e+00 5.00000000e-01</V>
<V ID="13">5.00000000e-01 1.00000000e+00 0.00000000e+00</V>
<V ID="14">5.00000000e-01 1.00000000e+00 1.00000000e+00</V>
<V ID="15">0.00000000e+00 1.00000000e+00 5.00000000e-01</V>
<V ID="16">0.00000000e+00 1.00000000e+00 0.00000000e+00</V>
<V ID="17">0.00000000e+00 1.00000000e+00 1.00000000e+00</V>
<V ID="18">1.00000000e+00 0.00000000e+00 5.00000000e-01</V>
<V ID="19">1.00000000e+00 0.00000000e+00 0.00000000e+00</V>
<V ID="20">1.00000000e+00 5.00000000e-01 5.00000000e-01</V>
<V ID="21">1.00000000e+00 5.00000000e-01 0.00000000e+00</V>
<V ID="22">1.00000000e+00 0.00000000e+00 1.00000000e+00</V>
<V ID="23">1.00000000e+00 5.00000000e-01 1.00000000e+00</V>
<V ID="24">1.00000000e+00 1.00000000e+00 5.00000000e-01</V>
<V ID="25">1.00000000e+00 1.00000000e+00 0.00000000e+00</V>
<V ID="26">1.00000000e+00 1.00000000e+00 1.00000000e+00</V>
</VERTEX>
<EDGE>
<E ID="0"> 0 1 </E>
<E ID="1"> 0 3 </E>
<E ID="2"> 0 4 </E>
<E ID="3"> 1 2 </E>
<E ID="4"> 1 4 </E>
<E ID="5"> 2 5 </E>
<E ID="6"> 3 2 </E>
<E ID="7"> 3 5 </E>
<E ID="8"> 4 5 </E>
<E ID="9"> 6 7 </E>
<E ID="10"> 6 0 </E>
<E ID="11"> 6 8 </E>
<E ID="12"> 7 1 </E>
<E ID="13"> 7 8 </E>
<E ID="14"> 8 4 </E>
<E ID="15"> 0 9 </E>
<E ID="16"> 4 9 </E>
<E ID="17"> 5 10 </E>
<E ID="18"> 3 10 </E>
<E ID="19"> 9 10 </E>
<E ID="20"> 6 11 </E>
<E ID="21"> 8 11 </E>
<E ID="22"> 11 9 </E>
<E ID="23"> 9 12 </E>
<E ID="24"> 4 12 </E>
<E ID="25"> 5 13 </E>
<E ID="26"> 10 13 </E>
<E ID="27"> 12 13 </E>
<E ID="28"> 11 14 </E>
<E ID="29"> 8 14 </E>
<E ID="30"> 14 12 </E>
<E ID="31"> 9 15 </E>
<E ID="32"> 12 15 </E>
<E ID="33"> 13 16 </E>
<E ID="34"> 10 16 </E>
<E ID="35"> 15 16 </E>
<E ID="36"> 11 17 </E>
<E ID="37"> 14 17 </E>
<E ID="38"> 17 15 </E>
<E ID="39"> 1 18 </E>
<E ID="40"> 1 20 </E>
<E ID="41"> 18 19 </E>
<E ID="42"> 18 20 </E>
<E ID="43"> 19 21 </E>
<E ID="44"> 2 19 </E>
<E ID="45"> 2 21 </E>
<E ID="46"> 20 21 </E>
<E ID="47"> 7 22 </E>
<E ID="48"> 7 23 </E>
<E ID="49"> 22 18 </E>
<E ID="50"> 22 23 </E>
<E ID="51"> 23 20 </E>
<E ID="52"> 20 4 </E>
<E ID="53"> 21 5 </E>
<E ID="54"> 23 8 </E>
<E ID="55"> 4 24 </E>
<E ID="56"> 20 24 </E>
<E ID="57"> 21 25 </E>
<E ID="58"> 5 25 </E>
<E ID="59"> 24 25 </E>
<E ID="60"> 8 26 </E>
<E ID="61"> 23 26 </E>
<E ID="62"> 26 24 </E>
<E ID="63"> 24 12 </E>
<E ID="64"> 25 13 </E>
<E ID="65"> 26 14 </E>
</EDGE>
<FACE>
<Q ID="0"> 0 3 6 1</Q>
<T ID="1"> 0 4 2</T>
<Q ID="2"> 3 5 8 4</Q>
<T ID="3"> 6 5 7</T>
<Q ID="4"> 1 7 8 2</Q>
<Q ID="5"> 9 12 0 10</Q>
<T ID="6"> 9 13 11</T>
<Q ID="7"> 12 4 14 13</Q>
<Q ID="8"> 10 2 14 11</Q>
<Q ID="9"> 15 1 18 19</Q>
<T ID="10"> 15 2 16</T>
<T ID="11"> 18 7 17</T>
<Q ID="12"> 19 17 8 16</Q>
<Q ID="13"> 20 10 15 22</Q>
<T ID="14"> 20 11 21</T>
<Q ID="15"> 22 16 14 21</Q>
<Q ID="16"> 23 19 26 27</Q>
<T ID="17"> 23 16 24</T>
<T ID="18"> 26 17 25</T>
<Q ID="19"> 27 25 8 24</Q>
<Q ID="20"> 28 22 23 30</Q>
<T ID="21"> 28 21 29</T>
<Q ID="22"> 30 24 14 29</Q>
<Q ID="23"> 32 35 33 27</Q>
<T ID="24"> 32 31 23</T>
<Q ID="25"> 35 34 19 31</Q>
<T ID="26"> 33 34 26</T>
<Q ID="27"> 37 38 32 30</Q>
<T ID="28"> 37 36 28</T>
<Q ID="29"> 38 31 22 36</Q>
<Q ID="30"> 39 41 44 3</Q>
<T ID="31"> 39 42 40</T>
<Q ID="32"> 41 43 46 42</Q>
<T ID="33"> 44 43 45</T>
<Q ID="34"> 3 45 46 40</Q>
<Q ID="35"> 47 49 39 12</Q>
<T ID="36"> 47 50 48</T>
<Q ID="37"> 49 42 51 50</Q>
<Q ID="38"> 12 40 51 48</Q>
<T ID="39"> 40 52 4</T>
<Q ID="40"> 46 53 8 52</Q>
<T ID="41"> 45 53 5</T>
<T ID="42"> 48 54 13</T>
<Q ID="43"> 51 52 14 54</Q>
<Q ID="44"> 56 59 57 46</Q>
<T ID="45"> 56 55 52</T>
<Q ID="46"> 59 58 8 55</Q>
<T ID="47"> 57 58 53</T>
<Q ID="48"> 61 62 56 51</Q>
<T ID="49"> 61 60 54</T>
<Q ID="50"> 62 55 14 60</Q>
<Q ID="51"> 63 27 64 59</Q>
<T ID="52"> 63 24 55</T>
<T ID="53"> 64 25 58</T>
<Q ID="54"> 65 30 63 62</Q>
<T ID="55"> 65 29 60</T>
</FACE>
<ELEMENT>
<R ID="0"> 0 1 2 3 4 </R>
<R ID="1"> 5 6 7 1 8 </R>
<R ID="2"> 9 10 4 11 12 </R>
<R ID="3"> 13 14 8 10 15 </R>
<R ID="4"> 16 17 12 18 19 </R>
<R ID="5"> 20 21 15 17 22 </R>
<R ID="6"> 23 24 25 26 16 </R>
<R ID="7"> 27 28 29 24 20 </R>
<R ID="8"> 30 31 32 33 34 </R>
<R ID="9"> 35 36 37 31 38 </R>
<R ID="10"> 34 39 40 41 2 </R>
<R ID="11"> 38 42 43 39 7 </R>
<R ID="12"> 44 45 46 47 40 </R>
<R ID="13"> 48 49 50 45 43 </R>
<R ID="14"> 51 52 19 53 46 </R>
<R ID="15"> 54 55 22 52 50 </R>
</ELEMENT>
<COMPOSITE>
<C ID="0"> R[0-15] </C>
<C ID="1"> F[0,3,5-6,9,11,13-14,18,21,23,25-30,32-33,35-37,41-42,44,47-49,51,53-55] </C>
</COMPOSITE>
<DOMAIN> C[0] </DOMAIN>
</GEOMETRY>
<EXPANSIONS>
<E COMPOSITE="C[0]" NUMMODES="3" TYPE="MODIFIED" FIELDS="u" />
</EXPANSIONS>
<CONDITIONS>
<PARAMETERS>
<P> Lambda = 1 </P>
</PARAMETERS>
<VARIABLES>
<V ID="0"> u </V>
</VARIABLES>
<BOUNDARYREGIONS>
<B ID="0"> C[1] </B>
</BOUNDARYREGIONS>
<BOUNDARYCONDITIONS>
<REGION REF="0">
<D VAR="u" VALUE="sin(PI/2*x)*sin(PI/2*y)*sin(PI/2*z)" />
</REGION>
</BOUNDARYCONDITIONS>
<FUNCTION NAME="Forcing">
<E VAR="u" VALUE="-(Lambda+3*PI*PI/4)*sin(PI/2*x)*sin(PI/2*y)*sin(PI/2*z)" />
</FUNCTION>
<FUNCTION NAME="ExactSolution">
<E VAR="u" VALUE="sin(PI/2*x)*sin(PI/2*y)*sin(PI/2*z)" />
</FUNCTION>
</CONDITIONS>
</NEKTAR>
......@@ -546,8 +546,8 @@ namespace Nektar
ASSERTL0(m_geom, "m_geom not define");
// get physical points defined in Geom
m_geom->FillGeom();
m_geom->FillGeom();
switch(m_geom->GetCoordim())
{
case 3:
......@@ -659,7 +659,7 @@ namespace Nektar
m_geom->FillGeom();
for(i = 0; i < m_geom->GetCoordDim(); ++i)
for(i = 0; i < m_geom->GetCoordim(); ++i)
{
coords[i] = m_geom->GetCoord(i,Lcoords);
}
......
......@@ -529,7 +529,7 @@ namespace Nektar
m_geom->FillGeom();
for(i = 0; i < m_geom->GetCoordDim(); ++i)
for(i = 0; i < m_geom->GetCoordim(); ++i)
{
coords[i] = m_geom->GetCoord(i,Lcoords);
}
......
......@@ -412,7 +412,7 @@ namespace Nektar
// m_geom->FillGeom(); // TODO: implement FillGeom()
for(i = 0; i < m_geom->GetCoordDim(); ++i)
for(i = 0; i < m_geom->GetCoordim(); ++i)
{
coords[i] = m_geom->GetCoord(i,Lcoords);
}
......@@ -623,7 +623,7 @@ namespace Nektar
{
// TODO: make sure 3D Laplacian is set up for Hex in three-dimensional in Standard Region.
// ASSERTL1(m_geom->GetCoordDim() == 2,"Standard Region Laplacian is only set up for Quads in two-dimensional");
ASSERTL1(m_geom->GetCoordDim() == 3,"Standard Region Laplacian is only set up for Hex in three-dimensional");
ASSERTL1(m_geom->GetCoordim() == 3,"Standard Region Laplacian is only set up for Hex in three-dimensional");
MatrixKey lap00key(StdRegions::eLaplacian00,
mkey.GetExpansionType(), *this);
MatrixKey lap01key(StdRegions::eLaplacian01,
......
......@@ -499,7 +499,7 @@ namespace Nektar
// m_geom->FillGeom(); // TODO: implement FillGeom()
for(i = 0; i < m_geom->GetCoordDim(); ++i)
for(i = 0; i < m_geom->GetCoordim(); ++i)
{
coords[i] = m_geom->GetCoord(i,Lcoords);
}
......
......@@ -29,16 +29,13 @@
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
// DEALINGS IN THE SOFTWARE.
//
// Description:
//
// Description: 3D geometry information.
//
////////////////////////////////////////////////////////////////////////////////
#include "pchSpatialDomains.h"
#include <SpatialDomains/Geometry3D.h>
namespace Nektar
{
namespace SpatialDomains
......@@ -51,39 +48,417 @@ namespace Nektar
Geometry(coordim)
{
ASSERTL0(m_coordim > 2,
"Coordinate dimension should be at least 3 for a 3D geometry");
"Coordinate dimension should be at least 3 for a 3D geometry.");
}
Geometry3D::~Geometry3D()
{
}
//---------------------------------------
// 3D Geometry Methods
//---------------------------------------
}; //end of namespace
}; //end of namespace
/**
* @brief Put all quadrature information into face/edge structure and
* backward transform.
*
* @see v_FillGeom()
*/
void Geometry3D::FillGeom()
{
Geometry3D::v_FillGeom();
}
void Geometry3D::GetLocCoords(
const Array<OneD, const NekDouble> &coords,
Array<OneD, NekDouble> &Lcoords)
{
v_GetLocCoords(coords, Lcoords);
}
/**
* @brief Given local collapsed coordinate Lcoord return the value of
* physical coordinate in direction i.
*/
NekDouble Geometry3D::GetCoord(
const int i, const Array<OneD, const NekDouble> &Lcoord)
{
return v_GetCoord(i, Lcoord);
}
//---------------------------------------
// Helper functions
//---------------------------------------
// $Log: Geometry3D.cpp,v $
// Revision 1.3 2008/08/26 02:28:07 ehan
// Added various virtual functions.
//
// Revision 1.2 2008/01/31 11:02:25 ehan
// Added constructor and destructor.
//
// Revision 1.1 2006/05/04 18:59:00 kirby
// *** empty log message ***
//
// Revision 1.12 2006/04/09 02:08:35 jfrazier
// Added precompiled header.
//
// Revision 1.11 2006/03/25 00:58:28 jfrazier
// Many changes dealing with fundamental structure and reading/writing.
//
// Revision 1.10 2006/03/12 07:42:02 sherwin
//
// Updated member names and StdRegions call. Still has not been compiled
//
// Revision 1.9 2006/02/19 01:37:33 jfrazier
// Initial attempt at bringing into conformance with the coding standard. Still more work to be done. Has not been compiled.
//
//
/**
* @brief Return the ID of edge i in this element.
*/
int Geometry3D::GetEid(int i) const
{
return v_GetEid(i);
}
/**
* @brief Return face i in this element.
*/
Geometry2DSharedPtr Geometry3D::GetFace(int i)
{
return v_GetFace(i);
}
/**
* @brief Return the orientation of face i in this element.
*/
StdRegions::FaceOrientation Geometry3D::GetFaceorient(const int i) const
{
return v_GetFaceorient(i);
}
/**
* @brief Return the ID of face i in this element.
*/
int Geometry3D::GetFid(int i) const
{
return v_GetFid(i);
}
/**
* @brief Return a reference to the physical space of co-ordinate
* dimension i.
*/
Array<OneD, NekDouble> &Geometry3D::UpdatePhys(const int i)
{
return v_UpdatePhys(i);
}
/**
* @brief Return the j-th basis of the i-th co-ordinate dimension.
*/
const LibUtilities::BasisSharedPtr Geometry3D::GetBasis(
const int i, const int j)
{
return v_GetBasis(i, j);
}
//---------------------------------------
// Element connection functions
//---------------------------------------
void Geometry3D::SetOwnData()
{
v_SetOwnData();
}
//---------------------------------------
// 3D Geometry Methods
//---------------------------------------
/**
* @brief 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.
*/
void Geometry3D::v_FillGeom()
{
if (m_state == ePtsFilled)
return;
int i,j,k;
for(i = 0; i < m_forient.size(); i++)
{
m_faces[i]->FillGeom();
int nFaceCoeffs = (*m_faces[i])[0]->GetNcoeffs();
Array<OneD, unsigned int> mapArray (nFaceCoeffs);
Array<OneD, int> signArray(nFaceCoeffs);
if (m_forient[i] < 4)
{
m_xmap[0]->GetFaceToElementMap(
i,m_forient[i],mapArray,signArray,
m_faces[i]->GetXmap(0)->GetEdgeNcoeffs(0),
m_faces[i]->GetXmap(0)->GetEdgeNcoeffs(1));
}
else
{
m_xmap[0]->GetFaceToElementMap(
i,m_forient[i],mapArray,signArray,
m_faces[i]->GetXmap(0)->GetEdgeNcoeffs(1),
m_faces[i]->GetXmap(0)->GetEdgeNcoeffs(0));
}
for(j = 0; j < m_coordim; j++)
{
const Array<OneD, const NekDouble> &coeffs =
(*m_faces[i])[j]->GetCoeffs();
for(k = 0; k < nFaceCoeffs; k++)
{
double v = signArray[k] * coeffs[k];
(m_xmap[j]->UpdateCoeffs())[mapArray[k]] = v;
}
}
}
for(i = 0; i < m_coordim; ++i)
{
m_xmap[i]->BwdTrans(m_xmap[i]->GetCoeffs (),
m_xmap[i]->UpdatePhys());
}
m_state = ePtsFilled;
}
void Geometry3D::v_GetLocCoords(
const Array<OneD, const NekDouble> &coords,
Array<OneD, NekDouble> &Lcoords)
{
ASSERTL0(false,
"This function must be implemented at a shape level.");
}
/**
* Generate the geometry factors for this element.
*/
void Geometry3D::v_GenGeomFactors(
const Array<OneD, const LibUtilities::BasisSharedPtr> &tbasis)
{
GeomType Gtype = eRegular;
GeomShapeType GSType = eQuadrilateral;
v_FillGeom();
// check to see if expansions are linear
for(int i = 0; i < m_coordim; ++i)
{
if (m_xmap[i]->GetBasisNumModes(0) != 2 ||
m_xmap[i]->GetBasisNumModes(1) != 2 ||
m_xmap[i]->GetBasisNumModes(2) != 2)
{
Gtype = eDeformed;
}
}
m_geomFactors = MemoryManager<GeomFactors3D>::AllocateSharedPtr(
Gtype, m_coordim, m_xmap, tbasis);
}
/**
* @brief Given local collapsed coordinate Lcoord return the value of
* physical coordinate in direction i.
*/
NekDouble Geometry3D::v_GetCoord(
const int i, const Array<OneD, const NekDouble> &Lcoord)
{
ASSERTL1(m_state == ePtsFilled, "Geometry is not in physical space.");
return m_xmap[i]->PhysEvaluate(Lcoord);
}
//---------------------------------------
// Helper functions
//---------------------------------------
/**
* @brief Return the co-ordinate mapping for dimension i.
*/
StdRegions::StdExpansion3DSharedPtr Geometry3D::v_GetXmap(const int i)
{
return m_xmap[i];
}
/**
* @brief Return the dimension of this element.
*/
int Geometry3D::v_GetShapeDim() const
{
return 3;
}
/**
* @brief Return the vertex ID of vertex i.
*/
int Geometry3D::v_GetVid(const int i) const
{
ASSERTL2(i >= 0 && i <= m_verts.num_elements()-1,
"Vertex ID must be between 0 and "+
boost::lexical_cast<string>(m_verts.num_elements()-1));
return m_verts[i]->GetVid();
}
/**
* @brief Return edge i of this element.
*/
const SegGeomSharedPtr Geometry3D::v_GetEdge(int i) const
{
ASSERTL2(i >= 0 && i <= m_edges.num_elements()-1,
"Edge ID must be between 0 and "+
boost::lexical_cast<string>(m_edges.num_elements()-1));
return m_edges[i];
}
/**
* @brief Return the orientation of edge i in this element.
*/
inline StdRegions::EdgeOrientation Geometry3D::v_GetEorient(
const int i) const
{
ASSERTL2(i >= 0 && i <= m_edges.num_elements()-1,
"Edge ID must be between 0 and "+
boost::lexical_cast<string>(m_edges.num_elements()-1));
return m_eorient[i];
}
/**
* @brief Return the ID of edge i in this element.
*/
int Geometry3D::v_GetEid(int i) const
{
ASSERTL2(i >= 0 && i <= m_edges.num_elements()-1,
"Edge ID must be between 0 and "+
boost::lexical_cast<string>(m_edges.num_elements()-1));
return m_edges[i]->GetEid();