OilfieldToolsNet/IFEM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Added adaptive simulation driver

Browse Source 
git-svn-id: http://svn.sintef.no/trondheim/IFEM/trunk@1185 e10b68d5-8a6e-419e-a041-bce267b0401d
This commit is contained in:
kmo
2011-09-22 14:54:12 +00:00
committed by Knut Morten Okstad
parent bcc84da438
commit b99db60741
8 changed files with 440 additions and 154 deletions
Download Patch File Download Diff File Expand all files Collapse all files

65
src/ASM/ASM2D.h Normal file
View File

@@ -0,0 +1,65 @@
// $Id$
//==============================================================================
//!
//! \file ASM2D.h
//!
//! \date Sep 20 2011
//!
//! \author Knut Morten Okstad / SINTEF
//!
//! \brief Abstract interface for 2D patches.
//!
//==============================================================================
#ifndef _ASM_2D_H
#define _ASM_2D_H
/*!
\brief Abstract interface for 2D spline patches.
\details This class contains an interface to methods common for structured and
unstructured 2D patches, such that these methods can be invoked without need
to type-cast the patch object to the actual class type.
*/
class ASM2D
{
protected:
//! \brief The constructor is protected to allow objects of sub-classes only.
ASM2D() {}
public:
//! \brief Empty destructor.
virtual ~ASM2D() {}
//! \brief Constrains all DOFs on a given boundary edge.
//! \param[in] dir Parameter direction defining the edge to constrain
//! \param[in] dof Which DOFs to constrain at each node on the edge
//! \param[in] code Inhomogeneous dirichlet condition code
virtual void constrainEdge(int dir, int dof = 123, int code = 0) = 0;
//! \brief Constrains a corner node identified by the two parameter indices.
//! \param[in] I Parameter index in u-direction
//! \param[in] J Parameter index in v-direction
//! \param[in] dof Which DOFs to constrain at the node
//! \param[in] code Inhomogeneous dirichlet condition code
//!
//! \details The sign of the two indices is used to define whether we want
//! the node at the beginning or the end of that parameter direction.
//! The magnitude of the indices are not used.
virtual void constrainCorner(int I, int J, int dof = 123, int code = 0) = 0;
//! \brief Constrains a node identified by two relative parameter values.
//! \param[in] xi Parameter in u-direction
//! \param[in] eta Parameter in v-direction
//! \param[in] dof Which DOFs to constrain at the node
//! \param[in] code Inhomogeneous dirichlet condition code
//!
//! \details The parameter values have to be in the domain [0.0,1.0], where
//! 0.0 means the beginning of the domain and 1.0 means the end. For values
//! in between, the actual index is taken as the integer value closest to
//! \a r*n, where \a r denotes the given relative parameter value,
//! and \a n is the number of nodes along that parameter direction.
virtual void constrainNode(double xi, double eta,
int dof = 123, int code = 0) = 0;
};
#endif

10
src/ASM/ASMs2D.h
View File

