From eb4f881a61ce195ccc524c2957c63b8b6b01f5fe Mon Sep 17 00:00:00 2001 From: akva Date: Thu, 8 Dec 2011 16:23:47 +0000 Subject: [PATCH] added: support for XML based input files in the Linear Elasticity applications git-svn-id: http://svn.sintef.no/trondheim/IFEM/trunk@1343 e10b68d5-8a6e-419e-a041-bce267b0401d --- Apps/LinearElasticity/SIMLinEl2D.C | 163 ++++++++++++++++++ Apps/LinearElasticity/SIMLinEl2D.h | 4 + Apps/LinearElasticity/SIMLinEl3D.C | 164 +++++++++++++++++++ Apps/LinearElasticity/SIMLinEl3D.h | 4 + Apps/LinearElasticity/SIMLinElBeamC1.C | 92 +++++++++++ Apps/LinearElasticity/SIMLinElBeamC1.h | 4 + Apps/LinearElasticity/SIMLinElKL.C | 136 +++++++++++++++ Apps/LinearElasticity/SIMLinElKL.h | 4 + Apps/LinearElasticity/Test/Hole2D-NURBS.xinp | 59 +++++++ 9 files changed, 630 insertions(+) create mode 100644 Apps/LinearElasticity/Test/Hole2D-NURBS.xinp diff --git a/Apps/LinearElasticity/SIMLinEl2D.C b/Apps/LinearElasticity/SIMLinEl2D.C index 2beeef79..d3c983b8 100644 --- a/Apps/LinearElasticity/SIMLinEl2D.C +++ b/Apps/LinearElasticity/SIMLinEl2D.C @@ -21,6 +21,8 @@ #include "AnaSol.h" #include +#include "tinyxml.h" + bool SIMLinEl2D::planeStrain = false; bool SIMLinEl2D::axiSymmetry = false; @@ -294,6 +296,167 @@ bool SIMLinEl2D::parse (char* keyWord, std::istream& is) } +bool SIMLinEl2D::parse(const TiXmlElement* elem) +{ + if (strcasecmp(elem->Value(),"linearelasticity")) + return SIM2D::parse(elem); + + Elasticity* elp = dynamic_cast(myProblem); + if (!elp) + return false; + + std::vector parsed = handlePriorityTags(elem); + const TiXmlElement* child = elem->FirstChildElement(); + while (child) { + if (find(parsed.begin(),parsed.end(),child) != parsed.end()) { + child = child->NextSiblingElement(); + continue; + } + if (!strcasecmp(child->Value(),"gravity")) { + double gx=0, gy=0; + if (child->Attribute("x")) + gx = atof(child->Attribute("x")); + if (child->Attribute("y")) + gy = atof(child->Attribute("y")); + elp->setGravity(gx,gy); + if (myPid == 0) + std::cout <<"\nGravitation vector: " << gx <<" "<< gy << std::endl; + } else if (!strcasecmp(child->Value(),"isotropic")) { + int code=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (code > 0) + setPropertyType(code,Property::MATERIAL,mVec.size()); + double E=1000.f, nu=0.3f, rho = 1.f; + if (child->Attribute("E")) + E = atof(child->Attribute("E")); + if (child->Attribute("rho")) + rho = atof(child->Attribute("rho")); + if (child->Attribute("nu")) + nu = atof(child->Attribute("nu")); + mVec.push_back(new LinIsotropic(E,nu,rho,!planeStrain,axiSymmetry)); + if (myPid == 0) + std::cout <<"\tMaterial code "<< code <<": " + << E <<" "<< nu <<" "<< rho << std::endl; + if (!mVec.empty()) + elp->setMaterial(mVec.front()); + } else if (!strcasecmp(child->Value(),"constantpressure")) { + int code=0, pdir=1; + double p=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (child->Attribute("dir")) + pdir = atoi(child->Attribute("dir")); + if (child->FirstChild() && child->FirstChild()->Value()) + p = atof(child->FirstChild()->Value()); + if (myPid == 0) + std::cout <<"\tPressure code "<< code <<" direction "<< pdir + <<": "<< p << std::endl; + setPropertyType(code,Property::NEUMANN); + myTracs[code] = new PressureField(p,pdir); + } else if (!strcasecmp(child->Value(),"anasol")) { + if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"hole")) { + double a=0, F0 = 0, nu = 0; + if (child->Attribute("a")) + a = atof(child->Attribute("a")); + if (child->Attribute("F0")) + F0 = atof(child->Attribute("F0")); + if (child->Attribute("nu")) + nu = atof(child->Attribute("nu")); + mySol = new AnaSol(new Hole(a,F0,nu)); + std::cout <<"\nAnalytical solution: Hole a="<< a <<" F0="<< F0 + <<" nu="<< nu << std::endl; + } else if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"Lshape")) { + double a=0, F0 = 0, nu = 0; + if (child->Attribute("a")) + a = atof(child->Attribute("a")); + if (child->Attribute("F0")) + F0 = atof(child->Attribute("F0")); + if (child->Attribute("nu")) + nu = atof(child->Attribute("nu")); + mySol = new AnaSol(new Lshape(a,F0,nu)); + std::cout <<"\nAnalytical solution: Lshape a="<< a <<" F0="<< F0 + <<" nu="<< nu << std::endl; + } else if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"cants")) { + double L=0, F0 = 0, H = 0; + if (child->Attribute("L")) + L = atof(child->Attribute("L")); + if (child->Attribute("F0")) + F0 = atof(child->Attribute("F0")); + if (child->Attribute("H")) + H = atof(child->Attribute("H")); + mySol = new AnaSol(new CanTS(L,H,F0)); + std::cout <<"\nAnalytical solution: CanTS L="<< L <<" H="<< H + <<" F0="<< F0 << std::endl; + } else if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"cantm")) { + double M0 = 0, H = 0; + if (child->Attribute("M0")) + M0 = atof(child->Attribute("M0")); + if (child->Attribute("H")) + H = atof(child->Attribute("H")); + mySol = new AnaSol(new CanTM(H,M0)); + std::cout <<"\nAnalytical solution: CanTM H="<< H + <<" M0="<< M0 << std::endl; + } else if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"curved")) { + double a=0, b=0, u0 = 0, E = 0; + if (child->Attribute("a")) + a = atof(child->Attribute("a")); + if (child->Attribute("b")) + b = atof(child->Attribute("b")); + if (child->Attribute("u0")) + u0 = atof(child->Attribute("u0")); + if (child->Attribute("E")) + E = atof(child->Attribute("E")); + mySol = new AnaSol(new CurvedBeam(u0,a,b,E)); + std::cout <<"\nAnalytical solution: Curved Beam a="<< a <<" b="<< b + <<" u0="<< u0 <<" E="<< E << std::endl; + } else if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"expression")) { + std::string variables, stress; + const TiXmlElement* var = child->FirstChildElement("variables"); + if (var && var->FirstChild() && var->FirstChild()->Value()) { + variables = var->FirstChild()->Value(); + if (variables[variables.size()-1] != ';') + variables += ";"; + } + const TiXmlElement* str = child->FirstChildElement("stress"); + if (str && str->FirstChild() && str->FirstChild()->Value()) + stress = str->FirstChild()->Value(); + mySol = new AnaSol(new STensorFuncExpr(stress,variables)); + std::cout <<"\nAnalytical solution:"; + if (!variables.empty()) + std::cout <<"\n\tVariables = "<< variables; + std::cout <<"\n\tStress = "<< stress << std::endl; + } else { + std::cerr <<" ** SIMLinEl2D::parse: Unknown analytical solution " + << (child->Attribute("type")?child->Attribute("type"):"") << std::endl; + } + + // Define the analytical boundary traction field + int code=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (code > 0 && mySol && mySol->getStressSol()) { + std::cout <<"Pressure code "<< code <<": Analytical traction"<< std::endl; + setPropertyType(code,Property::NEUMANN); + myTracs[code] = new TractionField(*mySol->getStressSol()); + } + } + else + SIM2D::parse(child); + + child = child->NextSiblingElement(); + } + + return true; +} + + bool SIMLinEl2D::initMaterial (size_t propInd) { Elasticity* elp = dynamic_cast(myProblem); diff --git a/Apps/LinearElasticity/SIMLinEl2D.h b/Apps/LinearElasticity/SIMLinEl2D.h index 19dad203..4a98b150 100644 --- a/Apps/LinearElasticity/SIMLinEl2D.h +++ b/Apps/LinearElasticity/SIMLinEl2D.h @@ -42,6 +42,10 @@ protected: //! \param is The file stream to read from virtual bool parse(char* keyWord, std::istream& is); + //! \brief Parses a data section from an XML element + //! \param[in] elem The XML element to parse + virtual bool parse(const TiXmlElement* elem); + //! \brief Initializes material properties for integration of interior terms. //! \param[in] propInd Physical property index virtual bool initMaterial(size_t propInd); diff --git a/Apps/LinearElasticity/SIMLinEl3D.C b/Apps/LinearElasticity/SIMLinEl3D.C index 641e9d16..1a360276 100644 --- a/Apps/LinearElasticity/SIMLinEl3D.C +++ b/Apps/LinearElasticity/SIMLinEl3D.C @@ -23,6 +23,8 @@ #include "AnaSol.h" #include +#include "tinyxml.h" + /*! \brief Local coordinate system for a cylinder along global z-axis. @@ -427,6 +429,168 @@ bool SIMLinEl3D::parse (char* keyWord, std::istream& is) } +bool SIMLinEl3D::parse(const TiXmlElement* elem) +{ + if (strcasecmp(elem->Value(),"linearelasticity")) + return SIM3D::parse(elem); + + Elasticity* elp = dynamic_cast(myProblem); + if (!elp) + return false; + + std::vector parsed = handlePriorityTags(elem); + const TiXmlElement* child = elem->FirstChildElement(); + while (child) { + if (find(parsed.begin(),parsed.end(),child) != parsed.end()) { + child = child->NextSiblingElement(); + continue; + } + if (!strcasecmp(child->Value(),"gravity")) { + double gx=0, gy=0, gz=0; + if (child->Attribute("x")) + gx = atof(child->Attribute("x")); + if (child->Attribute("y")) + gy = atof(child->Attribute("y")); + if (child->Attribute("z")) + gz = atof(child->Attribute("z")); + elp->setGravity(gx,gy,gz); + if (myPid == 0) + std::cout <<"\nGravitation vector: " << gx <<" "<< gy <<" " << gz << std::endl; + } else if (!strcasecmp(child->Value(),"isotropic")) { + int code=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (code > 0) + setPropertyType(code,Property::MATERIAL,mVec.size()); + double E=1000.f, nu=0.3f, rho = 1.f; + if (child->Attribute("E")) + E = atof(child->Attribute("E")); + if (child->Attribute("rho")) + rho = atof(child->Attribute("rho")); + if (child->Attribute("nu")) + nu = atof(child->Attribute("nu")); + mVec.push_back(new LinIsotropic(E,nu,rho)); + if (myPid == 0) + std::cout <<"\tMaterial code "<< code <<": " + << E <<" "<< nu <<" "<< rho << std::endl; + if (!mVec.empty()) + elp->setMaterial(mVec.front()); + } else if (!strcasecmp(child->Value(),"constantpressure") || + !strcasecmp(child->Value(),"linearpressure")) { + int code=0, pdir=1; + double p=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (child->Attribute("dir")) + pdir = atoi(child->Attribute("dir")); + if (child->FirstChild() && child->FirstChild()->Value()) + p = atof(child->FirstChild()->Value()); + if (myPid == 0) + std::cout <<"\tPressure code "<< code <<" direction "<< pdir + <<": "<< p << std::endl; + setPropertyType(code,Property::NEUMANN); + if (!strcasecmp(child->Value(),"linearpressure")) { + RealFunc* pfl = new ConstTimeFunc(new LinearFunc(p)); + myTracs[code] = new PressureField(pfl,pdir); + } else + myTracs[code] = new PressureField(p,pdir); + } else if (!strcasecmp(child->Value(),"anasol")) { + if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"hole")) { + double a=0, F0 = 0, nu = 0; + if (child->Attribute("a")) + a = atof(child->Attribute("a")); + if (child->Attribute("F0")) + F0 = atof(child->Attribute("F0")); + if (child->Attribute("nu")) + nu = atof(child->Attribute("nu")); + mySol = new AnaSol(new Hole(a,F0,nu,true)); + std::cout <<"\nAnalytical solution: Hole a="<< a <<" F0="<< F0 + <<" nu="<< nu << std::endl; + } else if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"Lshape")) { + double a=0, F0 = 0, nu = 0; + if (child->Attribute("a")) + a = atof(child->Attribute("a")); + if (child->Attribute("F0")) + F0 = atof(child->Attribute("F0")); + if (child->Attribute("nu")) + nu = atof(child->Attribute("nu")); + mySol = new AnaSol(new Lshape(a,F0,nu,true)); + std::cout <<"\nAnalytical solution: Lshape a="<< a <<" F0="<< F0 + <<" nu="<< nu << std::endl; + } else if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"cants")) { + double L=0, F0 = 0, H = 0; + if (child->Attribute("L")) + L = atof(child->Attribute("L")); + if (child->Attribute("F0")) + F0 = atof(child->Attribute("F0")); + if (child->Attribute("H")) + H = atof(child->Attribute("H")); + mySol = new AnaSol(new CanTS(L,H,F0,true)); + std::cout <<"\nAnalytical solution: CanTS L="<< L <<" H="<< H + <<" F0="<< F0 << std::endl; + } else if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"expression")) { + std::string variables, stress; + const TiXmlElement* var = child->FirstChildElement("variables"); + if (var && var->FirstChild() && var->FirstChild()->Value()) { + variables = var->FirstChild()->Value(); + if (variables[variables.size()-1] != ';') + variables += ";"; + } + const TiXmlElement* str = child->FirstChildElement("stress"); + if (str && str->FirstChild() && str->FirstChild()->Value()) + stress = str->FirstChild()->Value(); + mySol = new AnaSol(new STensorFuncExpr(stress,variables)); + std::cout <<"\nAnalytical solution:"; + if (!variables.empty()) + std::cout <<"\n\tVariables = "<< variables; + std::cout <<"\n\tStress = "<< stress << std::endl; + } else { + std::cerr <<" ** SIMLinEl3D::parse: Unknown analytical solution " + << (child->Attribute("type")?child->Attribute("type"):"") << std::endl; + } + + // Define the analytical boundary traction field + int code=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (code > 0 && mySol &&mySol->getStressSol()) { + std::cout <<"Pressure code "<< code <<": Analytical traction"<< std::endl; + setPropertyType(code,Property::NEUMANN); + myTracs[code] = new TractionField(*mySol->getStressSol()); + } + } else if (!strcasecmp(child->Value(),"localsystem")) { + bool ok=true; + if (child->FirstChild() && child->FirstChild()->Value()) { + if (!strcasecmp(child->FirstChild()->Value(),"cylindricz")) + elp->setLocalSystem(new CylinderCS); + else if (!strcasecmp(child->FirstChild()->Value(),"cylinder+sphere")) { + double H = 0; + if (child->Attribute("H")) + H = atof(child->Attribute("H")); + elp->setLocalSystem(new CylinderSphereCS(H)); + } + else + ok = false; + } + if (!ok) { + const char* type=NULL; + if (child->FirstChild() && child->FirstChild()->Value()) + type = child->FirstChild()->Value(); + std::cerr <<" *** SIMLinEl3D::parse: Unsupported coordinate system: " + << (type?type:"") << std::endl; + } + } + child = child->NextSiblingElement(); + } + + return true; +} + + bool SIMLinEl3D::initMaterial (size_t propInd) { Elasticity* elp = dynamic_cast(myProblem); diff --git a/Apps/LinearElasticity/SIMLinEl3D.h b/Apps/LinearElasticity/SIMLinEl3D.h index 438ddaca..e4124509 100644 --- a/Apps/LinearElasticity/SIMLinEl3D.h +++ b/Apps/LinearElasticity/SIMLinEl3D.h @@ -43,6 +43,10 @@ protected: //! \param is The file stream to read from virtual bool parse(char* keyWord, std::istream& is); + //! \brief Parses a data section from an XML element + //! \param[in] elem The XML element to parse + virtual bool parse(const TiXmlElement* elem); + //! \brief Initializes material properties for integration of interior terms. //! \param[in] propInd Physical property index virtual bool initMaterial(size_t propInd); diff --git a/Apps/LinearElasticity/SIMLinElBeamC1.C b/Apps/LinearElasticity/SIMLinElBeamC1.C index 9fa19d34..e6cf87fa 100644 --- a/Apps/LinearElasticity/SIMLinElBeamC1.C +++ b/Apps/LinearElasticity/SIMLinElBeamC1.C @@ -25,6 +25,8 @@ #include "Vec3Oper.h" #include +#include "tinyxml.h" + SIMLinElBeamC1::SIMLinElBeamC1 () : SIM1D(SIM::NONE) { @@ -155,6 +157,96 @@ bool SIMLinElBeamC1::parse (char* keyWord, std::istream& is) } +bool SIMLinElBeamC1::parse(const TiXmlElement* elem) +{ + if (strcasecmp(elem->Value(),"eulerbernoulli")) + return SIM1D::parse(elem); + + KirchhoffLovePlate* klp = dynamic_cast(myProblem); + if (!klp) + return false; + + std::vector parsed = handlePriorityTags(elem); + const TiXmlElement* child = elem->FirstChildElement(); + while (child) { + if (find(parsed.begin(),parsed.end(),child) != parsed.end()) { + child = child->NextSiblingElement(); + continue; + } + if (!strcasecmp(child->Value(),"gravity")) { + double g=0; + if (child->Attribute("g")) + g = atof(child->Attribute("g")); + klp->setGravity(g); + if (myPid == 0) + std::cout <<"\nGravitation constant: "<< g << std::endl; + } else if (!strcasecmp(child->Value(),"isotropic")) { + int code=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (code > 0) + setPropertyType(code,Property::MATERIAL,mVec.size()); + double E=1000.f, rho = 1.f, thk=0; + if (child->Attribute("E")) + E = atof(child->Attribute("E")); + if (child->Attribute("rho")) + rho = atof(child->Attribute("rho")); + if (child->Attribute("thickness")) + thk = atof(child->Attribute("thickness")); + mVec.push_back(new LinIsotropic(E,0.0,rho,true)); + tVec.push_back(thk); + if (myPid == 0) + std::cout <<"\tMaterial code "<< code <<": " + << E <<" " << rho <<" " << thk << std::endl; + if (!mVec.empty()) + klp->setMaterial(mVec.front()); + if (!tVec.empty() && tVec.front() != 0.0) + klp->setThickness(tVec.front()); + } else if (!strcasecmp(child->Value(),"pointload")) { + PointLoad load; + if (child->Attribute("patch")) + load.patch = atoi(child->Attribute("patch")); + if (child->Attribute("xi")) + load.xi = atof(child->Attribute("xi")); + if (child->FirstChild() && child->FirstChild()->Value()) + load.pload = atof(child->FirstChild()->Value()); + if (myPid == 0) + std::cout <<"\n\tPoint: P" << load.patch + <<" xi = " << load.xi + <<" load = " << load.pload; + myLoads.push_back(load); + } else if (!strcasecmp(child->Value(),"pressure")) { + int code=0; + double p=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (child->FirstChild() && child->FirstChild()->Value()) + p = atof(child->FirstChild()->Value()); + if (myPid == 0) + std::cout <<"\tPressure code "<< code <<": "; + char* function; + if (child->Attribute("function")) + function = strdup(child->Attribute("function")); + else + function = strdup(""); + char* function2 = function; + function = strtok(function," "); + myScalars[code] = const_cast(utl::parseRealFunc(const_cast(function),p)); + free(function2); + std::cout << std::endl; + if (code > 0) + setPropertyType(code,Property::BODYLOAD); + } + else + SIM1D::parse(child); + + child = child->NextSiblingElement(); + } + + return true; +} + + bool SIMLinElBeamC1::initMaterial (size_t propInd) { if (propInd >= mVec.size()) propInd = mVec.size()-1; diff --git a/Apps/LinearElasticity/SIMLinElBeamC1.h b/Apps/LinearElasticity/SIMLinElBeamC1.h index 978137b9..98ebb9c8 100644 --- a/Apps/LinearElasticity/SIMLinElBeamC1.h +++ b/Apps/LinearElasticity/SIMLinElBeamC1.h @@ -60,6 +60,10 @@ protected: //! \param is The file stream to read from virtual bool parse(char* keyWord, std::istream& is); + //! \brief Parses a data section from an XML element + //! \param[in] elem The XML element to parse + virtual bool parse(const TiXmlElement* elem); + //! \brief Initializes material properties for integration of interior terms. //! \param[in] propInd Physical property index virtual bool initMaterial(size_t propInd); diff --git a/Apps/LinearElasticity/SIMLinElKL.C b/Apps/LinearElasticity/SIMLinElKL.C index 4e052dbe..6f6aae60 100644 --- a/Apps/LinearElasticity/SIMLinElKL.C +++ b/Apps/LinearElasticity/SIMLinElKL.C @@ -25,6 +25,8 @@ #include "Vec3Oper.h" #include +#include "tinyxml.h" + SIMLinElKL::SIMLinElKL () : SIMLinEl2D(SIM::NONE) { @@ -162,6 +164,140 @@ bool SIMLinElKL::parse (char* keyWord, std::istream& is) } +bool SIMLinElKL::parse(const TiXmlElement* elem) +{ + if (strcasecmp(elem->Value(),"kirchhofflove")) + return SIMLinEl2D::parse(elem); + + KirchhoffLovePlate* klp = dynamic_cast(myProblem); + if (!klp) + return false; + + std::vector parsed = handlePriorityTags(elem); + const TiXmlElement* child = elem->FirstChildElement(); + while (child) { + if (find(parsed.begin(),parsed.end(),child) != parsed.end()) { + child = child->NextSiblingElement(); + continue; + } + if (!strcasecmp(child->Value(),"gravity")) { + double g=0; + if (child->Attribute("g")) + g = atof(child->Attribute("g")); + klp->setGravity(g); + if (myPid == 0) + std::cout <<"\nGravitation constant: "<< g << std::endl; + } else if (!strcasecmp(child->Value(),"isotropic")) { + int code=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (code > 0) + setPropertyType(code,Property::MATERIAL,mVec.size()); + double E=1000.f, nu=0.3f, rho = 1.f, thk=0; + if (child->Attribute("E")) + E = atof(child->Attribute("E")); + if (child->Attribute("rho")) + rho = atof(child->Attribute("rho")); + if (child->Attribute("nu")) + nu = atof(child->Attribute("nu")); + if (child->Attribute("thickness")) + thk = atof(child->Attribute("thickness")); + mVec.push_back(new LinIsotropic(E,nu,rho,true)); + tVec.push_back(thk); + if (myPid == 0) + std::cout <<"\tMaterial code "<< code <<": " + << E <<" "<< nu <<" "<< rho <<" " << thk << std::endl; + if (!mVec.empty()) + klp->setMaterial(mVec.front()); + if (!tVec.empty() && tVec.front() != 0.0) + klp->setThickness(tVec.front()); + } else if (!strcasecmp(child->Value(),"pointload")) { + PointLoad load; + if (child->Attribute("patch")) + load.patch = atoi(child->Attribute("patch")); + if (child->Attribute("xi")) + load.xi[0] = atof(child->Attribute("xi")); + if (child->Attribute("eta")) + load.xi[1] = atof(child->Attribute("eta")); + if (child->FirstChild() && child->FirstChild()->Value()) + load.pload = atof(child->FirstChild()->Value()); + if (myPid == 0) + std::cout <<"\n\tPoint: P"<< load.patch + <<" xi = "<< load.xi[0] <<" "<< load.xi[1] + <<" load = "<< load.pload; + myLoads.push_back(load); + } else if (!strcasecmp(child->Value(),"pressure")) { + int code=0; + double p=0; + if (child->Attribute("code")) + code = atoi(child->Attribute("code")); + if (child->FirstChild() && child->FirstChild()->Value()) + p = atof(child->FirstChild()->Value()); + if (myPid == 0) + std::cout <<"\tPressure code "<< code <<": "; + char* function; + if (child->Attribute("function")) + function = strdup(child->Attribute("function")); + else + function = strdup(""); + char* function2 = function; + function = strtok(function," "); + myScalars[code] = const_cast(utl::parseRealFunc(const_cast(function),p)); + free(function2); + std::cout << std::endl; + if (code > 0) + setPropertyType(code,Property::BODYLOAD); + } else if (!strcasecmp(child->Value(),"anasol")) { + if (child->Attribute("type") && + !strcasecmp(child->Attribute("type"),"navierplate")) { + double a=0, b=0, c=0, d=0, t=0, E=0, nu=0, pz=0, xi=0, eta=0; + if (child->Attribute("a")) + a = atof(child->Attribute("a")); + if (child->Attribute("b")) + b = atof(child->Attribute("b")); + if (child->Attribute("c")) + c = atof(child->Attribute("c")); + if (child->Attribute("d")) + d = atof(child->Attribute("d")); + if (child->Attribute("t")) + t = atof(child->Attribute("t")); + if (child->Attribute("E")) + E = atof(child->Attribute("E")); + if (child->Attribute("nu")) + nu = atof(child->Attribute("nu")); + if (child->Attribute("pz")) + pz = atof(child->Attribute("pz")); + if (child->Attribute("xi")) + xi = atof(child->Attribute("xi")); + if (child->Attribute("eta")) + eta = atof(child->Attribute("eta")); + std::cout <<"\nAnalytic solution: NavierPlate a="<< a <<" b="<< b + <<" t="<< t <<" E="<< E <<" nu="<< nu <<" pz="<< pz; + + if (xi != 0.f && eta != 0.f) { + std::cout <<" xi="<< xi <<" eta="<< eta; + if (c != 0.f && d != 0.f) { + std::cout <<" c="<< c <<" d="<< d; + mySol = new AnaSol(new NavierPlate(a,b,t,E,nu,pz,xi,eta,c,d)); + } else + mySol = new AnaSol(new NavierPlate(a,b,t,E,nu,pz,xi,eta)); + } else + mySol = new AnaSol(new NavierPlate(a,b,t,E,nu,pz)); + } else { + std::cerr <<" ** SIMLinElKL::parse: Unknown analytical solution " + << (child->Attribute("type")?child->Attribute("type"):"") << std::endl; + } + } + else + SIMLinEl2D::parse(child); + + child = child->NextSiblingElement(); + } + + return true; +} + + bool SIMLinElKL::initMaterial (size_t propInd) { if (propInd >= mVec.size()) propInd = mVec.size()-1; diff --git a/Apps/LinearElasticity/SIMLinElKL.h b/Apps/LinearElasticity/SIMLinElKL.h index 292cbc1c..89e2ee1b 100644 --- a/Apps/LinearElasticity/SIMLinElKL.h +++ b/Apps/LinearElasticity/SIMLinElKL.h @@ -59,6 +59,10 @@ protected: //! \param is The file stream to read from virtual bool parse(char* keyWord, std::istream& is); + //! \brief Parses a data section from an XML element + //! \param[in] elem The XML element to parse + virtual bool parse(const TiXmlElement* elem); + //! \brief Initializes material properties for integration of interior terms. //! \param[in] propInd Physical property index virtual bool initMaterial(size_t propInd); diff --git a/Apps/LinearElasticity/Test/Hole2D-NURBS.xinp b/Apps/LinearElasticity/Test/Hole2D-NURBS.xinp new file mode 100644 index 00000000..040e55cd --- /dev/null +++ b/Apps/LinearElasticity/Test/Hole2D-NURBS.xinp @@ -0,0 +1,59 @@ + + + + + hole2D.g2 + + + + + + + + + + + + + + + + + + + + + cg + ilu + petsc + 2 + 1.0e-12 + 1.0e-12 + 1.0e6 + 5000 + + + + + + nu=0.3; + F0=10; + a=1; + R=sqrt(x*x+y*y); + R2=if(above(R,a),a*a/(R*R),1); + R4=if(above(R,a),R2*R2,1); + th=atan2(y,x); + C2=cos(2*th); + C4=cos(4*th); + S2=sin(2*th); + S4=sin(4*th) + + F0*(1-R2*(1.5*C2+C4)+1.5*R4*C4)| + F0*(-R2*(0.5*S2+S4)+1.5*R4*S4)| + F0*(-R2*(0.5*C2-C4)-1.5*R4*C4)| + F0*nu*(1-2*R2*C2) + + + + +