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

Use utl::getValue to extract nonlinear solution parameters from the XML-tags.

Browse Source 
Use nGauss[1] when integrating energy norms also for nonlinear problems.

git-svn-id: http://svn.sintef.no/trondheim/IFEM/trunk@1552 e10b68d5-8a6e-419e-a041-bce267b0401d
This commit is contained in:
kmo 2012-03-27 12:44:45 +00:00 committed by Knut Morten Okstad
parent 13412d0a8a
commit a849f1fa17
1 changed files with 46 additions and 42 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

88
src/SIM/NonLinSIM.C
View File

@ -15,9 +15,8 @@
#include "SIMbase.h" #include "SIMbase.h"
#include "Profiler.h" #include "Profiler.h"
#include "Utilities.h" #include "Utilities.h"
#include <sstream>
#include "tinyxml.h" #include "tinyxml.h"
#include <sstream>
NonLinSIM::NonLinSIM (SIMbase* sim) : model(sim), nBlock(0) NonLinSIM::NonLinSIM (SIMbase* sim) : model(sim), nBlock(0)
@ -105,70 +104,71 @@ bool NonLinSIM::parse (char* keyWord, std::istream& is)
} }
void NonLinSIM::initSystem (int mType, size_t nGauss)
{
model->initSystem(mType,1,1);
model->setAssociatedRHS(0,0);
model->setQuadratureRule(nGauss);
}
bool NonLinSIM::parse (const TiXmlElement* elem) bool NonLinSIM::parse (const TiXmlElement* elem)
{ {
if (strcasecmp(elem->Value(),"nonlinearsolver")) if (strcasecmp(elem->Value(),"nonlinearsolver"))
return model->parse(elem); return model->parse(elem);
const char* value = 0;
const TiXmlElement* child = elem->FirstChildElement(); const TiXmlElement* child = elem->FirstChildElement();
while (child) { for (; child; child = child->NextSiblingElement())
if (!strcasecmp(child->Value(),"timestepping")) { if (!strcasecmp(child->Value(),"timestepping")) {
const TiXmlElement* step = child->FirstChildElement("step");
steps.clear(); steps.clear();
while (step) { const TiXmlElement* step = child->FirstChildElement("step");
for (; step; step = step->NextSiblingElement()) {
double start = 0.0, end = 0.0, dt = 0.0; double start = 0.0, end = 0.0, dt = 0.0;
std::pair<std::vector<double>,double> tstep; std::pair<std::vector<double>,double> timeStep;
utl::getAttribute(step,"start",start); utl::getAttribute(step,"start",start);
utl::getAttribute(step,"end",end); utl::getAttribute(step,"end",end);
if (steps.empty()) if (steps.empty()) startTime = start;
startTime = start; if (step->FirstChild()) {
if (step->FirstChild() && step->FirstChild()->Value()) { std::istringstream cline(step->FirstChild()->Value());
std::istringstream cline(child->FirstChild()->Value());
cline >> dt; cline >> dt;
if (dt > 1.0 && ceil(dt) == dt) { // number of steps specified if (dt > 1.0 && ceil(dt) == dt) {
dt = (end-steps.empty()?start:steps.back().second)/dt; // The number of steps are specified
tstep.first.push_back(dt); dt = (end - (steps.empty() ? start : steps.back().second))/dt;
} else while (!cline.fail() && !cline.bad()) { // steps specified timeStep.first.push_back(dt);
tstep.first.push_back(dt); }
else while (!cline.fail() && !cline.bad()) {
// The time step size(s) is/are specified
timeStep.first.push_back(dt);
cline >> dt; cline >> dt;
} }
tstep.second = end; timeStep.second = end;
steps.push_back(tstep); steps.push_back(timeStep);
} }
} }
stopTime = steps.back().second; stopTime = steps.back().second;
} else if (!strcasecmp(child->Value(),"maxits")) {
if (child->FirstChild() && child->FirstChild()->Value())
maxit = atoi(child->FirstChild()->Value());
} else if (!strcasecmp(child->Value(),"nupdate")) {
if (child->FirstChild() && child->FirstChild()->Value())
nupdat = atoi(child->FirstChild()->Value());
} else if (!strcasecmp(child->Value(),"rtol")) {
if (child->FirstChild() && child->FirstChild()->Value())
convTol = atof(child->FirstChild()->Value());
} else if (!strcasecmp(child->Value(),"dtol")) {
if (child->FirstChild() && child->FirstChild()->Value())
divgLim = atof(child->FirstChild()->Value());
} else if (!strcasecmp(child->Value(),"eta")) {
if (child->FirstChild() && child->FirstChild()->Value())
eta = atof(child->FirstChild()->Value());
} }
child = child->NextSiblingElement(); else if ((value = utl::getValue(child,"maxits")))
} maxit = atoi(value);
else if ((value = utl::getValue(child,"nupdate")))
nupdat = atoi(value);
else if ((value = utl::getValue(child,"rtol")))
convTol = atof(value);
else if ((value = utl::getValue(child,"dtol")))
divgLim = atof(value);
else if ((value = utl::getValue(child,"eta")))
eta = atof(value);
return true; return true;
} }
void NonLinSIM::initSystem (int mType, size_t nGauss)
{
model->initSystem(mType,1,1);
model->setAssociatedRHS(0,0);
model->setQuadratureRule(nGauss);
}
void NonLinSIM::init (SolvePrm& param, const RealArray& initVal) void NonLinSIM::init (SolvePrm& param, const RealArray& initVal)
{ {
param.initTime(startTime,stopTime,steps); param.initTime(startTime,stopTime,steps);
@ -226,6 +226,7 @@ bool NonLinSIM::solveStep (SolvePrm& param, SIM::SolutionMode mode,
return false; return false;
bool newTangent = true; bool newTangent = true;
model->setQuadratureRule(model->opt.nGauss[0]);
if (!model->assembleSystem(param.time,solution,newTangent)) if (!model->assembleSystem(param.time,solution,newTangent))
return false; return false;
@ -424,8 +425,11 @@ bool NonLinSIM::solutionNorms (const TimeDomain& time, const char* compName,
bool haveRFs = model->getCurrentReactions(RF,solution.front()); bool haveRFs = model->getCurrentReactions(RF,solution.front());
if (energyNorm) if (energyNorm)
{
model->setQuadratureRule(model->opt.nGauss[1]);
if (!model->solutionNorms(time,solution,gNorm)) if (!model->solutionNorms(time,solution,gNorm))
gNorm.clear(); gNorm.clear();
}
if (myPid == 0) if (myPid == 0)
{ {