43 #include "PanzerAdaptersSTK_config.hpp"
44 #ifdef PANZER_HAVE_TEKO
48 namespace panzer_stk {
53 ParameterListCallbackBlocked::ParameterListCallbackBlocked(
57 : connManager_(connManager), blocked_ugi_(blocked_ugi), aux_blocked_ugi_(aux_blocked_ugi)
61 ParameterListCallbackBlocked::request(
const Teko::RequestMesg & rm)
69 for(itr=inputPL->
begin();itr!=inputPL->
end();++itr) {
70 std::string * str_ptr = 0;
72 setFieldByKey(itr->first,field,*outputPL);
78 bool ParameterListCallbackBlocked::handlesRequest(
const Teko::RequestMesg & rm)
82 if(rm.getName()==
"Parameter List") {
83 bool isHandled =
true;
86 if(pl->
isType<std::string>(
"x-coordinates")) {
87 field = pl->
get<std::string>(
"x-coordinates");
92 if(pl->
isType<std::string>(
"y-coordinates")) {
94 if(field != pl->
get<std::string>(
"y-coordinates")) {
98 if(pl->
isType<std::string>(
"z-coordinates")) {
100 if(field != pl->
get<std::string>(
"z-coordinates")) {
104 if(pl->
isType<std::string>(
"Coordinates")){
105 field = pl->
get<std::string>(
"Coordinates");
107 #ifdef PANZER_HAVE_EPETRA_STACK
108 if(pl->
isType<std::string>(
"Coordinates-Epetra")){
109 field = pl->
get<std::string>(
"Coordinates-Epetra");
118 void ParameterListCallbackBlocked::preRequest(
const Teko::RequestMesg & rm)
122 const std::string&
field(getHandledField(*rm.getParameterList()));
127 std::vector<Teuchos::RCP<panzer::GlobalIndexer>>
128 fieldDOFMngrs = blocked_ugi_->getFieldDOFManagers();
129 for (
int b(0); b < static_cast<int>(fieldDOFMngrs.size()); ++b)
131 for (
int f(0); f < fieldDOFMngrs[b]->getNumFields(); ++f)
133 if (fieldDOFMngrs[b]->getFieldString(f) ==
field)
141 aux_blocked_ugi_->getFieldBlock(aux_blocked_ugi_->getFieldNum(field));
143 block = blocked_ugi_->getFieldBlock(blocked_ugi_->getFieldNum(field));
146 #ifdef PANZER_HAVE_EPETRA_STACK
147 if (rm.getParameterList()->isType<std::string>(
"Coordinates-Epetra")) {
148 buildArrayToVectorEpetra(block, field, useAux);
149 buildCoordinatesEpetra(field, useAux);
153 buildArrayToVectorTpetra(block, field, useAux);
154 buildCoordinatesTpetra(field, useAux);
158 void ParameterListCallbackBlocked::setFieldByKey(
const std::string & key,
const std::string & field,
Teuchos::ParameterList & pl)
const
161 if(key==
"x-coordinates") {
162 double * x =
const_cast<double *
>(&getCoordinateByField(0,field)[0]);
163 pl.
set<
double*>(key,x);
165 else if(key==
"y-coordinates") {
166 double * y =
const_cast<double *
>(&getCoordinateByField(1,field)[0]);
167 pl.
set<
double*>(key,y);
169 else if(key==
"z-coordinates") {
170 double * z =
const_cast<double *
>(&getCoordinateByField(2,field)[0]);
171 pl.
set<
double*>(key,z);
172 }
else if(key ==
"Coordinates") {
174 #ifdef PANZER_HAVE_EPETRA_STACK
175 }
else if(key ==
"Coordinates-Epetra") {
182 "ParameterListCallback cannot handle key=\"" << key <<
"\"");
185 void ParameterListCallbackBlocked::buildArrayToVectorTpetra(
int block,
const std::string & field,
const bool useAux)
187 if(arrayToVectorTpetra_[field]==Teuchos::null) {
190 ugi = aux_blocked_ugi_->getFieldDOFManagers()[block];
192 ugi = blocked_ugi_->getFieldDOFManagers()[block];
197 #ifdef PANZER_HAVE_EPETRA_STACK
198 void ParameterListCallbackBlocked::buildArrayToVectorEpetra(
int block,
const std::string & field,
const bool useAux)
200 if(arrayToVectorEpetra_[field]==Teuchos::null) {
203 ugi = aux_blocked_ugi_->getFieldDOFManagers()[block];
205 ugi = blocked_ugi_->getFieldDOFManagers()[block];
211 void ParameterListCallbackBlocked::buildCoordinatesTpetra(
const std::string & field,
const bool useAux)
213 std::map<std::string,Kokkos::DynRankView<double,PHX::Device> > data;
217 std::vector<std::string> elementBlocks;
219 aux_blocked_ugi_->getElementBlockIds(elementBlocks);
221 blocked_ugi_->getElementBlockIds(elementBlocks);
222 for(std::size_t i=0;i<elementBlocks.size();++i) {
223 std::string blockId = elementBlocks[i];
224 std::vector<std::size_t> localCellIds;
227 Kokkos::DynRankView<double,PHX::Device> & fieldData = data[blockId];
228 fieldData = Kokkos::DynRankView<double,PHX::Device>(
"fieldData",connManager_->getElementBlock(blockId).size(),fieldPattern->
numberIds());
232 connManager_->getDofCoords(blockId,*fieldPattern,localCellIds,fieldData);
236 out.setOutputToRootOnly(-1);
237 out <<
"WARNING: In ParameterListCallback::buildCoordinates(), the Intrepid2::FieldPattern in "
238 <<
"block \"" << blockId <<
"\" does not support interpolatory coordinates. "
239 <<
"This may be fine if coordinates are not actually needed. However if they are then bad things "
240 <<
"will happen. Enjoy!" << std::endl;
246 coordsVecTp_ = arrayToVectorTpetra_[
field]->template getDataVector<double>(
field,data);
248 switch(coordsVecTp_->getNumVectors()) {
250 zcoords_[
field].resize(coordsVecTp_->getLocalLength());
251 coordsVecTp_->getVector(2)->get1dCopy(Teuchos::arrayViewFromVector(zcoords_[field]));
254 ycoords_[
field].resize(coordsVecTp_->getLocalLength());
255 coordsVecTp_->getVector(1)->get1dCopy(Teuchos::arrayViewFromVector(ycoords_[field]));
258 xcoords_[
field].resize(coordsVecTp_->getLocalLength());
259 coordsVecTp_->getVector(0)->get1dCopy(Teuchos::arrayViewFromVector(xcoords_[field]));
263 "ParameterListCallback::buildCoordinates: Constructed multivector has nonphysical dimensions.");
268 #ifdef PANZER_HAVE_EPETRA_STACK
269 void ParameterListCallbackBlocked::buildCoordinatesEpetra(
const std::string & field,
const bool useAux)
271 std::map<std::string,Kokkos::DynRankView<double,PHX::Device> > data;
275 std::vector<std::string> elementBlocks;
277 aux_blocked_ugi_->getElementBlockIds(elementBlocks);
279 blocked_ugi_->getElementBlockIds(elementBlocks);
280 for(std::size_t i=0;i<elementBlocks.size();++i) {
281 std::string blockId = elementBlocks[i];
282 std::vector<std::size_t> localCellIds;
285 Kokkos::DynRankView<double,PHX::Device> & fieldData = data[blockId];
286 fieldData = Kokkos::DynRankView<double,PHX::Device>(
"fieldData",connManager_->getElementBlock(blockId).size(),fieldPattern->
numberIds());
290 connManager_->getDofCoords(blockId,*fieldPattern,localCellIds,fieldData);
294 out.setOutputToRootOnly(-1);
295 out <<
"WARNING: In ParameterListCallback::buildCoordinates(), the Intrepid2::FieldPattern in "
296 <<
"block \"" << blockId <<
"\" does not support interpolatory coordinates. "
297 <<
"This may be fine if coordinates are not actually needed. However if they are then bad things "
298 <<
"will happen. Enjoy!" << std::endl;
304 coordsVecEp_ = arrayToVectorEpetra_[
field]->template getDataVector<double>(
field,data);
308 std::string ParameterListCallbackBlocked::
312 if(pl.
isType<std::string>(
"x-coordinates"))
313 return pl.
get<std::string>(
"x-coordinates");
314 else if(pl.
isType<std::string>(
"Coordinates"))
315 return pl.
get<std::string>(
"Coordinates");
316 #ifdef PANZER_HAVE_EPETRA_STACK
317 else if(pl.
isType<std::string>(
"Coordinates-Epetra"))
318 return pl.
get<std::string>(
"Coordinates-Epetra");
321 #ifdef PANZER_HAVE_EPETRA_STACK
322 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"Neither x-coordinates nor Coordinates or Coordinates-Epetra field provided.");
328 const std::vector<double> & ParameterListCallbackBlocked::
329 getCoordinateByField(
int dim,
const std::string & field)
const
335 const std::map<std::string,std::vector<double> > * coord =
nullptr;
336 if(dim==0) coord = &xcoords_;
337 else if(dim==1) coord = &ycoords_;
338 else if(dim==2) coord = &zcoords_;
341 "ParameterListCallbackBlocked::getCoordinateByField: dim is not in range of [0,2] for field \""
345 std::map<std::string,std::vector<double> >::const_iterator itr;
346 itr = coord->find(field);
349 "ParameterListCallbackBlocked::getCoordinateByField: Coordinates for field \"" + field +
350 "\" dimension " << dim <<
" have not been built!");
356 ::getFieldPattern(
const std::string & fieldName,
const bool useAux)
const
358 std::vector<std::string> elementBlocks;
360 aux_blocked_ugi_->getElementBlockIds(elementBlocks);
362 blocked_ugi_->getElementBlockIds(elementBlocks);
364 for(std::size_t e=0;e<elementBlocks.size();e++) {
365 std::string blockId = elementBlocks[e];
367 if(blocked_ugi_->fieldInBlock(fieldName,blockId))
368 return Teuchos::rcp_dynamic_cast<const panzer::Intrepid2FieldPattern>(blocked_ugi_->getFieldPattern(blockId,fieldName),
true);
370 if(aux_blocked_ugi_->fieldInBlock(fieldName,blockId))
371 return Teuchos::rcp_dynamic_cast<const panzer::Intrepid2FieldPattern>(aux_blocked_ugi_->getFieldPattern(blockId,fieldName),
true);
374 return Teuchos::null;
ConstIterator end() const
T & get(const std::string &name, T def_value)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
T & getValue(T *ptr) const
virtual int numberIds() const
bool supportsInterpolatoryCoordinates() const
Does this field pattern support interpolatory coordinates?
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
params_t::ConstIterator ConstIterator
ConstIterator begin() const
const ParameterEntry & entry(ConstIterator i) const
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field
A field to which we'll contribute, or in which we'll store, the result of computing this integral...
bool isType(const std::string &name) const
#define TEUCHOS_ASSERT(assertion_test)