@@ -15,6 +15,7 @@
#define _ASM_S2D_H
#include "ASMstruct.h"
#include "ASM2D.h"
namespace Go {
class SplineSurface;
@@ -28,7 +29,7 @@ namespace Go {
\details This class contains methods common for structured 2D spline patches.
*/
class ASMs2D : public ASMstruct
class ASMs2D : public ASMstruct, public ASM2D
{
//! \brief Struct for nodal point data.
struct IJ
@@ -134,7 +135,7 @@ public:
//! \param[in] dir Parameter direction defining the edge to constrain
//! \param[in] dof Which DOFs to constrain at each node on the edge
//! \param[in] code Inhomogeneous dirichlet condition code
void constrainEdge(int dir, int dof = 123, int code = 0);
virtual void constrainEdge(int dir, int dof = 123, int code = 0);
//! \brief Constrains a corner node identified by the two parameter indices.
//! \param[in] I Parameter index in u-direction
@@ -145,7 +146,7 @@ public:
//! \details The sign of the two indices is used to define whether we want
//! the node at the beginning or the end of that parameter direction.
//! The magnitude of the indices are not used.
void constrainCorner(int I, int J, int dof = 123, int code = 0);
virtual void constrainCorner(int I, int J, int dof = 123, int code = 0);
//! \brief Constrains a node identified by two relative parameter values.
//! \param[in] xi Parameter in u-direction
//! \param[in] eta Parameter in v-direction
@@ -157,7 +158,8 @@ public:
//! in between, the actual index is taken as the integer value closest to
//! \a r*n, where \a r denotes the given relative parameter value,
//! and \a n is the number of nodes along that parameter direction.
void constrainNode(double xi, double eta, int dof = 123, int code = 0);
virtual void constrainNode(double xi, double eta,
int dof = 123, int code = 0);
//! \brief Connects all matching nodes on two adjacent boundary edges.
//! \param[in] edge Local edge index of this patch, in range [1,4]

31
src/ASM/LR/ASMu2D.C
View File

@@ -18,7 +18,6 @@
#include "LRSpline/LRSplineSurface.h"
#include "LRSpline/Element.h"
#include "LRSpline/Basisfunction.h"
#include <fstream>
#include "ASMu2D.h"
#include "TimeDomain.h"
@@ -32,6 +31,8 @@
#include "Profiler.h"
#include "Vec3Oper.h"
#include <ctype.h>
#include <fstream>
ASMu2D::ASMu2D (const char* fName, unsigned char n_s, unsigned char n_f)
: ASMunstruct(2,n_s,n_f), lrspline(0), tensorspline(0)
@@ -120,19 +121,23 @@ bool ASMu2D::write (std::ostream& os, int) const
}
void ASMu2D::clear ()
void ASMu2D::clear (bool retainGeometry)
{
if (!retainGeometry) {
// Erase spline data
if (lrspline) delete lrspline;
if (tensorspline) delete tensorspline;
lrspline = 0;
tensorspline = 0;
geo = 0;
}
// Erase the FE data
ASMbase::clear();
// Erase the FE data
this->ASMbase::clear(retainGeometry);
}
bool ASMu2D::cornerRefine (int minBasisfunctions)
{
if (!lrspline ) return false;
@@ -289,6 +294,20 @@ bool ASMu2D::raiseOrder (int ru, int rv)
}
bool ASMu2D::refine (const std::vector<int>& elements, const char* fName)
{
if (!lrspline) return false;
lrspline->refineElement(elements);
if (fName)
{
std::ofstream meshFile(fName);
lrspline->writePostscriptMesh(meshFile);
}
return true;
}
bool ASMu2D::generateFEMTopology ()
{
if (!lrspline) return false;
@@ -533,9 +552,9 @@ void ASMu2D::constrainCorner (int I, int J, int dof, int code)
// Hopefully we don't have to constrain non-corner singlenodes inside patches
#if 0
void ASMu2D::constrainNode (double xi, double eta, int dof, int code)
{
#if 0
if (xi < 0.0 || xi > 1.0) return;
if (eta < 0.0 || eta > 1.0) return;
@@ -547,8 +566,8 @@ void ASMu2D::constrainNode (double xi, double eta, int dof, int code)
if (eta > 0.0) node += n1*int(0.5+(n2-1)*eta);
this->prescribe(node,dof,code);
}
#endif
}
#define DERR -999.99

14
src/ASM/LR/ASMu2D.h
View File

@@ -15,6 +15,7 @@
#define _ASM_U2D_H
#include "ASMunstruct.h"
#include "ASM2D.h"
namespace Go {
class SplineSurface;
@@ -27,11 +28,11 @@ namespace LR {
/*!
\brief Driver for assembly of structured 2D spline FE models.
\details This class contains methods common for structured 2D spline patches.
\brief Driver for assembly of unstructured 2D spline FE models.
\details This class contains methods common for 2D LR-spline patches.
*/
class ASMu2D : public ASMunstruct
class ASMu2D : public ASMunstruct, public ASM2D
{
//! \brief Struct for nodal point data.
struct IJ
@@ -98,7 +99,9 @@ public:
virtual bool generateFEMTopology();
//! \brief Clears the contents of the patch, making it empty.
virtual void clear();
//! \param[in] retainGeometry If \e true, the spline geometry is not cleared.
//! This is used to reinitialize the patch after it has been refined.
virtual void clear(bool retainGeometry = false);
//! \brief Returns the global coordinates for the given node.
//! \param[in] inod 1-based node index local to current patch
@@ -155,7 +158,8 @@ public:
//! \param[in] ru Number of times to raise the order in u-direction
//! \param[in] rv Number of times to raise the order in v-direction
bool raiseOrder(int ru, int rv);
//! \brief Refines the specified list of elements.
bool refine(const std::vector<int>& elements, const char* fName = 0);
// Various methods for preprocessing of boundary conditions and patch topology

122
src/SIM/AdaptiveSIM.C Normal file
View File

@@ -0,0 +1,122 @@
// $Id$
//==============================================================================
//!
//! \file AdaptiveSIM.h
//!
//! \date Sep 22 2011
//!
//! \author Knut Morten Okstad / SINTEF
//!
//! \brief Adaptive solution driver for isogeometric FEM simulators.
//!
//==============================================================================
#include "AdaptiveSIM.h"
#include "ASMunstruct.h"
#include "SIMbase.h"
#include "SIMenums.h"
#include "Utilities.h"
#include <sstream>
AdaptiveSIM::AdaptiveSIM (SIMbase* sim) : model(sim)
{
// Default adaptation parameters
nStep = 10;
stopTol = 1.0;
beta = 25.0;
}
AdaptiveSIM::~AdaptiveSIM ()
{
if (model) delete model;
}
bool AdaptiveSIM::parse (char* keyWord, std::istream& is)
{
if (!strncasecmp(keyWord,"ADAPTIVE",8))
{
std::istringstream cline(utl::readLine(is));
cline >> nStep >> stopTol >> beta;
if (cline.fail() || cline.bad()) return false;
}
else
return model->parse(keyWord,is);
return true;
}
bool AdaptiveSIM::solveStep (const char* inputfile, SystemMatrix::Type solver,
int iStep)
{
if (iStep > 1)
{
// Re-generate the FE model after the refinement
ASMunstruct::resetNumbering();
if (!model->read(inputfile) || !model->preprocess())
return false;
}
// Assemble the linear FE equation system
model->setMode(SIM::STATIC,true);
model->initSystem(iStep == 1 ? SystemMatrix::DENSE : solver, 1,1);
model->setAssociatedRHS(0,0);
if (!model->assembleSystem())
return false;
// Solve the linear system of equations
if (!model->solveSystem(linsol,1))
return false;
// Evaluate solution norms
model->setMode(SIM::RECOVERY);
return model->solutionNorms(Vectors(1,linsol),eNorm,gNorm);
}
static bool larger (double a, double b) { return a > b; }
bool AdaptiveSIM::adaptMesh (int iStep)
{
printNorms(gNorm,std::cout);
double eta = gNorm.size() > 3 ? 100.0*gNorm(4)/gNorm(3) : 0.0;
if (eta < stopTol || iStep > nStep) return false;
// Find the list of elements to refine (the beta % with the highest error)
std::vector<int> elements;
Vector errors(eNorm.getRow(4));
size_t ipivot = ceil(errors.size()*beta/100.0);
if (ipivot < 1 || ipivot > errors.size()) return false;
std::partial_sort(errors.begin(),errors.begin()+ipivot,errors.end(),larger);
double pivot = errors(ipivot);
std::cout <<"\nRefining "<< ipivot <<" elements with errors in range ["
<< pivot <<","<< errors.front() <<"]"<< std::endl;
elements.reserve(ipivot);
for (size_t e = 1; e <= errors.size(); e++)
if (eNorm(4,e) >= pivot)
elements.push_back(e-1);
char fname[10] = "mesh_.eps";
fname[4] = '0' + iStep;
return model->refine(elements,fname);
}
void AdaptiveSIM::printNorms (const Vector& norms, std::ostream& os)
{
os <<"Energy norm |u^h| = a(u^h,u^h)^0.5 : "<< norms(1)
<<"\nExternal energy ((h,u^h)+(t,u^h)^0.5 : "<< norms(2);
if (norms.size() > 2)
os <<"\nExact norm |u| = a(u,u)^0.5 : "<< norms(3);
if (norms.size() > 3)
os <<"\nExact error a(e,e)^0.5, e=u-u^h : "<< norms(4)
<<"\nExact relative error (%) : "<< 100.0*norms(4)/norms(3);
std::cout << std::endl;
}

78
src/SIM/AdaptiveSIM.h Normal file
View File

@@ -0,0 +1,78 @@
// $Id$
//==============================================================================
//!
//! \file AdaptiveSIM.h
//!
//! \date Sep 22 2011
//!
//! \author Knut Morten Okstad / SINTEF
//!
//! \brief Adaptive solution driver for isogeometric FEM simulators.
//!
//==============================================================================
#ifndef _ADAPTIVE_SIM_H
#define _ADAPTIVE_SIM_H
#include "SIMinput.h"
#include "SystemMatrix.h"
class SIMbase;
/*!
\brief Nonlinear solution driver for isogeometric FEM simulators.
\details This class contains data and methods for computing the nonlinear
solution to a FE problem based on splines/NURBS basis functions, through
Newton-Raphson iterations.
*/
class AdaptiveSIM : public SIMinput
{
public:
//! \brief The constructor initialized default solution parameters.
//! \param sim Pointer to the spline FE model
AdaptiveSIM(SIMbase* sim = 0);
//! \brief The destructor frees the dynamically allocated FE model object.
virtual ~AdaptiveSIM();
//! \brief Solves the nonlinear equations by Newton-Raphson iterations.
//! \param param Solution algorithm parameters
//! \param[in] mode Solution mode to use for this step
//! \param[in] compName Solution name to be used in the norm output
//! \param[in] energyNorm If \e true, integrate energy norm of the solution
//! \param[in] zero_tolerance Truncate norm values small than this to zero
//! \param[in] outPrec Number of digits after the decimal point in norm print
bool solveStep(const char* inputfile, SystemMatrix::Type solver, int iStep);
//! \brief Computes and prints some solution norm quantities.
//! \param[in] time Parameters for nonlinear/time-dependent simulations
//! \param[in] compName Solution name to be used in the norm output
//! \param[in] energyNorm If \e true, integrate energy norm of the solution
//! \param[in] zero_tolerance Truncate norm values small than this to zero
//! \param[in] outPrec Number of digits after the decimal point in norm print
bool adaptMesh(int iStep);
//! \brief Prints out the global norms to given stream
static void printNorms(const Vector& norms, std::ostream& os);
protected:
//! \brief Parses a data section from an input stream.
//! \param[in] keyWord Keyword of current data section to read
//! \param is The file stream to read from
virtual bool parse(char* keyWord, std::istream& is);
private:
SIMbase* model; //!< The isogeometric FE model
int nStep;
double stopTol;
double beta;
Vector linsol; //!< Linear solution vector
Vector gNorm; //!< Global norms
Matrix eNorm; //!< Element norms
};
#endif

256
src/SIM/SIM2D.C
View File

@@ -11,13 +11,13 @@
//!
//==============================================================================
#ifdef HAS_LRSPLINE
#include "LR/ASMu2D.h"
#endif
#include "SIM2D.h"
#include "ASMs2Dmx.h"
#include "ASMs2DmxLag.h"
#include "ASMs2DSpec.h"
#ifdef HAS_LRSPLINE
#include "LR/ASMu2D.h"
#endif
#include "Functions.h"
#include "Utilities.h"
#include <fstream>
@@ -38,59 +38,39 @@ bool SIM2D::parse (char* keyWord, std::istream& is)
char* cline = 0;
if (!strncasecmp(keyWord,"PATCHES",7))
{
ASMbase* pch = 0;
int npatch = atoi(keyWord+7);
std::cout <<"\nNumber of patches: "<< npatch << std::endl;
for (int i = 0; i < npatch && (cline = utl::readLine(is)); i++)
if (myModel.empty())
{
cline = strtok(cline," ");
switch (discretization) {
case Lagrange:
if (nf[1] > 0)
pch = new ASMs2DmxLag(cline,2,nf[0],nf[1]);
else
pch = new ASMs2DLag(cline,2,nf[0]);
break;
case Spectral:
pch = new ASMs2DSpec(cline,2,nf[0]);
break;
#ifdef HAS_LRSPLINE
case LRSpline:
pch = new ASMu2D(cline,2,nf[0]);
break;
#endif
default:
if (nf[1] > 0)
pch = new ASMs2Dmx(cline,2,nf[0],nf[1]);
else
pch = new ASMs2D(cline,2,nf[0]);
}
if (pch->empty() || this->getLocalPatchIndex(i+1) < 1)
delete pch;
else
myModel.push_back(pch);
}
std::cout <<"\nNumber of patches: "<< npatch << std::endl;
for (int i = 0; i < npatch && (cline = utl::readLine(is)); i++)
this->readPatch(strtok(cline," "),i);
if ((int)myModel.size() < npatch)
{
std::cerr <<" *** SIM2D::parse: Expected "<< npatch
<<" patches but could read only "<< myModel.size()
<< std::endl;
return false;
if ((int)myModel.size() < npatch)
{
std::cerr <<" *** SIM2D::parse: Expected "<< npatch
<<" patches but could read only "<< myModel.size()
<< std::endl;
return false;
}
}
else // just read through the npatch next lines without doing anything
for (int i = 0; i < npatch && utl::readLine(is); i++);
}
else if (!strncasecmp(keyWord,"PATCHFILE",9))
{
size_t i = 9; while (i < strlen(keyWord) && isspace(keyWord[i])) i++;
std::cout <<"\nReading data file "<< keyWord+i << std::endl;
std::ifstream isp(keyWord+i);
this->readPatches(isp);
if (myModel.empty())
{
std::cerr <<" *** SIM2D::parse: No patches read"<< std::endl;
return false;
size_t i = 9; while (i < strlen(keyWord) && isspace(keyWord[i])) i++;
std::cout <<"\nReading data file "<< keyWord+i << std::endl;
std::ifstream isp(keyWord+i);
this->readPatches(isp);
if (myModel.empty())
{
std::cerr <<" *** SIM2D::parse: No patches read"<< std::endl;
return false;
}
}
}
@@ -465,76 +445,40 @@ bool SIM2D::addConstraint (int patch, int lndx, int ldim, int dirs, int code)
<< (ldim == 0 ? " V" : " E") << lndx <<" in direction(s) "<< dirs;
if (code) std::cout <<" code = "<< code <<" ";
if(discretization == LRSpline ) {
#ifdef HAS_LRSPLINE
ASMu2D* pch = static_cast<ASMu2D*>(myModel[patch-1]);
switch (ldim)
{
case 0: // Vertex constraints
switch (lndx)
{
case 1: pch->constrainCorner(-1,-1,dirs,code); break;
case 2: pch->constrainCorner( 1,-1,dirs,code); break;
case 3: pch->constrainCorner(-1, 1,dirs,code); break;
case 4: pch->constrainCorner( 1, 1,dirs,code); break;
default: std::cout << std::endl;
return constrError("vertex index ",lndx);
}
break;
case 1: // Edge constraints
switch (lndx)
{
case 1: pch->constrainEdge(-1,dirs,code); break;
case 2: pch->constrainEdge( 1,dirs,code); break;
case 3: pch->constrainEdge(-2,dirs,code); break;
case 4: pch->constrainEdge( 2,dirs,code); break;
default: std::cout << std::endl;
return constrError("edge index ",lndx);
}
break;
default:
std::cout << std::endl;
return constrError("local dimension switch ",ldim);
}
#else
std::cerr << "Error: No LR-spline library detected\n";
return constrError("discretization: LR-spline ", discretization);
#endif
} else {
ASMs2D* pch = static_cast<ASMs2D*>(myModel[patch-1]);
switch (ldim)
{
case 0: // Vertex constraints
switch (lndx)
{
case 1: pch->constrainCorner(-1,-1,dirs,code); break;
case 2: pch->constrainCorner( 1,-1,dirs,code); break;
case 3: pch->constrainCorner(-1, 1,dirs,code); break;
case 4: pch->constrainCorner( 1, 1,dirs,code); break;
default: std::cout << std::endl;
return constrError("vertex index ",lndx);
}
break;
case 1: // Edge constraints
switch (lndx)
{
case 1: pch->constrainEdge(-1,dirs,code); break;
case 2: pch->constrainEdge( 1,dirs,code); break;
case 3: pch->constrainEdge(-2,dirs,code); break;
case 4: pch->constrainEdge( 2,dirs,code); break;
default: std::cout << std::endl;
return constrError("edge index ",lndx);
}
break;
default:
std::cout << std::endl;
return constrError("local dimension switch ",ldim);
}
}
// Must dynamic cast here, since ASM2D is not derived from ASMbase
ASM2D* pch = dynamic_cast<ASM2D*>(myModel[patch-1]);
if (!pch) return constrError("2D patch ",patch);
switch (ldim)
{
case 0: // Vertex constraints
switch (lndx)
{
case 1: pch->constrainCorner(-1,-1,dirs,code); break;
case 2: pch->constrainCorner( 1,-1,dirs,code); break;
case 3: pch->constrainCorner(-1, 1,dirs,code); break;
case 4: pch->constrainCorner( 1, 1,dirs,code); break;
default: std::cout << std::endl;
return constrError("vertex index ",lndx);
}
break;
case 1: // Edge constraints
switch (lndx)
{
case 1: pch->constrainEdge(-1,dirs,code); break;
case 2: pch->constrainEdge( 1,dirs,code); break;
case 3: pch->constrainEdge(-2,dirs,code); break;
case 4: pch->constrainEdge( 2,dirs,code); break;
default: std::cout << std::endl;
return constrError("edge index ",lndx);
}
break;
default:
std::cout << std::endl;
return constrError("local dimension switch ",ldim);
}
return true;
}
@@ -548,6 +492,37 @@ void SIM2D::setQuadratureRule (size_t ng)
}
void SIM2D::readPatch (const char* patchFile, int pchInd)
{
ASMbase* pch = 0;
switch (discretization) {
case Lagrange:
if (nf[1] > 0)
pch = new ASMs2DmxLag(patchFile,2,nf[0],nf[1]);
else
pch = new ASMs2DLag(patchFile,2,nf[0]);
break;
case Spectral:
pch = new ASMs2DSpec(patchFile,2,nf[0]);
break;
#ifdef HAS_LRSPLINE
case LRSpline:
pch = new ASMu2D(patchFile,2,nf[0]);
break;
#endif
default:
if (nf[1] > 0)
pch = new ASMs2Dmx(patchFile,2,nf[0],nf[1]);
else
pch = new ASMs2D(patchFile,2,nf[0]);
}
if (pch->empty() || this->getLocalPatchIndex(pchInd+1) < 1)
delete pch;
else
myModel.push_back(pch);
}
void SIM2D::readPatches (std::istream& isp)
{
ASMbase* pch = 0;
@@ -555,25 +530,25 @@ void SIM2D::readPatches (std::istream& isp)
{
std::cout <<"Reading patch "<< patchNo << std::endl;
switch (discretization) {
case Lagrange:
if (nf[1] > 0)
pch = new ASMs2DmxLag(isp,2,nf[0],nf[1]);
else
pch = new ASMs2DLag(isp,2,nf[0]);
break;
case Spectral:
pch = new ASMs2DSpec(isp,2,nf[0]);
break;
case LRSpline:
#ifdef HAS_LRSPLINE
pch = new ASMu2D(isp,2,nf[0]);
#endif
break;
default:
if (nf[1] > 0)
pch = new ASMs2Dmx(isp,2,nf[0],nf[1]);
else
pch = new ASMs2D(isp,2,nf[0]);
case Lagrange:
if (nf[1] > 0)
pch = new ASMs2DmxLag(isp,2,nf[0],nf[1]);
else
pch = new ASMs2DLag(isp,2,nf[0]);
break;
case Spectral:
pch = new ASMs2DSpec(isp,2,nf[0]);
break;
case LRSpline:
#ifdef HAS_LRSPLINE
pch = new ASMu2D(isp,2,nf[0]);
break;
#endif
default:
if (nf[1] > 0)
pch = new ASMs2Dmx(isp,2,nf[0],nf[1]);
else
pch = new ASMs2D(isp,2,nf[0]);
}
if (pch->empty() || this->getLocalPatchIndex(patchNo) < 1)
delete pch;
@@ -581,3 +556,14 @@ void SIM2D::readPatches (std::istream& isp)
myModel.push_back(pch);
}
}
bool SIM2D::refine (const std::vector<int>& elements, const char* fName)
{
for (size_t i = 0; i < myModel.size(); i++)
if (!myModel.empty())
if (!static_cast<ASMu2D*>(myModel[i])->refine(elements,fName))
return false;
return true;
}

18
src/SIM/SIM2D.h
View File

@@ -37,16 +37,26 @@ public:
//! \param[in] ng Number of Gauss points in each parameter direction
virtual void setQuadratureRule(size_t ng);
//! \brief Read patches from given stream
//! \param[in] isp The stream to read from
void readPatches(std::istream& isp);
protected:
//! \brief Parses a data section from an input stream.
//! \param[in] keyWord Keyword of current data section to read
//! \param is The file stream to read from
virtual bool parse(char* keyWord, std::istream& is);
//! \brief Reads a patch from given input file.
//! \param[in] patchFile Name of file to read from
//! \param[in] pchInd 0-based index of the patch to read
void readPatch(const char* patchFile, int pchInd);
//! \brief Reads patches from given input stream.
//! \param[in] isp The file stream to read from
void readPatches(std::istream& isp);
//! \brief Refines a list of elements.
//! \param[in] elements 1-based indices of the elements to refine
//! \param[in] fName Optional mesh output file (Encapsulated PostScript)
virtual bool refine(const std::vector<int>& elements,
const char* fName = 0);
//! \brief Preprocesses a user-defined Dirichlet boundary property.
//! \param[in] patch 1-based index of the patch to receive the property
//! \param[in] lndx Local index of the boundary item to receive the property