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 Data Layout") ),
00043 val_grad_qp(p.get<std::string>("Gradient QP Variable Name"),
00044 p.get< Teuchos::RCP<PHX::DataLayout> >("QP 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
00064 template<typename EvalT, typename Traits>
00065 void FEInterpolation<EvalT, Traits>::
00066 evaluateFields(typename Traits::EvalData cell_data)
00067 {
00068
00069 const int nodes_per_cell = val_node.fieldTag().dataLayout().size();
00070 const int qp_per_cell = val_qp.fieldTag().dataLayout().size();
00071 std::vector<MyCell>::iterator cell_it = cell_data.begin;
00072
00073
00074 for (std::size_t cell = 0; cell < cell_data.num_cells; ++cell) {
00075
00076 std::vector< std::vector<double> >& phi =
00077 cell_it->getBasisFunctions();
00078 std::vector< std::vector< MyVector<double> > >& grad_phi =
00079 cell_it->getBasisFunctionGradients();
00080
00081 int node_offset = cell * nodes_per_cell;
00082 int qp_offset = cell * qp_per_cell;
00083
00084
00085 for (int qp = 0; qp < qp_per_cell; ++qp) {
00086
00087 val_qp[qp_offset + qp] = 0.0;
00088 val_grad_qp[qp_offset + qp] = MyVector<ScalarT>(0.0, 0.0, 0.0);
00089
00090
00091 for (int node = 0; node < nodes_per_cell; ++node) {
00092
00093 val_qp[qp_offset + qp] += phi[qp][node] * val_node[node_offset + node];
00094
00095 val_grad_qp[qp_offset + qp] +=
00096 grad_phi[qp][node] * val_node[node_offset + node];
00097 }
00098 }
00099
00100 ++cell_it;
00101
00102 }
00103
00104 }
00105
00106