• Main Page
  • Related Pages
  • Namespaces
  • Classes
  • Files
  • File List
  • File Members

ElasticityProblem.cpp

Go to the documentation of this file.
00001 //*****************************************************************//
00002 //    Albany 2.0:  Copyright 2012 Sandia Corporation               //
00003 //    This Software is released under the BSD license detailed     //
00004 //    in the file "license.txt" in the top-level Albany directory  //
00005 //*****************************************************************//
00006 #include "ElasticityProblem.hpp"
00007 #include "Albany_Utils.hpp"
00008 #include "Albany_ProblemUtils.hpp"
00009 
00010 Albany::ElasticityProblem::
00011 ElasticityProblem(const Teuchos::RCP<Teuchos::ParameterList>& params_,
00012       const Teuchos::RCP<ParamLib>& paramLib_,
00013       const int numDim_) :
00014   Albany::AbstractProblem(params_, paramLib_, numDim_),
00015   haveSource(false),
00016   numDim(numDim_)
00017 {
00018   std::string& method = params->get("Name", "Elasticity ");
00019   *out << "Problem Name = " << method << std::endl;
00020 
00021   haveSource =  params->isSublist("Source Functions");
00022 
00023   matModel = params->sublist("Material Model").get("Model Name", "LinearElasticity");
00024 
00025   computeError = params->get<bool>("Compute Error", false);
00026 
00027   if (computeError)
00028     this->setNumEquations(2*neq);
00029 
00030 // the following function returns the problem information required for setting the rigid body modes (RBMs) for elasticity problems
00031 //written by IK, Feb. 2012
00032 
00033   int numScalar = 0;
00034   int nullSpaceDim = 0;
00035   if (numDim == 1) {nullSpaceDim = 0; }
00036   else {
00037     if (numDim == 2) {nullSpaceDim = 3; }
00038     if (numDim == 3) {nullSpaceDim = 6; }
00039   }
00040 
00041   if (computeError)
00042     rigidBodyModes->setParameters(2*numDim, numDim, numScalar, nullSpaceDim);
00043   else
00044     rigidBodyModes->setParameters(numDim, numDim, numScalar, nullSpaceDim);
00045 
00046 }
00047 
00048 Albany::ElasticityProblem::
00049 ~ElasticityProblem()
00050 {
00051 }
00052 
00053 void
00054 Albany::ElasticityProblem::
00055 buildProblem(
00056   Teuchos::ArrayRCP<Teuchos::RCP<Albany::MeshSpecsStruct> >  meshSpecs,
00057   Albany::StateManager& stateMgr)
00058 {
00059   /* Construct All Phalanx Evaluators */
00060   TEUCHOS_TEST_FOR_EXCEPTION(meshSpecs.size()!=1,std::logic_error,"Problem supports one Material Block");
00061 
00062   fm.resize(1);
00063 
00064   fm[0]  = Teuchos::rcp(new PHX::FieldManager<PHAL::AlbanyTraits>);
00065   buildEvaluators(*fm[0], *meshSpecs[0], stateMgr, BUILD_RESID_FM, 
00066       Teuchos::null);
00067 
00068   if(meshSpecs[0]->nsNames.size() > 0) // Build a nodeset evaluator if nodesets are present
00069 
00070     constructDirichletEvaluators(*meshSpecs[0]);
00071 
00072   if(meshSpecs[0]->ssNames.size() > 0) // Build a sideset evaluator if sidesets are present
00073 
00074     constructNeumannEvaluators(meshSpecs[0]);
00075 
00076 }
00077 
00078 Teuchos::Array<Teuchos::RCP<const PHX::FieldTag> >
00079 Albany::ElasticityProblem::
00080 buildEvaluators(
00081   PHX::FieldManager<PHAL::AlbanyTraits>& fm0,
00082   const Albany::MeshSpecsStruct& meshSpecs,
00083   Albany::StateManager& stateMgr,
00084   Albany::FieldManagerChoice fmchoice,
00085   const Teuchos::RCP<Teuchos::ParameterList>& responseList)
00086 {
00087   // Call constructeEvaluators<EvalT>(*rfm[0], *meshSpecs[0], stateMgr);
00088   // for each EvalT in PHAL::AlbanyTraits::BEvalTypes
00089   ConstructEvaluatorsOp<ElasticityProblem> op(
00090     *this, fm0, meshSpecs, stateMgr, fmchoice, responseList);
00091   boost::mpl::for_each<PHAL::AlbanyTraits::BEvalTypes>(op);
00092   return *op.tags;
00093 }
00094 
00095 // Dirichlet BCs
00096 void
00097 Albany::ElasticityProblem::constructDirichletEvaluators(
00098         const Albany::MeshSpecsStruct& meshSpecs)
00099 {
00100   // Construct Dirichlet evaluators for all nodesets and names
00101   std::vector<std::string> dirichletNames(neq);
00102   dirichletNames[0] = "X";
00103   if (neq>1) dirichletNames[1] = "Y";
00104   if (neq>2) dirichletNames[2] = "Z";
00105   Albany::BCUtils<Albany::DirichletTraits> dirUtils;
00106   dfm = dirUtils.constructBCEvaluators(meshSpecs.nsNames, dirichletNames,
00107                                        this->params, this->paramLib);
00108 }
00109 
00110 // Neumann BCs
00111 void
00112 Albany::ElasticityProblem::constructNeumannEvaluators(
00113         const Teuchos::RCP<Albany::MeshSpecsStruct>& meshSpecs)
00114 {
00115    // Note: we only enter this function if sidesets are defined in the mesh file
00116    // i.e. meshSpecs.ssNames.size() > 0
00117 
00118    Albany::BCUtils<Albany::NeumannTraits> neuUtils;
00119 
00120    // Check to make sure that Neumann BCs are given in the input file
00121 
00122    if(!neuUtils.haveBCSpecified(this->params))
00123 
00124       return;
00125 
00126    // Construct BC evaluators for all side sets and names
00127    // Note that the string index sets up the equation offset, so ordering is important
00128    std::vector<std::string> neumannNames(neq + 1);
00129    Teuchos::Array<Teuchos::Array<int> > offsets;
00130    offsets.resize(neq + 1);
00131 
00132    neumannNames[0] = "sig_x";
00133    offsets[0].resize(1);
00134    offsets[0][0] = 0;
00135    offsets[neq].resize(neq);
00136    offsets[neq][0] = 0;
00137 
00138    if (neq>1){ 
00139       neumannNames[1] = "sig_y";
00140       offsets[1].resize(1);
00141       offsets[1][0] = 1;
00142       offsets[neq][1] = 1;
00143    }
00144 
00145    if (neq>2){
00146      neumannNames[2] = "sig_z";
00147       offsets[2].resize(1);
00148       offsets[2][0] = 2;
00149       offsets[neq][2] = 2;
00150    }
00151 
00152    neumannNames[neq] = "all";
00153 
00154    // Construct BC evaluators for all possible names of conditions
00155    // Should only specify flux vector components (dudx, dudy, dudz), or dudn, not both
00156    std::vector<std::string> condNames(3); //dudx, dudy, dudz, dudn, P
00157    Teuchos::ArrayRCP<std::string> dof_names(1);
00158      dof_names[0] = "Displacement";
00159 
00160    // Note that sidesets are only supported for two and 3D currently
00161    if(numDim == 2)
00162     condNames[0] = "(t_x, t_y)";
00163    else if(numDim == 3)
00164     condNames[0] = "(t_x, t_y, t_z)";
00165    else
00166     TEUCHOS_TEST_FOR_EXCEPTION(true, Teuchos::Exceptions::InvalidParameter,
00167        std::endl << "Error: Sidesets only supported in 2 and 3D." << std::endl);
00168 
00169    condNames[1] = "dudn";
00170    condNames[2] = "P";
00171 
00172    nfm.resize(1); // Elasticity problem only has one element block
00173 
00174    nfm[0] = neuUtils.constructBCEvaluators(meshSpecs, neumannNames, dof_names, true, 0,
00175                                           condNames, offsets, dl,
00176                                           this->params, this->paramLib);
00177 
00178 }
00179 
00180 
00181 Teuchos::RCP<const Teuchos::ParameterList>
00182 Albany::ElasticityProblem::getValidProblemParameters() const
00183 {
00184   Teuchos::RCP<Teuchos::ParameterList> validPL =
00185     this->getGenericProblemParams("ValidElasticityProblemParams");
00186 
00187   validPL->sublist("Elastic Modulus", false, "");
00188   validPL->sublist("Poissons Ratio", false, "");
00189   validPL->sublist("Material Model", false, "");
00190 
00191   validPL->set<bool>("Compute Error", false, "");
00192 
00193   if (matModel == "CapExplicit"|| matModel == "CapImplicit")
00194   {
00195   validPL->set<double>("A", false, "");
00196   validPL->set<double>("B", false, "");
00197   validPL->set<double>("C", false, "");
00198   validPL->set<double>("theta", false, "");
00199   validPL->set<double>("R", false, "");
00200   validPL->set<double>("kappa0", false, "");
00201   validPL->set<double>("W", false, "");
00202   validPL->set<double>("D1", false, "");
00203   validPL->set<double>("D2", false, "");
00204   validPL->set<double>("calpha", false, "");
00205   validPL->set<double>("psi", false, "");
00206   validPL->set<double>("N", false, "");
00207   validPL->set<double>("L", false, "");
00208   validPL->set<double>("phi", false, "");
00209   validPL->set<double>("Q", false, "");
00210   }
00211 
00212   if (matModel == "GursonSD")
00213   {
00214   validPL->set<double>("f0", false, "");
00215   validPL->set<double>("Y0", false, "");
00216   validPL->set<double>("kw", false, "");
00217   validPL->set<double>("N", false, "");
00218   validPL->set<double>("q1", false, "");
00219   validPL->set<double>("q2", false, "");
00220   validPL->set<double>("q3", false, "");
00221   validPL->set<double>("eN", false, "");
00222   validPL->set<double>("sN", false, "");
00223   validPL->set<double>("fN", false, "");
00224   validPL->set<double>("fc", false, "");
00225   validPL->set<double>("ff", false, "");
00226   validPL->set<double>("flag", false, "");
00227   }  
00228 
00229   return validPL;
00230 }
00231 
00232 void
00233 Albany::ElasticityProblem::getAllocatedStates(
00234    Teuchos::ArrayRCP<Teuchos::ArrayRCP<Teuchos::RCP<Intrepid::FieldContainer<RealType> > > > oldState_,
00235    Teuchos::ArrayRCP<Teuchos::ArrayRCP<Teuchos::RCP<Intrepid::FieldContainer<RealType> > > > newState_
00236    ) const
00237 {
00238   oldState_ = oldState;
00239   newState_ = newState;
00240 }

Generated on Wed Mar 26 2014 18:36:38 for Albany: a Trilinos-based PDE code by  doxygen 1.7.1