Some smaller adaptivity updates. Added virtual method initBodyLoad and made two other methods in SIMbase virtual as well.

git-svn-id: http://svn.sintef.no/trondheim/IFEM/trunk@1187 e10b68d5-8a6e-419e-a041-bce267b0401d
This commit is contained in:
kmo 2011-09-23 08:06:32 +00:00 committed by Knut Morten Okstad
parent 61be192a55
commit e469c6caee
5 changed files with 37 additions and 23 deletions

View File

@ -17,6 +17,7 @@
#include "SIMenums.h"
#include "Utilities.h"
#include <sstream>
#include <cstdio>
AdaptiveSIM::AdaptiveSIM (SIMbase* sim) : model(sim)
@ -52,9 +53,11 @@ bool AdaptiveSIM::parse (char* keyWord, std::istream& is)
bool AdaptiveSIM::solveStep (const char* inputfile, SystemMatrix::Type solver,
int iStep)
{
std::cout <<"\nAdaptive step "<< iStep << std::endl;
if (iStep > 1)
{
// Re-generate the FE model after the refinement
model->clearProperties();
ASMunstruct::resetNumbering();
if (!model->read(inputfile) || !model->preprocess())
return false;
@ -73,13 +76,11 @@ bool AdaptiveSIM::solveStep (const char* inputfile, SystemMatrix::Type solver,
// Evaluate solution norms
model->setMode(SIM::RECOVERY);
return model->solutionNorms(Vectors(1,linsol),eNorm,gNorm);
return (model->solutionNorms(Vectors(1,linsol),eNorm,gNorm) &&
model->dumpResults(linsol,0.0,std::cout,true,6));
}
static bool larger (double a, double b) { return a > b; }
bool AdaptiveSIM::adaptMesh (int iStep)
{
printNorms(gNorm,std::cout);
@ -92,7 +93,8 @@ bool AdaptiveSIM::adaptMesh (int iStep)
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);
std::partial_sort(errors.begin(),errors.begin()+ipivot,
errors.end(),std::greater<double>());
double pivot = errors(ipivot);
std::cout <<"\nRefining "<< ipivot <<" elements with errors in range ["
@ -103,8 +105,8 @@ bool AdaptiveSIM::adaptMesh (int iStep)
if (eNorm(4,e) >= pivot)
elements.push_back(e-1);
char fname[10] = "mesh_.eps";
fname[4] = '0' + iStep;
char fname[12];
sprintf(fname,"mesh_%02d.eps",iStep);
return model->refine(elements,fname);
}
@ -118,5 +120,5 @@ void AdaptiveSIM::printNorms (const Vector& norms, std::ostream& os)
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;
os << std::endl;
}

View File

