00001
00002
00003
00004
00005
00006
00007 #include "Teuchos_TestForException.hpp"
00008 #include "Phalanx_DataLayout.hpp"
00009
00010 #include "Intrepid_FunctionSpaceTools.hpp"
00011
00012 namespace PHAL {
00013
00014
00015 template<typename EvalT, typename Traits>
00016 NSNeutronEqResid<EvalT, Traits>::
00017 NSNeutronEqResid(const Teuchos::ParameterList& p) :
00018 wBF (p.get<std::string> ("Weighted BF Name"),
00019 p.get<Teuchos::RCP<PHX::DataLayout> >("Node QP Scalar Data Layout") ),
00020 Neutron (p.get<std::string> ("QP Variable Name"),
00021 p.get<Teuchos::RCP<PHX::DataLayout> >("QP Scalar Data Layout") ),
00022 NeutronDiff (p.get<std::string> ("Neutron Diffusion Name"),
00023 p.get<Teuchos::RCP<PHX::DataLayout> >("QP Scalar Data Layout") ),
00024 wGradBF (p.get<std::string> ("Weighted Gradient BF Name"),
00025 p.get<Teuchos::RCP<PHX::DataLayout> >("Node QP Vector Data Layout") ),
00026 NGrad (p.get<std::string> ("Gradient QP Variable Name"),
00027 p.get<Teuchos::RCP<PHX::DataLayout> >("QP Vector Data Layout") ),
00028 Absorp (p.get<std::string> ("Neutron Absorption Name"),
00029 p.get<Teuchos::RCP<PHX::DataLayout> >("QP Scalar Data Layout") ),
00030 Fission (p.get<std::string> ("Neutron Fission Name"),
00031 p.get<Teuchos::RCP<PHX::DataLayout> >("QP Scalar Data Layout") ),
00032 nu (p.get<std::string> ("Neutrons per Fission Name"),
00033 p.get<Teuchos::RCP<PHX::DataLayout> >("QP Scalar Data Layout") ),
00034 NResidual (p.get<std::string> ("Residual Name"),
00035 p.get<Teuchos::RCP<PHX::DataLayout> >("Node Scalar Data Layout") ),
00036 haveNeutSource (p.get<bool>("Have Neutron Source"))
00037 {
00038
00039 this->addDependentField(wBF);
00040 this->addDependentField(wGradBF);
00041 this->addDependentField(Neutron);
00042 this->addDependentField(NGrad);
00043 this->addDependentField(NeutronDiff);
00044 this->addDependentField(Absorp);
00045 this->addDependentField(Fission);
00046 this->addDependentField(nu);
00047
00048 if (haveNeutSource) {
00049 Source = PHX::MDField<ScalarT,Cell,QuadPoint>(
00050 p.get<std::string>("Source Name"),
00051 p.get<Teuchos::RCP<PHX::DataLayout> >("QP Scalar Data Layout") );
00052 this->addDependentField(Source);
00053 }
00054
00055 this->addEvaluatedField(NResidual);
00056
00057 Teuchos::RCP<PHX::DataLayout> vector_dl =
00058 p.get< Teuchos::RCP<PHX::DataLayout> >("QP Vector Data Layout");
00059 std::vector<PHX::DataLayout::size_type> dims;
00060 vector_dl->dimensions(dims);
00061 numQPs = dims[1];
00062 numDims = dims[2];
00063
00064
00065 flux.resize(dims[0], numQPs, numDims);
00066 abscoeff.resize(dims[0], numQPs);
00067
00068 this->setName("NSNeutronEqResid"+PHX::TypeString<EvalT>::value);
00069 }
00070
00071
00072 template<typename EvalT, typename Traits>
00073 void NSNeutronEqResid<EvalT, Traits>::
00074 postRegistrationSetup(typename Traits::SetupData d,
00075 PHX::FieldManager<Traits>& fm)
00076 {
00077 this->utils.setFieldData(wBF,fm);
00078 this->utils.setFieldData(wGradBF,fm);
00079 this->utils.setFieldData(Neutron,fm);
00080 this->utils.setFieldData(NGrad,fm);
00081 this->utils.setFieldData(NeutronDiff,fm);
00082 this->utils.setFieldData(Absorp,fm);
00083 this->utils.setFieldData(Fission,fm);
00084 this->utils.setFieldData(nu,fm);
00085 if (haveNeutSource) this->utils.setFieldData(Source,fm);
00086
00087 this->utils.setFieldData(NResidual,fm);
00088 }
00089
00090
00091 template<typename EvalT, typename Traits>
00092 void NSNeutronEqResid<EvalT, Traits>::
00093 evaluateFields(typename Traits::EvalData workset)
00094 {
00095 typedef Intrepid::FunctionSpaceTools FST;
00096
00097 FST::scalarMultiplyDataData<ScalarT> (flux, NeutronDiff, NGrad);
00098
00099 FST::integrate<ScalarT>(NResidual, flux, wGradBF, Intrepid::COMP_CPP, false);
00100
00101 for (std::size_t cell=0; cell < workset.numCells; ++cell) {
00102 for (std::size_t qp=0; qp < numQPs; ++qp) {
00103 abscoeff(cell,qp) =
00104 (Absorp(cell,qp) - nu(cell,qp)*Fission(cell,qp)) * Neutron(cell,qp);
00105 if (haveNeutSource) abscoeff(cell,qp) -= Source(cell,qp);
00106 }
00107 }
00108
00109 FST::integrate<ScalarT>(NResidual, abscoeff, wBF, Intrepid::COMP_CPP, true);
00110
00111 }
00112
00113
00114 }
00115
