00001
00002
00003
00004
00005
00006
00007 #include "Teuchos_TestForException.hpp"
00008 #include "Phalanx_DataLayout.hpp"
00009
00010
00011
00012
00013 namespace PHAL {
00014
00015 template<typename EvalT, typename Traits>
00016 SeparableScatterScalarResponseBase<EvalT, Traits>::
00017 SeparableScatterScalarResponseBase(const Teuchos::ParameterList& p,
00018 const Teuchos::RCP<Albany::Layouts>& dl)
00019 {
00020 setup(p, dl);
00021 }
00022
00023 template<typename EvalT, typename Traits>
00024 void SeparableScatterScalarResponseBase<EvalT, Traits>::
00025 postRegistrationSetup(typename Traits::SetupData d,
00026 PHX::FieldManager<Traits>& fm)
00027 {
00028 this->utils.setFieldData(local_response,fm);
00029 }
00030
00031 template<typename EvalT, typename Traits>
00032 void
00033 SeparableScatterScalarResponseBase<EvalT, Traits>::
00034 setup(const Teuchos::ParameterList& p, const Teuchos::RCP<Albany::Layouts>& dl)
00035 {
00036 bool stand_alone = p.get<bool>("Stand-alone Evaluator");
00037
00038
00039 PHX::Tag<ScalarT> local_response_tag =
00040 p.get<PHX::Tag<ScalarT> >("Local Response Field Tag");
00041 local_response = PHX::MDField<ScalarT>(local_response_tag);
00042 if (stand_alone)
00043 this->addDependentField(local_response);
00044 else
00045 this->addEvaluatedField(local_response);
00046 }
00047
00048
00049
00050
00051
00052 template<typename Traits>
00053 SeparableScatterScalarResponse<PHAL::AlbanyTraits::Jacobian, Traits>::
00054 SeparableScatterScalarResponse(const Teuchos::ParameterList& p,
00055 const Teuchos::RCP<Albany::Layouts>& dl)
00056 {
00057 this->setup(p,dl);
00058 }
00059
00060 template<typename Traits>
00061 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::Jacobian, Traits>::
00062 preEvaluate(typename Traits::PreEvalData workset)
00063 {
00064
00065 Teuchos::RCP<Epetra_MultiVector> dgdx = workset.dgdx;
00066 Teuchos::RCP<Epetra_MultiVector> overlapped_dgdx = workset.overlapped_dgdx;
00067 if (dgdx != Teuchos::null) {
00068 dgdx->PutScalar(0.0);
00069 overlapped_dgdx->PutScalar(0.0);
00070 }
00071
00072 Teuchos::RCP<Epetra_MultiVector> dgdxdot = workset.dgdxdot;
00073 Teuchos::RCP<Epetra_MultiVector> overlapped_dgdxdot =
00074 workset.overlapped_dgdxdot;
00075 if (dgdxdot != Teuchos::null) {
00076 dgdxdot->PutScalar(0.0);
00077 overlapped_dgdxdot->PutScalar(0.0);
00078 }
00079 }
00080
00081 template<typename Traits>
00082 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::Jacobian, Traits>::
00083 evaluateFields(typename Traits::EvalData workset)
00084 {
00085
00086 Teuchos::RCP<Epetra_MultiVector> dgdx = workset.overlapped_dgdx;
00087 Teuchos::RCP<Epetra_MultiVector> dgdxdot = workset.overlapped_dgdxdot;
00088 Teuchos::RCP<Epetra_MultiVector> dg;
00089 if (dgdx != Teuchos::null)
00090 dg = dgdx;
00091 else
00092 dg = dgdxdot;
00093
00094
00095 for (std::size_t cell=0; cell < workset.numCells; ++cell) {
00096 const Teuchos::ArrayRCP<Teuchos::ArrayRCP<int> >& nodeID =
00097 workset.wsElNodeEqID[cell];
00098
00099
00100 for (std::size_t res = 0; res < this->global_response.size(); res++) {
00101 ScalarT& val = this->local_response(cell, res);
00102
00103
00104 for (unsigned int node_dof=0; node_dof<numNodes; node_dof++) {
00105 int neq = nodeID[node_dof].size();
00106
00107
00108 for (unsigned int eq_dof=0; eq_dof<neq; eq_dof++) {
00109
00110
00111 int deriv = neq * node_dof + eq_dof;
00112
00113
00114 int dof = nodeID[node_dof][eq_dof];
00115
00116
00117 dg->SumIntoMyValue(dof, res, val.dx(deriv));
00118
00119 }
00120 }
00121 }
00122 }
00123 }
00124
00125 template<typename Traits>
00126 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::Jacobian, Traits>::
00127 postEvaluate(typename Traits::PostEvalData workset)
00128 {
00129
00130 Teuchos::RCP<Epetra_Vector> g = workset.g;
00131 if (g != Teuchos::null)
00132 for (std::size_t res = 0; res < this->global_response.size(); res++) {
00133 (*g)[res] = this->global_response[res].val();
00134 }
00135
00136
00137 Teuchos::RCP<Epetra_MultiVector> dgdx = workset.dgdx;
00138 Teuchos::RCP<Epetra_MultiVector> overlapped_dgdx = workset.overlapped_dgdx;
00139 if (dgdx != Teuchos::null)
00140 dgdx->Export(*overlapped_dgdx, *workset.x_importer, Add);
00141
00142 Teuchos::RCP<Epetra_MultiVector> dgdxdot = workset.dgdxdot;
00143 Teuchos::RCP<Epetra_MultiVector> overlapped_dgdxdot =
00144 workset.overlapped_dgdxdot;
00145 if (dgdxdot != Teuchos::null)
00146 dgdxdot->Export(*overlapped_dgdxdot, *workset.x_importer, Add);
00147 }
00148
00149
00150
00151
00152 #ifdef ALBANY_SG_MP
00153 template<typename Traits>
00154 SeparableScatterScalarResponse<PHAL::AlbanyTraits::SGJacobian, Traits>::
00155 SeparableScatterScalarResponse(const Teuchos::ParameterList& p,
00156 const Teuchos::RCP<Albany::Layouts>& dl)
00157 {
00158 this->setup(p,dl);
00159 }
00160
00161 template<typename Traits>
00162 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::SGJacobian, Traits>::
00163 preEvaluate(typename Traits::PreEvalData workset)
00164 {
00165
00166 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdx_sg =
00167 workset.sg_dgdx;
00168 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> overlapped_dgdx_sg =
00169 workset.overlapped_sg_dgdx;
00170 if (dgdx_sg != Teuchos::null) {
00171 dgdx_sg->init(0.0);
00172 overlapped_dgdx_sg->init(0.0);
00173 }
00174
00175 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdxdot_sg =
00176 workset.sg_dgdxdot;
00177 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> overlapped_dgdxdot_sg =
00178 workset.overlapped_sg_dgdxdot;
00179 if (dgdxdot_sg != Teuchos::null) {
00180 dgdxdot_sg->init(0.0);
00181 overlapped_dgdxdot_sg->init(0.0);
00182 }
00183 }
00184
00185 template<typename Traits>
00186 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::SGJacobian, Traits>::
00187 evaluateFields(typename Traits::EvalData workset)
00188 {
00189
00190 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdx_sg =
00191 workset.overlapped_sg_dgdx;
00192 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdxdot_sg =
00193 workset.overlapped_sg_dgdxdot;
00194 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dg_sg;
00195 if (dgdx_sg != Teuchos::null)
00196 dg_sg = dgdx_sg;
00197 else
00198 dg_sg = dgdxdot_sg;
00199
00200
00201 for (std::size_t cell=0; cell < workset.numCells; ++cell) {
00202 const Teuchos::ArrayRCP<Teuchos::ArrayRCP<int> >& nodeID =
00203 workset.wsElNodeEqID[cell];
00204
00205
00206 for (std::size_t res = 0; res < this->global_response.size(); res++) {
00207 ScalarT& val = this->local_response(cell, res);
00208
00209
00210 for (unsigned int node_dof=0; node_dof<numNodes; node_dof++) {
00211 int neq = nodeID[node_dof].size();
00212
00213
00214 for (unsigned int eq_dof=0; eq_dof<neq; eq_dof++) {
00215
00216
00217 int deriv = neq * node_dof + eq_dof;
00218
00219
00220 int dof = nodeID[node_dof][eq_dof];
00221
00222
00223 for (int block=0; block<dg_sg->size(); block++)
00224 (*dg_sg)[block].SumIntoMyValue(dof, res,
00225 val.dx(deriv).coeff(block));
00226
00227 }
00228 }
00229 }
00230 }
00231 }
00232
00233 template<typename Traits>
00234 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::SGJacobian, Traits>::
00235 postEvaluate(typename Traits::PostEvalData workset)
00236 {
00237
00238 Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > g_sg = workset.sg_g;
00239 if (g_sg != Teuchos::null) {
00240 for (std::size_t res = 0; res < this->global_response.size(); res++) {
00241 ScalarT& val = this->global_response[res];
00242 for (int block=0; block<g_sg->size(); block++)
00243 (*g_sg)[block][res] = val.val().coeff(block);
00244 }
00245 }
00246
00247
00248 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdx_sg =
00249 workset.sg_dgdx;
00250 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> overlapped_dgdx_sg =
00251 workset.overlapped_sg_dgdx;
00252 if (dgdx_sg != Teuchos::null)
00253 for (int block=0; block<dgdx_sg->size(); block++)
00254 (*dgdx_sg)[block].Export((*overlapped_dgdx_sg)[block],
00255 *workset.x_importer, Add);
00256
00257 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdxdot_sg =
00258 workset.sg_dgdxdot;
00259 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> overlapped_dgdxdot_sg =
00260 workset.overlapped_sg_dgdxdot;
00261 if (dgdxdot_sg != Teuchos::null)
00262 for (int block=0; block<dgdxdot_sg->size(); block++)
00263 (*dgdxdot_sg)[block].Export((*overlapped_dgdxdot_sg)[block],
00264 *workset.x_importer, Add);
00265 }
00266
00267
00268
00269
00270
00271 template<typename Traits>
00272 SeparableScatterScalarResponse<PHAL::AlbanyTraits::MPJacobian, Traits>::
00273 SeparableScatterScalarResponse(const Teuchos::ParameterList& p,
00274 const Teuchos::RCP<Albany::Layouts>& dl)
00275 {
00276 this->setup(p,dl);
00277 }
00278
00279 template<typename Traits>
00280 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::MPJacobian, Traits>::
00281 preEvaluate(typename Traits::PreEvalData workset)
00282 {
00283
00284 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> dgdx_mp =
00285 workset.mp_dgdx;
00286 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> overlapped_dgdx_mp =
00287 workset.overlapped_mp_dgdx;
00288 if (dgdx_mp != Teuchos::null) {
00289 dgdx_mp->init(0.0);
00290 overlapped_dgdx_mp->init(0.0);
00291 }
00292
00293 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> dgdxdot_mp =
00294 workset.mp_dgdxdot;
00295 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> overlapped_dgdxdot_mp =
00296 workset.overlapped_mp_dgdxdot;
00297 if (dgdxdot_mp != Teuchos::null) {
00298 dgdxdot_mp->init(0.0);
00299 overlapped_dgdxdot_mp->init(0.0);
00300 }
00301 }
00302
00303 template<typename Traits>
00304 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::MPJacobian, Traits>::
00305 evaluateFields(typename Traits::EvalData workset)
00306 {
00307
00308 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> dgdx_mp =
00309 workset.overlapped_mp_dgdx;
00310 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> dgdxdot_mp =
00311 workset.overlapped_mp_dgdxdot;
00312 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> dg_mp;
00313 if (dgdx_mp != Teuchos::null)
00314 dg_mp = dgdx_mp;
00315 else
00316 dg_mp = dgdxdot_mp;
00317
00318
00319 for (std::size_t cell=0; cell < workset.numCells; ++cell) {
00320 const Teuchos::ArrayRCP<Teuchos::ArrayRCP<int> >& nodeID =
00321 workset.wsElNodeEqID[cell];
00322
00323
00324 for (std::size_t res = 0; res < this->global_response.size(); res++) {
00325 ScalarT& val = this->local_response(cell, res);
00326
00327
00328 for (unsigned int node_dof=0; node_dof<numNodes; node_dof++) {
00329 int neq = nodeID[node_dof].size();
00330
00331
00332 for (unsigned int eq_dof=0; eq_dof<neq; eq_dof++) {
00333
00334
00335 int deriv = neq * node_dof + eq_dof;
00336
00337
00338 int dof = nodeID[node_dof][eq_dof];
00339
00340
00341 for (int block=0; block<dg_mp->size(); block++)
00342 (*dg_mp)[block].SumIntoMyValue(dof, res,
00343 val.dx(deriv).coeff(block));
00344
00345 }
00346 }
00347 }
00348 }
00349 }
00350
00351 template<typename Traits>
00352 void SeparableScatterScalarResponse<PHAL::AlbanyTraits::MPJacobian, Traits>::
00353 postEvaluate(typename Traits::PostEvalData workset)
00354 {
00355
00356 Teuchos::RCP<Stokhos::ProductEpetraVector> g_mp = workset.mp_g;
00357 if (g_mp != Teuchos::null) {
00358 for (std::size_t res = 0; res < this->global_response.size(); res++) {
00359 ScalarT& val = this->global_response[res];
00360 for (int block=0; block<g_mp->size(); block++)
00361 (*g_mp)[block][res] = val.val().coeff(block);
00362 }
00363 }
00364
00365
00366 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> dgdx_mp =
00367 workset.mp_dgdx;
00368 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> overlapped_dgdx_mp =
00369 workset.overlapped_mp_dgdx;
00370 if (dgdx_mp != Teuchos::null)
00371 for (int block=0; block<dgdx_mp->size(); block++)
00372 (*dgdx_mp)[block].Export((*overlapped_dgdx_mp)[block],
00373 *workset.x_importer, Add);
00374
00375 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> dgdxdot_mp =
00376 workset.mp_dgdxdot;
00377 Teuchos::RCP<Stokhos::ProductEpetraMultiVector> overlapped_dgdxdot_mp =
00378 workset.overlapped_mp_dgdxdot;
00379 if (dgdxdot_mp != Teuchos::null)
00380 for (int block=0; block<dgdxdot_mp->size(); block++)
00381 (*dgdxdot_mp)[block].Export((*overlapped_dgdxdot_mp)[block],
00382 *workset.x_importer, Add);
00383 }
00384 #endif //ALBANY_SG_MP
00385
00386 }
00387