@ -28,6 +28,7 @@ struct Property
{
UNDEFINED,
MATERIAL,
BODYLOAD,
DIRICHLET,
DIRICHLET_INHOM,
NEUMANN

View File

@ -81,6 +81,7 @@ SIMbase::~SIMbase ()
for (FEModelVec::iterator i1 = myModel.begin(); i1 != myModel.end(); i1++)
delete *i1;
myModel.clear();
this->clearProperties();
#ifdef SP_DEBUG
@ -91,6 +92,8 @@ SIMbase::~SIMbase ()
void SIMbase::clearProperties ()
{
for (FEModelVec::iterator i1 = myModel.begin(); i1 != myModel.end(); i1++)
(*i1)->clear(true); // retain the geometry only
for (SclFuncMap::iterator i2 = myScalars.begin(); i2 != myScalars.end(); i2++)
delete i2->second;
for (VecFuncMap::iterator i3 = myVectors.begin(); i3 != myVectors.end(); i3++)
@ -173,6 +176,9 @@ bool SIMbase::parse (char* keyWord, std::istream& is)
return !cline.fail() && !cline.bad();
}
else if (!strncasecmp(keyWord,"ADAPTIVE",8))
;
else if (isalpha(keyWord[0]))
std::cerr <<" *** SIMbase::parse: Unknown keyword: "<< keyWord << std::endl;
@ -221,13 +227,13 @@ int SIMbase::getLocalPatchIndex (int patchNo) const
}
bool SIMbase::preprocess (const std::vector<int>& ignoredPatches, bool fixDup)
bool SIMbase::preprocess (const std::vector<int>& ignored, bool fixDup)
{
if (!this->createFEMmodel()) return false;
// Erase all patches that should be ignored in the analysis
std::vector<int>::const_iterator it;
for (it = ignoredPatches.begin(); it != ignoredPatches.end(); it++)
for (it = ignored.begin(); it != ignored.end(); it++)
if (*it > 0 && (size_t)*it <= myModel.size())
myModel[*it-1]->clear();
@ -524,6 +530,7 @@ bool SIMbase::assembleSystem (const TimeDomain& time, const Vectors& prevSol,
if (msgLevel > 1)
std::cout <<"\nAssembling interior matrix terms for P"<< j
<< std::endl;
this->initBodyLoad(j);
this->extractPatchSolution(myModel[j-1],prevSol);
ok = myModel[j-1]->integrate(*myProblem,*myEqSys,time);
lp = j;
@ -539,6 +546,7 @@ bool SIMbase::assembleSystem (const TimeDomain& time, const Vectors& prevSol,
if (msgLevel > 1)
std::cout <<"\nAssembling interior matrix terms for P"<< i+1
<< std::endl;
this->initBodyLoad(i+1);
this->extractPatchSolution(myModel[i],prevSol);
ok = myModel[i]->integrate(*myProblem,*myEqSys,time);
lp = i+1;

View File

@ -100,10 +100,10 @@ public:
const char* = 0) { return false; }
//! \brief Performs some pre-processing tasks on the FE model.
//! \param[in] ignoredPatches Indices of patches to ignore in the analysis
//! \param[in] ignored Indices of patches to ignore in the analysis
//! \param[in] fixDup Merge duplicated FE nodes on patch interfaces?
bool preprocess(const std::vector<int>& ignoredPatches = std::vector<int>(),
bool fixDup = false);
virtual bool preprocess(const std::vector<int>& ignored = std::vector<int>(),
bool fixDup = false);
//! \brief Allocates the system matrices of the FE problem to be solved.
//! \param[in] mType The matrix format to use
@ -152,7 +152,7 @@ public:
int getNoPatches() const { return nGlPatches; }
//! \brief Returns the number of registered result points.
size_t getNoResultPoints() const { return myPoints.size(); }
//! \brief Return the visualization dump interval
//! \brief Returns the visualization dump interval.
int getDumpInterval() const { return vizIncr; }
//! \brief Initializes time-dependent in-homogeneous Dirichlet coefficients.
@ -490,7 +490,8 @@ protected:
//! \brief Initializes material properties for integration of interior terms.
virtual bool initMaterial(size_t) { return true; }
//! \brief Initializes the body load properties for current patch.
virtual bool initBodyLoad(size_t) { return true; }
//! \brief Initializes for integration of Neumann terms for a given property.
virtual bool initNeumann(size_t) { return true; }
@ -498,8 +499,8 @@ protected:
//! \details This method helps with encapsulating PETSc in libIFEM.
void readLinSolParams(std::istream& is, int npar);
//! \brief Finalizes the system assembly.
bool finalizeAssembly(bool newLHSmatrix);
//! \brief Finalizes the global equation system assembly.
virtual bool finalizeAssembly(bool newLHSmatrix);
public:
//! \brief Enum defining the available discretization methods.
@ -529,16 +530,17 @@ protected:
TracFuncMap myTracs; //!< Traction property fields
IntegrandBase* myProblem; //!< Problem-specific data and methods
AnaSol* mySol; //!< Analytical/Exact solution
VTF* myVtf; //!< VTF-file for result visualization
ResPointVec myPoints; //!< User-defined result sampling points
// Post-processing attributes
ResPointVec myPoints; //!< User-defined result sampling points
VTF* myVtf; //!< VTF-file for result visualization
int vizIncr; //!< Number increments between each result output
int format; //!< Output file format
// Parallel computing attributes
int nGlPatches; //!< Number of global patches
std::vector<int> myPatches; //!< Global patch numbers for current processor
// Visualization attributes
int vizIncr;
int format;
// Equation solver attributes
AlgEqSystem* myEqSys; //!< The actual linear equation system

View File

@ -26,10 +26,11 @@ class SIMinput
protected:
//! \brief The default constructor initializes \a myPid and \a nProc.
SIMinput();
public:
//! \brief Empty destructor.
virtual ~SIMinput() {}
public:
//! \brief Reads model data from the specified input file \a *fileName.
bool read(const char* fileName);
//! \brief Parses a data section from an input stream.