00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032 #include "Teuchos_TestForException.hpp"
00033 #include "Phalanx_DataLayout.hpp"
00034
00035
00036 template<typename EvalT, typename Traits>
00037 FEInterpolation<EvalT, Traits>::
00038 FEInterpolation(const Teuchos::ParameterList& p) :
00039 val_node(p.get<std::string>("Node Variable Name"),
00040 p.get< Teuchos::RCP<PHX::DataLayout> >("Node Data Layout") ),
00041 val_qp(p.get<std::string>("QP Variable Name"),
00042 p.get< Teuchos::RCP<PHX::DataLayout> >("QP Scalar Data Layout") ),
00043 val_grad_qp(p.get<std::string>("Gradient QP Variable Name"),
00044 p.get< Teuchos::RCP<PHX::DataLayout> >("QP Vector Data Layout") )
00045 {
00046 this->addDependentField(val_node);
00047 this->addEvaluatedField(val_qp);
00048 this->addEvaluatedField(val_grad_qp);
00049
00050 this->setName("FEInterpolation");
00051 }
00052
00053
00054 template<typename EvalT, typename Traits>
00055 void FEInterpolation<EvalT, Traits>::
00056 postRegistrationSetup(PHX::FieldManager<Traits>& fm)
00057 {
00058 this->utils.setFieldData(val_node,fm);
00059 this->utils.setFieldData(val_qp,fm);
00060 this->utils.setFieldData(val_grad_qp,fm);
00061
00062
00063 typename std::vector< typename PHX::template MDField<ScalarT,PHX::NaturalOrder,Cell,Node>::size_type > dims;
00064 val_node.dimensions(dims);
00065 num_nodes = dims[1];
00066
00067 val_grad_qp.dimensions(dims);
00068 num_qp = dims[1];
00069 num_dim = dims[2];
00070 }
00071
00072
00073 template<typename EvalT, typename Traits>
00074 void FEInterpolation<EvalT, Traits>::
00075 evaluateFields(typename Traits::EvalData cell_data)
00076 {
00077 using namespace PHX;
00078
00079 std::vector<MyCell>::iterator cell_it = cell_data.begin;
00080
00081
00082 for (std::size_t cell = 0; cell < cell_data.num_cells; ++cell) {
00083
00084 Array<double,NaturalOrder,QuadPoint,Node>& phi =
00085 cell_it->getBasisFunctions();
00086
00087 Array<double,NaturalOrder,QuadPoint,Node,Dim>& grad_phi =
00088 cell_it->getBasisFunctionGradients();
00089
00090
00091 for (std::size_t qp = 0; qp < num_qp; ++qp) {
00092
00093 val_qp(cell,qp) = 0.0;
00094
00095 for (std::size_t dim = 0; dim < num_dim; ++dim)
00096 val_grad_qp(cell,qp,dim) = 0.0;
00097
00098
00099 for (std::size_t node = 0; node < num_nodes; ++node) {
00100
00101 val_qp(cell,qp) += phi(qp,node) * val_node(cell,node);
00102
00103 for (std::size_t dim = 0; dim < num_dim; ++dim)
00104 val_grad_qp(cell,qp,dim) +=
00105 grad_phi(qp,node,dim) * val_node(cell,node);
00106
00107 }
00108 }
00109
00110 ++cell_it;
00111
00112 }
00113
00114 }
00115
00116