Commit dd125064 authored by Douglas Serson's avatar Douglas Serson
Browse files

Merge branch 'master' into feature/CFS-RealGas

parents 78c5dd0a 5051a446
......@@ -22,17 +22,19 @@ v5.0.0
- Simplify RawType typedefs (!840)
- Remove unused files from BasicUtils (!841)
- Remove checks for old boost versions which are no longer supported (!841)
- Refactor ParseUtils to be more consistent (!843)
- Refactor ParseUtils to be more consistent (!843, !896)
- Added support for using the distance to a specific region (e.g. outlet) in the
function definitions for the Absorption Forcing (!769)
- Improve performance of DisContField2D::v_ExtractTracePhys (!824)
- Fix small bug in Jacobian Energy (!857)
- fix variable name overriding in file functions (!870)
- Adds CFI CAD engine back-end (!864)
- Adds CFI Mesh IO support (!864)
- Cleanup of CAD system data structures (!864)
- Fix mac OSX on buildbots (!876)
- Fix error from (!826) (!876)
- Fix minor bug in ARPACK thirdparty build cmake (!874)
- Switch MeshGraph to use factory pattern and add HDF5 geometry support (!900)
**NekMesh**:
- Add feature to read basic 2D geo files as CAD (!731)
......@@ -48,19 +50,34 @@ v5.0.0
- Additional curve types in GEO reader: BSpline, Circle, Ellipse (!800)
- Fix default command line argument value (!823)
- Add projection meshing module which can curve linear meshes with CAD (!826)
- XML meshes now write with provenance information, including information about
their source, for debugging purposes (!872)
- Force 3-node loops to avoid degenerate 1-triangle faces (!875)
- Smooth BL normals in 2D when normals intersect or cause invalid macro BL
elements (!877)
- Revert triangle code to ThirdParty library (!883)
- Fix coinciding nodes issue with very fine meshes (!883)
- Skip CFI groups of bodies and non-numbered nodes (!891)
- Add ability to space out 2D BL nodes to better fit local target Delta (!890)
- Fix automatic peralign call in 2D periodic meshing (!888)
**FieldConvert**:
- Add input module for Semtex field files (!777)
- Fixed interppoints module (!760)
- Fix OutputTecplot in 2DH1D (!818)
- Move StreamFunction utility to a FieldConvert module (!809)
- Extend wss module to compressible flows (!810)
- Allow explicitly setting bool options of FieldConvert modules as false (!811)
- Enable output to multiple files (!844)
- Allow using xml file without expansion tag in FieldConvert (!849)
- Add Lambda 2 vortex detection criteria (!882)
- Add module for modifying/adding fields from expressions (!889, !903)
- Add module for evaluating the mean of variables on the domain (!894)
**IncNavierStokesSolver**
- Replace steady-state check based on difference of norms by check based on
norm of the difference, to be consistent with the compressible solver (!832)
- Updated SVV to allow for the DGKernel extension (!851)
**CompressibleFlowSolver**
- Add 3D regression tests (!567)
......@@ -70,6 +87,10 @@ v5.0.0
seg, quad and hex elements (!771, !862)
- Introduce equations of state to account for real gas effects (!880)
**APESolver:**
- Added two new boundary conditions to the APE system: RiemannInvariantBC
and WhiteNoise (!782)
**Documentation**:
- Added the developer-guide repository as a submodule (!751)
......@@ -79,16 +100,25 @@ v4.4.2
- Fix evaluation of points (e.g. HistoryPoints, Interpolation to pts) close to
the interface of two elements (!836)
- Fix deadlock in Hdf5 with homogeneous expansions (!858)
- Fix a few memory leaks in polylib (!863)
- Fix a crash when Interpolator is called on an empty field (!869)
- Fix petsc compile without MPI (!873)
- Fix calculation of BLPoints (!892)
- Fix deadlock in DiffusionLDG (!885)
- Fix uninitialised coefficients in DirectFull solver (!898)
**NekMesh**
- Fix missing periodic boundary meshing and boundary layer mesh adjustment
configurations in 2D (!859)
- Fix 2D BL splitting where out-of-plane nodes would be created (!887)
**Documentation**:
- Fix sign of the viscous term in the velocity correction scheme equations in
the user guide (!856)
**FieldConvert**
- Allow passing input name with trailing separator (!879)
v4.4.1
------
**Library**
......@@ -120,6 +150,7 @@ v4.4.1
IMEXGear, CNAB, 2nd order IMEX-DIRK, 3rd order IMEX-DIRK (!854)
- Fix bug due to subtractive cancellation in polylib routines (!778)
**FieldConvert:**
- Fix issue with field ordering in the interppointdatatofld module (!754)
- Fix issue with FieldConvert when range flag used (!761)
......
......@@ -245,6 +245,7 @@ INCLUDE (ThirdPartyArpack)
INCLUDE (ThirdPartyMPI)
INCLUDE (ThirdPartyVTK)
INCLUDE (ThirdPartyOCE)
INCLUDE (ThirdPartyTriangle)
INCLUDE (ThirdPartyTetGen)
INCLUDE (ThirdPartyCCM)
INCLUDE (FindCFI)
......
########################################################################
#
# ThirdParty configuration for Nektar++
#
# TRIANGLE
#
########################################################################
IF(NEKTAR_USE_MESHGEN)
SET(BUILD_TRIANGLE ON)
OPTION(THIRDPARTY_BUILD_TRIANGLE
"Build Triangle library from ThirdParty." ${BUILD_TRIANGLE})
IF (THIRDPARTY_BUILD_TRIANGLE)
INCLUDE(ExternalProject)
EXTERNALPROJECT_ADD(
triangle-1.6
PREFIX ${TPSRC}
URL ${TPURL}/triangle.zip
URL_MD5 357cb7107f51f3f89940c47435d4fa49
STAMP_DIR ${TPBUILD}/stamp
DOWNLOAD_DIR ${TPSRC}
SOURCE_DIR ${TPSRC}/triangle-1.6
BINARY_DIR ${TPBUILD}/triangle-1.6
TMP_DIR ${TPBUILD}/triangle-1.6-tmp
INSTALL_DIR ${TPDIST}
CONFIGURE_COMMAND ${CMAKE_COMMAND}
-G ${CMAKE_GENERATOR}
-DCMAKE_C_COMPILER:FILEPATH=${CMAKE_C_COMPILER}
-DCMAKE_CXX_COMPILER:FILEPATH=${CMAKE_CXX_COMPILER}
-DCMAKE_INSTALL_PREFIX:PATH=${TPDIST}
${TPSRC}/triangle-1.6
)
THIRDPARTY_LIBRARY(TRIANGLE_LIBRARY STATIC triangle
DESCRIPTION "Triangle library")
SET(TRIANGLE_INCLUDE_DIR ${TPDIST}/include CACHE FILEPATH
"Triangle include" FORCE)
MESSAGE(STATUS "Build Triangle: ${TRIANGLE_LIBRARY}")
SET(TRIANGLE_CONFIG_INCLUDE_DIR ${TPINC})
ELSE()
ADD_CUSTOM_TARGET(triangle-1.6 ALL)
MESSAGE(STATUS "Found Triangle: ${TRIANGLE_LIBRARY}")
SET(TRIANGLE_CONFIG_INCLUDE_DIR ${TRIANGLE_INCLUDE_DIR})
ENDIF (THIRDPARTY_BUILD_TRIANGLE)
MARK_AS_ADVANCED(TRIANGLE_LIBRARY)
MARK_AS_ADVANCED(TRIANGLE_INCLUDE_DIR)
INCLUDE_DIRECTORIES(${TRIANGLE_INCLUDE_DIR})
ENDIF(NEKTAR_USE_MESHGEN)
......@@ -51,11 +51,57 @@ Currently, only \inltt{APEUpwind} supported.
\subsection{Functions}
\begin{itemize}
\item \inltt{BaseFlow} Baseflow $(\overline{u}_i, \overline{p}, \overline{\rho})$ defined by the variables \inltt{u0,v0,w0,p0,rho0}
\item \inltt{Source} Source term $\overline{c}^2 q_c$
\item \inltt{InitialConditions}
\end{itemize}
\subsection{Boundary Conditions}
In addition to plain Dirichlet and Neumann boundary conditions, the APESolver features a wall boundray condition, a non-reflecting boundary and a white noise boundary condition.
\begin{itemize}
\item Rigid (Slip-) Wall Boundary Condition
\begin{lstlisting}[style=XmlStyle]
<BOUNDARYCONDITIONS>
<REGION REF="0">
<D VAR="p" USERDEFINEDTYPE="Wall" VALUE="0" />
<D VAR="u" USERDEFINEDTYPE="Wall" VALUE="0" />
<D VAR="v" USERDEFINEDTYPE="Wall" VALUE="0" />
<D VAR="w" USERDEFINEDTYPE="Wall" VALUE="0" />
</REGION>
</BOUNDARYCONDITIONS>
\end{lstlisting}
This BC imposes zero wall-normal perturbation velocity in a way that is more robust than using a Dirichlet boundary condition directly.
\item Non-Reflecting Boundary Condition
\begin{lstlisting}[style=XmlStyle]
<BOUNDARYCONDITIONS>
<REGION REF="0">
<D VAR="p" USERDEFINEDTYPE="RiemannInvariantBC"/>
<D VAR="u" USERDEFINEDTYPE="RiemannInvariantBC"/>
<D VAR="v" USERDEFINEDTYPE="RiemannInvariantBC"/>
<D VAR="w" USERDEFINEDTYPE="RiemannInvariantBC"/>
</REGION>
</BOUNDARYCONDITIONS>
\end{lstlisting}
The Riemann-Invariant BC approximates a non-reflecting (r.g. Farfield) boundary condition by setting incoming invariants to zero.
\item White Noise Boundary Condition
\begin{lstlisting}[style=XmlStyle]
<BOUNDARYCONDITIONS>
<REGION REF="0">
<D VAR="p" USERDEFINEDTYPE="Wall" VALUE="10" />
<D VAR="u" USERDEFINEDTYPE="Wall" VALUE="10" />
<D VAR="v" USERDEFINEDTYPE="Wall" VALUE="10" />
<D VAR="w" USERDEFINEDTYPE="Wall" VALUE="10" />
</REGION>
</BOUNDARYCONDITIONS>
\end{lstlisting}
The white noise BC imposes a stochastic, uniform pressure at the boundary. The implementation uses a Mersenne-Twister pseudo random number generator to generate white Gaussian noise.
The standard deviation $\sigma$ of the pressure is specified by the \inltt{VALUE} attribute.
\end{itemize}
\section{Examples}
\subsection{Aeroacoustic Wave Propagation}
In this section we explain how to set up a simple simulation of aeroacoustics in
......@@ -108,10 +154,6 @@ eventual source terms using the following functions:
<E VAR="rho0" VALUE="Rho0" />
</FUNCTION>
<FUNCTION NAME="Source">
<E VAR="S" VALUE="0" />
</FUNCTION>
<!-- Gaussian pulse located at the origin -->
<FUNCTION NAME="InitialConditions">
<E VAR="p" VALUE="100*exp(-32*(((x)*(x))+((y)*(y))))" />
......
......@@ -155,7 +155,9 @@ Specifically, FieldConvert has these additional functionalities
\begin{enumerate}
\item \inltt{C0Projection}: Computes the C0 projection of a given output file;
\item \inltt{QCriterion}: Computes the Q-Criterion for a given output file;
\item \inltt{L2Criterion}: Computes the Lambda 2 Criterion for a given output file;
\item \inltt{addcompositeid}: Adds the composite ID of an element as an additional field;
\item \inltt{fieldfromstring}: Modifies or adds a new field from an expression involving the existing fields;
\item \inltt{addFld}: Sum two .fld files;
\item \inltt{combineAvg}: Combine two \nekpp binary output (.chk or .fld) field file containing averages of fields (and
possibly also Reynolds stresses) into single file;
......@@ -173,6 +175,7 @@ possibly also Reynolds stresses) into single file;
\item \inltt{isocontour}: Extract an isocontour of ``fieldid'' variable and at value ``fieldvalue''. Optionally ``fieldstr'' can be specified for a string definition or ``smooth'' for smoothing;
\item \inltt{jacobianenergy}: Shows high frequency energy of Jacobian;
\item \inltt{qualitymetric}: Evaluate a quality metric of the underlying mesh to show mesh quality;
\item \inltt{mean}: Evaluate the mean of variables on the domain;
\item \inltt{meanmode}: Extract mean mode (plane zero) of 3DH1D expansions;
\item \inltt{pointdatatofld}: Given discrete data at quadrature points
project them onto an expansion basis and output fld file;
......@@ -250,6 +253,21 @@ to visualise the result either in Tecplot, Paraview or VisIt.
%
%
\subsection{Calculate $\lambda_2$: \textit{L2Criterion} module}
To perform the $\lambda_2$ vortex detection calculation and obtain an output
data containing the values of the $\lambda_2$ eigenvalue, the user can run
%
\begin{lstlisting}[style=BashInputStyle]
FieldConvert -m L2Criterion test.xml test.fld test-L2Crit.fld
\end{lstlisting}
%
where the file \inltt{test-L2Crit.fld} can be processed in a similar
way as described in section \ref{s:utilities:fieldconvert:sub:convert}
to visualise the result either in Tecplot, Paraview or VisIt.
%
%
%
\subsection{Add composite ID: \textit{addcompositeid} module}
When dealing with a geometry that has many surfaces, we need to identify the
composites to assign boundary conditions. To assist in this, FieldConvert has a
......@@ -265,7 +283,28 @@ variable that contains the composite ID. To assist in boundary identification,
the input file \inlsh{mesh.xml} should be a surface XML file that can be
obtained through the \nm \inltt{extract} module (see section
\ref{s:utilities:nekmesh:extract}).
%
%
%
\subsection{Add new field: \textit{fieldfromstring} module}
To modify or create a new field using an expression involving the existing fields, one can
use the \inltt{fieldfromstring} module of FieldConvert
%
\begin{lstlisting}[style=BashInputStyle]
FieldConvert -m fieldfromstring:fieldstr="x+y+u":fieldname="result" \
file1.xml file2.fld file3.fld
\end{lstlisting}
%
In this case \inltt{fieldstr} is a required parameter describing a function of
the coordinates and the existing variables, and \inltt{fieldname} is an optional
parameter defining the name of the new or modified field (the default is newfield).
\inltt{file3.fld} is the output containing both the original and the new fields,
and can be processed in a similar way as described
in section \ref{s:utilities:fieldconvert:sub:convert} to visualise
the result either in Tecplot, Paraview or VisIt.
%
%
%
\subsection{Sum two .fld files: \textit{addFld} module}
To sum two .fld files one can use the \inltt{addFld} module of FieldConvert
%
......@@ -768,7 +807,20 @@ Two quality metrics are implemented that produce scalar fields $Q$:
%
%
%
\subsection{Evaluate the mean of variables on the domain: \textit{mean} module}
To evaluate the mean of variables on the domain one can use the
\inltt{mean} module of FieldConvert
%
\begin{lstlisting}[style=BashInputStyle]
FieldConvert -m mean file1.xml file2.fld out.stdout
\end{lstlisting}
%
This module does not create an output file which is reinforced by the
out.stdout option. The integral and mean for each field variable are
then printed to the stdout.
%
%
%
\subsection{Extract mean mode of 3DH1D expansion: \textit{meanmode} module}
To obtain a 2D expansion containing the mean mode (plane zero in Fourier space) of a
......
SUBDIRS(LibUtilities StdRegions LocalRegions Collections MultiRegions)
SUBDIRS(LibUtilities StdRegions SpatialDomains LocalRegions Collections MultiRegions)
ADD_NEKTAR_EXECUTABLE(PartitionAnalyse
COMPONENT demos DEPENDS LibUtilities SOURCES PartitionAnalyse.cpp)
ADD_NEKTAR_EXECUTABLE(FoundationDemo
COMPONENT demos DEPENDS LibUtilities SOURCES FoundationDemo.cpp)
ADD_NEKTAR_EXECUTABLE(NodalDemo
......
......@@ -9,7 +9,7 @@ using namespace Nektar;
static double solutionpoly(double x, int order);
static double solutionfourier(double x, int order, double a, double b);
// This routine projects a polynomial or trigonmetric functions which
// This routine projects a polynomial or trigonmetric functions which
// has energy in all mdoes of the expansions and report an error.
int main(int argc, char *argv[])
......@@ -38,7 +38,7 @@ int main(int argc, char *argv[])
fprintf(stderr,"\t Chebyshev = 13\n");
fprintf(stderr,"\t Monomial = 14\n");
fprintf(stderr,"\t FourierSingleMode = 15\n");
fprintf(stderr,"Note type = 1,2,4,5,7,8 are for higher dimensional basis\n");
exit(1);
......@@ -81,8 +81,8 @@ int main(int argc, char *argv[])
const double dZero=0.0;
SpatialDomains::PointGeomSharedPtr vert1 = MemoryManager<SpatialDomains::PointGeom>::AllocateSharedPtr(one,zero,x[0],dZero,dZero);
SpatialDomains::PointGeomSharedPtr vert2 = MemoryManager<SpatialDomains::PointGeom>::AllocateSharedPtr(one,zero,x[1],dZero,dZero);
SpatialDomains::SegGeomSharedPtr geom = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(zero,vert1,vert2);
geom->SetOwnData();
SpatialDomains::PointGeomSharedPtr verts[2] = {vert1,vert2};
SpatialDomains::SegGeomSharedPtr geom = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(zero,1,verts);
const LibUtilities::PointsKey Pkey(nq,Qtype);
const LibUtilities::BasisKey Bkey(btype,order,Pkey);
......@@ -90,7 +90,7 @@ int main(int argc, char *argv[])
//-----------------------------------------------
//----------------------------------------------
// Define solution to be projected
// Define solution to be projected
Array<OneD,NekDouble> xc(nq);
E->GetCoords(xc);
......@@ -116,14 +116,14 @@ int main(int argc, char *argv[])
Array<OneD, NekDouble> phys (nq);
E->FwdTrans(sol, coeffs);
//---------------------------------------------
//-------------------------------------------
// Backward Transform Solution to get projected values
E->BwdTrans(coeffs, phys);
//-------------------------------------------
//-------------------------------------------
//--------------------------------------------
// Calculate L_inf error
// Calculate L_inf error
cout << "L infinity error: " << E->Linf(phys, sol) << endl;
cout << "L 2 error: " << E->L2 (phys, sol) << endl;
//--------------------------------------------
......
......@@ -21,7 +21,7 @@ NekDouble Quad_sol(NekDouble x, NekDouble y, int order1, int order2,
int main(int argc, char *argv[])
{
int i;
int order1,order2, nq1,nq2;
LibUtilities::PointsType Qtype1,Qtype2;
LibUtilities::BasisType btype1 = LibUtilities::eOrtho_A;
......@@ -31,24 +31,23 @@ int main(int argc, char *argv[])
StdRegions::StdExpansion2D *E = NULL;
Array<OneD, NekDouble> sol;
Array<OneD, NekDouble> coords(8);
StdRegions::Orientation edgeDir = StdRegions::eForwards;
if((argc != 16)&&(argc != 14))
{
// arg[0] arg[1] arg[2] arg[3] arg[4] arg[5] arg[6] arg[7] arg[8] arg[9] arg[10] arg[11] arg[12] arg[13]
fprintf(stderr,"Usage: Project2D RegionShape Type1 Type2 order1 order2 nq1 nq2 x1, y1, x2, y2, x3, y3 \n");
fprintf(stderr,"Example : ./LocProject2D-g 2 4 5 3 3 4 4 .1 -.5 .6 .1 .3 .2 \n" );
fprintf(stderr,"Where RegionShape is an integer value which "
"dictates the region shape:\n");
fprintf(stderr,"\t Triangle = 3\n");
fprintf(stderr,"\t Quadrilateral = 4\n");
fprintf(stderr,"Where type is an integer value which "
"dictates the basis as:\n");
fprintf(stderr,"\t Ortho_A = 1\n");
fprintf(stderr,"\t Ortho_B = 2\n");
fprintf(stderr,"\t Modified_A = 4\n");
......@@ -61,24 +60,24 @@ int main(int argc, char *argv[])
fprintf(stderr,"\t Monomial = 14\n");
fprintf(stderr,"\t Nodal tri (Electro) = 15\n");
fprintf(stderr,"\t Nodal tri (Fekete) = 16\n");
fprintf(stderr,"Note type = 3,6,7,8 are for three-dimensional basis\n");
fprintf(stderr,"The last series of values are the coordinates\n");
exit(1);
}
regionshape = (LibUtilities::ShapeType) atoi(argv[1]);
// Check to see if 2D region
if((regionshape != LibUtilities::eTriangle)&&(regionshape != LibUtilities::eQuadrilateral))
{
NEKERROR(ErrorUtil::efatal,"This shape is not a 2D region");
}
int btype1_val = atoi(argv[2]);
int btype2_val = atoi(argv[3]);
if(( btype1_val <= 14)&&( btype2_val <= 14))
{
btype1 = (LibUtilities::BasisType) btype1_val;
......@@ -88,7 +87,7 @@ int main(int argc, char *argv[])
{
btype1 = LibUtilities::eOrtho_A;
btype2 = LibUtilities::eOrtho_B;
if(btype1_val == 15)
{
NodalType = LibUtilities::eNodalTriElec;
......@@ -98,8 +97,8 @@ int main(int argc, char *argv[])
NodalType = LibUtilities::eNodalTriFekete;
}
}
// Check to see that correct Expansions are used
switch(regionshape)
{
......@@ -109,7 +108,7 @@ int main(int argc, char *argv[])
NEKERROR(ErrorUtil::efatal,
"Basis 1 cannot be of type Ortho_B or Modified_B");
}
break;
case LibUtilities::eQuadrilateral:
if((btype1 == LibUtilities::eOrtho_B)||(btype1 == LibUtilities::eOrtho_C)||
......@@ -118,7 +117,7 @@ int main(int argc, char *argv[])
NEKERROR(ErrorUtil::efatal,
"Basis 1 is for 2 or 3D expansions");
}
if((btype2 == LibUtilities::eOrtho_B)||(btype2 == LibUtilities::eOrtho_C)||
(btype2 == LibUtilities::eModified_B)||(btype2 == LibUtilities::eModified_C))
{
......@@ -129,15 +128,15 @@ int main(int argc, char *argv[])
default:
NEKERROR(ErrorUtil::efatal, "Not a valid 2D expansion.");
}
order1 = atoi(argv[4]);
order2 = atoi(argv[5]);
nq1 = atoi(argv[6]);
nq2 = atoi(argv[7]);
sol = Array<OneD, NekDouble>(nq1*nq2);
if(btype1 != LibUtilities::eFourier)
{
Qtype1 = LibUtilities::eGaussLobattoLegendre;
......@@ -146,7 +145,7 @@ int main(int argc, char *argv[])
{
Qtype1 = LibUtilities::eFourierEvenlySpaced;
}
if(btype2 != LibUtilities::eFourier)
{
if (regionshape == LibUtilities::eTriangle) {
......@@ -161,10 +160,10 @@ int main(int argc, char *argv[])
{
Qtype2 = LibUtilities::eFourierEvenlySpaced;
}
//-----------------------------------------------
// Define a 2D expansion based on basis definition
switch(regionshape)
{
case LibUtilities::eTriangle:
......@@ -175,7 +174,7 @@ int main(int argc, char *argv[])
coords[3] = atof(argv[11]);
coords[4] = atof(argv[12]);
coords[5] = atof(argv[13]);
// Set up coordinates
SpatialDomains::PointGeomSharedPtr verts[3];
const int zero = 0;
......@@ -185,21 +184,18 @@ int main(int argc, char *argv[])
verts[0] = MemoryManager<SpatialDomains::PointGeom>::AllocateSharedPtr(two,zero,coords[0],coords[1],dZero);
verts[1] = MemoryManager<SpatialDomains::PointGeom>::AllocateSharedPtr(two,one,coords[2],coords[3],dZero);
verts[2] = MemoryManager<SpatialDomains::PointGeom>::AllocateSharedPtr(two,two,coords[4],coords[5],dZero);
SpatialDomains::PointGeomSharedPtr verta[2] = {verts[0],verts[1]};
SpatialDomains::PointGeomSharedPtr vertb[2] = {verts[1],verts[2]};
SpatialDomains::PointGeomSharedPtr vertc[2] = {verts[2],verts[0]};
// Set up Edges
SpatialDomains::SegGeomSharedPtr edges[3];
edges[0] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(zero,verts[0],verts[1]);
edges[1] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(one,verts[1],verts[2]);
edges[2] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(two,verts[2],verts[0]);
StdRegions::Orientation eorient[3];
eorient[0] = edgeDir;
eorient[1] = edgeDir;
eorient[2] = edgeDir;
SpatialDomains::TriGeomSharedPtr geom = MemoryManager<SpatialDomains::TriGeom>::AllocateSharedPtr(zero,verts,edges,eorient);
geom->SetOwnData();
edges[0] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(zero,2,verta);
edges[1] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(one,2,vertb);
edges[2] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(two,2,vertc);
SpatialDomains::TriGeomSharedPtr geom = MemoryManager<SpatialDomains::TriGeom>::AllocateSharedPtr(zero,edges);
const LibUtilities::PointsKey Pkey1(nq1,Qtype1);
const LibUtilities::PointsKey Pkey2(nq2,Qtype2);
const LibUtilities::BasisKey Bkey1(btype1,order1,Pkey1);
......@@ -213,11 +209,11 @@ int main(int argc, char *argv[])
{
E = new LocalRegions::TriExp(Bkey1,Bkey2,geom);
}
Array<OneD,NekDouble> x = Array<OneD,NekDouble>(nq1*nq2);
Array<OneD,NekDouble> y = Array<OneD,NekDouble>(nq1*nq2);
E->GetCoords(x,y);
//----------------------------------------------
// Define solution to be projected
for(i = 0; i < nq1*nq2; ++i)
......@@ -225,7 +221,7 @@ int main(int argc, char *argv[])
sol[i] = Tri_sol(x[i],y[i],order1,order2);
}
//----------------------------------------------
}
break;
case LibUtilities::eQuadrilateral:
......@@ -239,7 +235,7 @@ int main(int argc, char *argv[])
coords[5] = atof(argv[13]);
coords[6] = atof(argv[14]);
coords[7] = atof(argv[15]);
// Set up coordinates
const int zero=0;
const int one=1;
......@@ -251,43 +247,41 @@ int main(int argc, char *argv[])
verts[1] = MemoryManager<SpatialDomains::PointGeom>::AllocateSharedPtr(two,one,coords[2],coords[3],dZero);
verts[2] = MemoryManager<SpatialDomains::PointGeom>::AllocateSharedPtr(two,two,coords[4],coords[5],dZero);
verts[3] = MemoryManager<SpatialDomains::PointGeom>::AllocateSharedPtr(two,three,coords[6],coords[7],dZero);
SpatialDomains::PointGeomSharedPtr verta[2] = {verts[0],verts[1]};
SpatialDomains::PointGeomSharedPtr vertb[2] = {verts[1],verts[2]};
SpatialDomains::PointGeomSharedPtr vertc[2] = {verts[2],verts[3]};
SpatialDomains::PointGeomSharedPtr vertd[2] = {verts[3],verts[0]};
// Set up Edges
SpatialDomains::SegGeomSharedPtr edges[4];
edges[0] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(zero,verts[0],verts[1]);
edges[1] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(one,verts[1],verts[2]);
edges[2] = MemoryManager<SpatialDomains::SegGeom>::AllocateSharedPtr(two,verts[2],verts[3]);
edges[3]