Commit 335d1b73 authored by Chris Cantwell's avatar Chris Cantwell
Browse files

Merge branch 'master' into feature/SpechpSVV

Conflicts:
	library/StdRegions/StdExpansion.cpp
parents 0d8ca1b7 c12b0358
......@@ -4,7 +4,7 @@ options are: None(CMAKE_CXX_FLAGS or CMAKE_C_FLAGS used) Debug Release
RelWithDebInfo MinSizeRel.")
PROJECT(Nektar++)
# Helps organize projects in IDEs.
SET_PROPERTY(GLOBAL PROPERTY USE_FOLDERS ON)
......@@ -68,7 +68,7 @@ OPTION(NEKTAR_FULL_DEBUG "Enable Full Debugging." OFF)
MARK_AS_ADVANCED(NEKTAR_FULL_DEBUG)
IF (${CMAKE_COMPILER_IS_GNUCXX})
OPTION(NEKTAR_ENABLE_PROFILE "Uses -pg compiler flag" ON)
OPTION(NEKTAR_ENABLE_PROFILE "Uses -pg compiler flag" OFF)
MARK_AS_ADVANCED(NEKTAR_ENABLE_PROFILE)
ENDIF (${CMAKE_COMPILER_IS_GNUCXX})
......@@ -160,6 +160,7 @@ INCLUDE (ThirdPartyFFTW)
INCLUDE (ThirdPartyArpack)
INCLUDE (ThirdPartyMPI)
INCLUDE (ThirdPartyVTK)
INCLUDE (ThirdPartyQT4)
IF( NEKTAR_USE_MKL )
INCLUDE (FindMKL)
......
SET(METIS_SEARCH_PATHS
${CMAKE_SOURCE_DIR}/ThirdParty/modmetis-4.0/
${CMAKE_SOURCE_DIR}/ThirdParty/modmetis-4.0/build/
${CMAKE_SOURCE_DIR}/../ThirdParty/modmetis-4.0/
${CMAKE_SOURCE_DIR}/../ThirdParty/modmetis-4.0/build
${CMAKE_SOURCE_DIR}/ThirdParty/metis-5.0.2/
${CMAKE_SOURCE_DIR}/ThirdParty/metis-5.0.2/build/
${CMAKE_SOURCE_DIR}/../ThirdParty/metis-5.0.2/
${CMAKE_SOURCE_DIR}/../ThirdParty/metis-5.0.2/build
${CMAKE_SOURCE_DIR}/ThirdParty/dist/lib
${CMAKE_SOURCE_DIR}/../ThirdParty/dist/lib)
FIND_LIBRARY(METIS_LIB NAMES modmetis PATHS ${METIS_SEARCH_PATHS})
FIND_LIBRARY(METIS_LIB NAMES metis PATHS ${METIS_SEARCH_PATHS})
SET(METIS_FOUND FALSE)
IF (METIS_LIB)
SET(METIS_FOUND TRUE)
GET_FILENAME_COMPONENT(METIS_PATH ${METIS_LIB} PATH)
INCLUDE_DIRECTORIES(${METIS_PATH}/../include)
MARK_AS_ADVANCED(METIS_PATH)
MARK_AS_ADVANCED(METIS_LIB)
ENDIF (METIS_LIB)
......
......@@ -221,7 +221,7 @@ MACRO(ADD_NEKTAR_LIBRARY name component type)
# NIST Sparse BLAS only static, so link into Nektar libraries directly.
TARGET_LINK_LIBRARIES( ${name} ${NIST_SPARSE_BLAS} ${METIS_LIB})
ADD_DEPENDENCIES(${name} spblastk0.9b modmetis-4.0 boost tinyxml zlib)
ADD_DEPENDENCIES(${name} spblastk0.9b modmetis-5.0.2 boost tinyxml zlib)
SET_PROPERTY(TARGET ${name} PROPERTY FOLDER ${component})
IF (NEKTAR_USE_MPI)
TARGET_LINK_LIBRARIES( ${name} ${GSMPI_LIBRARY} ${XXT_LIBRARY})
......
......@@ -4,17 +4,18 @@ SET(THIRDPARTY_BUILD_METIS ON CACHE BOOL
IF (THIRDPARTY_BUILD_METIS)
INCLUDE(ExternalProject)
EXTERNALPROJECT_ADD(
modmetis-4.0
modmetis-5.0.2
PREFIX ${TPSRC}
URL ${TPURL}/modmetis-4.0.tar.bz2
URL_MD5 "4267355e04dcc2ef36f9889a98c923a5"
URL ${TPURL}/modmetis-5.0.2.tar.bz2
URL_MD5 "ffbdc6a50283934389a0b3b0c32b62c0"
DOWNLOAD_DIR ${TPSRC}
CONFIGURE_COMMAND ${CMAKE_COMMAND} -DCMAKE_INSTALL_PREFIX:PATH=${TPSRC}/dist ${TPSRC}/src/modmetis-4.0
CONFIGURE_COMMAND ${CMAKE_COMMAND} -DCMAKE_INSTALL_PREFIX:PATH=${TPSRC}/dist -DCMAKE_C_FLAGS:STRING=-fPIC -DGKLIB_PATH:PATH=${TPSRC}/src/modmetis-5.0.2/GKlib ${TPSRC}/src/modmetis-5.0.2
)
SET(METIS_LIB modmetis CACHE FILEPATH
SET(METIS_LIB metis CACHE FILEPATH
"METIS library" FORCE)
MARK_AS_ADVANCED(METIS_LIB)
LINK_DIRECTORIES(${TPSRC}/dist/lib)
INCLUDE_DIRECTORIES(${TPSRC}/dist/include)
ELSE (THIRDPARTY_BUILD_METIS)
INCLUDE (FindMetis)
ENDIF (THIRDPARTY_BUILD_METIS)
......
# QT4
OPTION(NEKTAR_USE_QT4
"Use QT4 for graphical tools and utilities." OFF)
IF( NEKTAR_USE_QT4 )
FIND_PACKAGE(Qt4 COMPONENTS QtCore QtGui REQUIRED)
INCLUDE(${QT_USE_FILE})
ADD_DEFINITIONS(${QT_DEFINITIONS})
ENDIF( NEKTAR_USE_QT4 )
......@@ -8,6 +8,7 @@ library/CPackRPM.cmake
library/Demos/StdRegions/ExtraDemos
library/Demos/MultiRegions/ExtraDemos
solvers/ADR2DManifoldSolver
solvers/CardiacEPSolver/Utilities
solvers/FitzHughNagumoSolver
solvers/ImageWarpingSolver
solvers/PulseWaveSolver
......
......@@ -110,6 +110,7 @@ ADD_NEKTAR_TEST(Helmholtz2D_HDG_P7_Modes_AllBCs)
ADD_NEKTAR_TEST(Helmholtz3D_CG_Hex)
ADD_NEKTAR_TEST(Helmholtz3D_CG_Hex_AllBCs)
ADD_NEKTAR_TEST(Helmholtz3D_CG_Hex_AllBCs_iter_ml)
ADD_NEKTAR_TEST(Helmholtz3D_CG_Hex_AllBCs_iter_sc_cont)
ADD_NEKTAR_TEST(Helmholtz3D_CG_Tet)
ADD_NEKTAR_TEST(Helmholtz3D_CG_Prism)
ADD_NEKTAR_TEST(Helmholtz3D_CG_Prism_iter_ml)
......
<?xml version="1.0" encoding="utf-8"?>
<test>
<description>Helmholtz 3D CG, hexes, mixed BCs, iterative SC, contiguous storage strategy</description>
<executable>Helmholtz3D</executable>
<parameters>-I GlobalSysSoln=IterativeStaticCond -I LocalMatrixStorageStrategy=Contiguous Helmholtz3D_Hex_AllBCs_P6.xml</parameters>
<files>
<file description="Session File">Helmholtz3D_Hex_AllBCs_P6.xml</file>
</files>
<metrics>
<metric type="L2" id="1">
<value tolerance="1e-8">0.000416575</value>
</metric>
<metric type="Linf" id="2">
<value tolerance="1e-8">0.000871589</value>
</metric>
</metrics>
</test>
......@@ -341,7 +341,13 @@ namespace Nektar
try
{
// Attempt partitioning using METIS.
Metis::PartGraphVKway(nGraphVerts, xadj, adjncy, vwgt, vsize, npart, vol, part);
// Check METIS produced a valid partition and fix if not.
CheckPartitions(part);
// Distribute partitioning to all processes.
for (i = 1; i < m_comm->GetSize(); ++i)
{
m_comm->Send(i, part);
......@@ -376,6 +382,39 @@ namespace Nektar
}
}
void MeshPartition::CheckPartitions(Array<OneD, int> &pPart)
{
unsigned int i = 0;
unsigned int cnt = 0;
const unsigned int npart = m_comm->GetSize();
bool valid = true;
// Check that every process has at least one element assigned
for (i = 0; i < npart; ++i)
{
cnt = std::count(pPart.begin(), pPart.end(), i);
if (cnt == 0)
{
valid = false;
}
}
// If METIS produced an invalid partition, repartition naively.
// Elements are assigned to processes in a round-robin fashion.
// It is assumed that METIS failure only occurs when the number of
// elements is approx. the number of processes, so this approach
// should not be too inefficient communication-wise.
if (!valid)
{
for (i = 0; i < pPart.num_elements(); ++i)
{
pPart[i] = i % npart;
}
}
}
void MeshPartition::OutputPartition(
LibUtilities::SessionReaderSharedPtr& pSession,
BoostSubGraph& pGraph,
......
......@@ -180,6 +180,7 @@ namespace Nektar
void PartitionGraph(BoostSubGraph& pGraph,
BoostSubGraph& pLocalPartition);
void OutputPartition(SessionReaderSharedPtr& pSession, BoostSubGraph& pGraph, TiXmlElement* pGeometry);
void CheckPartitions(Array<OneD, int> &pPart);
};
typedef boost::shared_ptr<MeshPartition> MeshPartitionSharedPtr;
......
......@@ -39,99 +39,63 @@
#include <LibUtilities/BasicConst/NektarUnivTypeDefs.hpp>
#include <LibUtilities/BasicUtils/SharedArray.hpp>
#include "metis.h"
namespace Metis
{
extern "C"
{
// -- Sparse MAtrix Reordering (equivalent to onmetis)
void METIS_NodeND (int *nVerts, int *xadj, int *adjncy, int *numflag, int *options,
int *perm, int *iperm);
void AS_METIS_NodeND(int *nVerts, int *xadj, int *adjncy, int *numflag, int *options,
int *perm, int *iperm, int *map, int *mdswitch);
void METIS_PartMeshNodal(int *nElmts, int *nVerts, int *elType, int *numflag,
int *nparts, int *edgecut, int *epart, int *npart);
void METIS_PartGraphVKway(int *nVerts, int *xadj, int *adjcy, int *vertWgt, int *vertSize,
int *wgtFlag, int *numflag, int *nparts, int *options, int *volume, int *part);
}
//#ifdef NEKTAR_USING_METIS
inline static void onmetis(int *nVerts, int *xadj, int *adjncy, int *numflag, int *options,
int *perm, int *iperm)
{
METIS_NodeND(nVerts,xadj,adjncy,numflag,options,perm,iperm);
}
inline static void as_onmetis(int *nVerts, int *xadj, int *adjncy, int *numflag, int *options,
int *perm, int *iperm, int *map, int *mdswitch)
{
AS_METIS_NodeND(nVerts,xadj,adjncy,numflag,options,perm,iperm,map,mdswitch) ;
void AS_METIS_NodeND(int *nVerts, int *xadj, int *adjncy, int *vwgt,
int *options, int *perm, int *iperm, int *map,
int *mdswitch);
}
inline static void onmetis(int nVerts, Nektar::Array<Nektar::OneD, int> xadj, Nektar::Array<Nektar::OneD, int> adjncy,
Nektar::Array<Nektar::OneD, int> perm, Nektar::Array<Nektar::OneD, int> iperm)
inline static void as_onmetis(
int nVerts,
Nektar::Array<Nektar::OneD, int> xadj,
Nektar::Array<Nektar::OneD, int> adjncy,
Nektar::Array<Nektar::OneD, int> perm,
Nektar::Array<Nektar::OneD, int> iperm,
Nektar::Array<Nektar::OneD, int> map,
int mdswitch = 1)
{
ASSERTL1(xadj.num_elements() == nVerts+1,"Array xadj out of bounds");
ASSERTL1(perm.num_elements() == nVerts,"Array perm out of bounds");
ASSERTL1(iperm.num_elements() == nVerts,"Array iperm out of bounds");
int numflag = 0;
int options[8];
options[0]=0;
METIS_NodeND(&nVerts,&xadj[0],&adjncy[0],&numflag,options,&perm[0],&iperm[0]);
AS_METIS_NodeND(&nVerts, &xadj[0], &adjncy[0], NULL, NULL, &perm[0],
&iperm[0], &map[0], &mdswitch);
}
inline static void as_onmetis(int nVerts, Nektar::Array<Nektar::OneD, int> xadj, Nektar::Array<Nektar::OneD, int> adjncy,
Nektar::Array<Nektar::OneD, int> perm, Nektar::Array<Nektar::OneD, int> iperm, Nektar::Array<Nektar::OneD, int> map,
int mdswitch)
{
ASSERTL1(xadj.num_elements() == nVerts+1,"Array xadj out of bounds");
ASSERTL1(perm.num_elements() == nVerts,"Array perm out of bounds");
ASSERTL1(iperm.num_elements() == nVerts,"Array iperm out of bounds");
int numflag = 0;
int options[8];
options[0]=0;
AS_METIS_NodeND(&nVerts,&xadj[0],&adjncy[0],&numflag,options,&perm[0],&iperm[0],&map[0],&mdswitch);
}
inline static void MeshPartition(int nElmts, int nVerts, Nektar::Array<Nektar::OneD, int>& mesh, int type, int nparts,
Nektar::Array<Nektar::OneD, int>& edgePart, Nektar::Array<Nektar::OneD, int>& nodePart)
{
int numflag = 0;
METIS_PartMeshNodal(&nElmts, &nVerts, &mesh[0], &type, &numflag, &nparts, &edgePart[0], &nodePart[0]);
}
inline static void PartGraphVKway( int& nVerts,
Nektar::Array<Nektar::OneD, int>& xadj,
Nektar::Array<Nektar::OneD, int>& adjcy,
Nektar::Array<Nektar::OneD, int>& vertWgt,
Nektar::Array<Nektar::OneD, int>& vertSize,
int& nparts,
int& volume,
Nektar::Array<Nektar::OneD, int>& part)
inline static void PartGraphVKway(
int& nVerts,
Nektar::Array<Nektar::OneD, int>& xadj,
Nektar::Array<Nektar::OneD, int>& adjcy,
Nektar::Array<Nektar::OneD, int>& vertWgt,
Nektar::Array<Nektar::OneD, int>& vertSize,
int& nparts,
int& volume,
Nektar::Array<Nektar::OneD, int>& part)
{
int wgtflag = 0;
int *wgts = 0;
int *sizes = 0;
int *vwgt = 0;
int *vsize = 0;
if (vertWgt.num_elements() > 0)
{
wgtflag += 1;
wgts = &vertWgt[0];
vwgt = &vertWgt[0];
}
if (vertSize.num_elements() > 0)
{
wgtflag += 2;
sizes = &vertSize[0];
vsize = &vertSize[0];
}
int numflag = 0;
int options[5];
options[0]=0;
METIS_PartGraphVKway(&nVerts, &xadj[0], &adjcy[0], wgts, sizes, &wgtflag,
&numflag, &nparts, options, &volume, &part[0]);
int ncon = 1;
int options[METIS_NOPTIONS];
METIS_SetDefaultOptions(options);
METIS_PartGraphKway(&nVerts, &ncon, &xadj[0], &adjcy[0], vwgt, vsize,
0, &nparts, 0, 0, options, &volume, &part[0]);
}
//#endif //NEKTAR_USING_METIS
}
#endif //NEKTAR_LIB_UTILITIES_BASICUTILS_METIS_HPP
......@@ -684,7 +684,9 @@ namespace Nektar
boost::algorithm::trim(valueStr);
const parser_id parserID = location->value.id();
#if defined(NEKTAR_DEBUG) || defined(NEKTAR_FULLDEBUG)
const int num_children = location->children.size();
#endif
if (parserID == AnalyticExpression::constantID)
{
......
......@@ -62,6 +62,32 @@ namespace Nektar
const Array<OneD, const Array<OneD, NekDouble> > normals
= GetEdgeNormal(edge);
if (m_requireNeg.size() == 0)
{
m_requireNeg.resize(GetNedges());
for (int i = 0; i < GetNedges(); ++i)
{
m_requireNeg[i] = false;
if (m_negatedNormals[i])
{
m_requireNeg[i] = true;
}
Expansion1DSharedPtr edgeExp = boost::dynamic_pointer_cast<
Expansion1D>(m_edgeExp[i].lock());
if (edgeExp->GetRightAdjacentElementExp())
{
if (edgeExp->GetRightAdjacentElementExp()->GetGeom2D()
->GetGlobalID() == GetGeom2D()->GetGlobalID())
{
m_requireNeg[i] = true;
}
}
}
}
// We allow the case of mixed polynomial order by supporting only
// those modes on the edge common to both adjoining elements. This
// is enforced here by taking the minimum size and padding with
......@@ -79,6 +105,7 @@ namespace Nektar
boost::dynamic_pointer_cast<
LocalRegions::Expansion1D>(EdgeExp);
if (m_negatedNormals[edge])
{
Vmath::Neg(nquad_e,EdgeExp->UpdatePhys(),1);
......@@ -104,6 +131,32 @@ namespace Nektar
{
int i;
if (m_requireNeg.size() == 0)
{
m_requireNeg.resize(GetNedges());
for (i = 0; i < GetNedges(); ++i)
{
m_requireNeg[i] = false;
if (m_negatedNormals[i])
{
m_requireNeg[i] = true;
}
Expansion1DSharedPtr edgeExp = boost::dynamic_pointer_cast<
Expansion1D>(m_edgeExp[i].lock());
if (edgeExp->GetRightAdjacentElementExp())
{
if (edgeExp->GetRightAdjacentElementExp()->GetGeom2D()
->GetGlobalID() == GetGeom2D()->GetGlobalID())
{
m_requireNeg[i] = true;
}
}
}
}
StdRegions::Orientation edgedir = GetEorient(edge);
unsigned short num_mod0 = EdgeExp->GetBasis(0)->GetNumModes();
......@@ -126,18 +179,10 @@ namespace Nektar
boost::dynamic_pointer_cast<
LocalRegions::Expansion1D>(EdgeExp);
if (m_negatedNormals[edge])
if (m_requireNeg[edge])
{
Vmath::Neg(n_coeffs,EdgeExp->UpdateCoeffs(),1);
}
else if (locExp->GetRightAdjacentElementEdge() != -1)
{
if (locExp->GetRightAdjacentElementExp()->GetGeom2D()->GetGlobalID()
== GetGeom2D()->GetGlobalID())
{
Vmath::Neg(n_coeffs,EdgeExp->UpdateCoeffs(),1);
}
}
Array<OneD, NekDouble> coeff(n_coeffs,0.0);
LibUtilities::BasisType btype = ((LibUtilities::BasisType) 1); //1-->Ortho_A
......@@ -162,18 +207,10 @@ namespace Nektar
boost::dynamic_pointer_cast<
LocalRegions::Expansion1D>(EdgeExp);
if (m_negatedNormals[edge])
if (m_requireNeg[edge])
{
Vmath::Neg(n_coeffs,EdgeExp->UpdateCoeffs(),1);
}
else if (locExp->GetRightAdjacentElementEdge() != -1)
{
if (locExp->GetRightAdjacentElementExp()->GetGeom2D()->GetGlobalID()
== GetGeom2D()->GetGlobalID())
{
Vmath::Neg(order_e,EdgeExp->UpdateCoeffs(),1);
}
}
}
// add data to outarray if forward edge normal is outwards
......
......@@ -134,6 +134,7 @@ namespace Nektar
private:
std::vector<ExpansionWeakPtr> m_edgeExp;
std::vector<bool> m_requireNeg;
Expansion3DWeakPtr m_elementLeft;
Expansion3DWeakPtr m_elementRight;
......
......@@ -42,7 +42,7 @@ namespace Nektar
{
namespace LocalRegions
{
Expansion3D::Expansion3D(){}
Expansion3D::Expansion3D() : m_requireNeg() {}
void Expansion3D::AddHDGHelmholtzTraceTerms(
const NekDouble tau,
......@@ -914,6 +914,40 @@ namespace Nektar
StdRegions::IndexMapValuesSharedPtr map =
StdExpansion::GetIndexMap(ikey);
/*
* Coming into this routine, the velocity V will have been
* multiplied by the trace normals to give the input vector Vn. By
* convention, these normals are inwards facing for elements which
* have FaceExp as their right-adjacent face. This conditional
* statement therefore determines whether the normals must be
* negated, since the integral being performed here requires an
* outwards facing normal.
*/
if (m_requireNeg.size() == 0)
{
m_requireNeg.resize(GetNfaces());
for (i = 0; i < GetNfaces(); ++i)
{
m_requireNeg[i] = false;
if (m_negatedNormals[i])
{
m_requireNeg[i] = true;
}
Expansion2DSharedPtr faceExp = m_faceExp[i].lock();
if (faceExp->GetRightAdjacentElementExp())
{
if (faceExp->GetRightAdjacentElementExp()->GetGeom3D()->GetGlobalID()
== GetGeom3D()->GetGlobalID())
{
m_requireNeg[i] = true;
}
}
}
}
int order_e = (*map).num_elements(); // Order of the element
int n_coeffs = FaceExp->GetCoeffs().num_elements(); // Order of the trace
......@@ -924,33 +958,21 @@ namespace Nektar
else
{
FaceExp->IProductWRTBase(Fn,FaceExp->UpdateCoeffs());
LocalRegions::Expansion2DSharedPtr locExp =
boost::dynamic_pointer_cast<
LocalRegions::Expansion2D>(FaceExp);
/*
* Coming into this routine, the velocity V will have been
* multiplied by the trace normals to give the input vector
* Vn. By convention, these normals are inwards facing for
* elements which have FaceExp as their right-adjacent face.
* This conditional statement therefore determines whether the
* normals must be negated, since the integral being performed
* here requires an outwards facing normal.
*/
if (locExp->GetRightAdjacentElementFace() != -1)
}
if (m_requireNeg[face])
{
for(i = 0; i < order_e; ++i)
{
if (locExp->GetRightAdjacentElementExp()->GetGeom3D()->GetGlobalID()
== GetGeom3D()->GetGlobalID())
{
Vmath::Neg(order_e,FaceExp->UpdateCoeffs(),1);
}
outarray[(*map)[i].index] -= (*map)[i].sign*FaceExp->GetCoeff(i);
}
}
for(i = 0; i < order_e; ++i)
else
{
outarray[(*map)[i].index] += (*map)[i].sign*FaceExp->GetCoeff(i);
for(i = 0; i < order_e; ++i)
{
outarray[(*map)[i].index] += (*map)[i].sign*FaceExp->GetCoeff(i);
}
}
}
......
......@@ -106,6 +106,7 @@ namespace Nektar
// Only use this class for member functions
std::vector<Expansion2DWeakPtr> m_faceExp;
std::vector<bool> m_requireNeg;
};
// type defines for use of PrismExp in a boost vector
......
This diff is collapsed.
......@@ -162,7 +162,7 @@ namespace Nektar
LOCAL_REGIONS_EXPORT virtual
StdRegions::ExpansionType v_DetExpansionType() const;
LOCAL_REGIONS_EXPORT virtual const
LOCAL_REGIONS_EXPORT virtual const
SpatialDomains::GeomFactorsSharedPtr& v_GetMetricInfo() const;
LOCAL_REGIONS_EXPORT virtual const
......
......@@ -279,25 +279,22 @@ namespace Nektar
const Array<OneD, const NekDouble>& inarray,
Array<OneD, NekDouble>& outarray)
{
int nquad0 = m_base[0]->GetNumPoints();
int nquad1 = m_base[1]->GetNumPoints();
int nquad2 = m_base[2]->GetNumPoints();
Array<OneD, const NekDouble> jac = m_metricinfo->GetJac();
Array<OneD,NekDouble> tmp(nquad0*nquad1*nquad2);
const int nquad0 = m_base[0]->GetNumPoints();
const int nquad1 = m_base[1]->GetNumPoints();