Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_Response_Probe_impl.hpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Panzer: A partial differential equation assembly
5// engine for strongly coupled complex multiphysics systems
6// Copyright (2011) Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact Roger P. Pawlowski (rppawlo@sandia.gov) and
39// Eric C. Cyr (eccyr@sandia.gov)
40// ***********************************************************************
41// @HEADER
42
43#ifndef __Panzer_Response_Probe_impl_hpp__
44#define __Panzer_Response_Probe_impl_hpp__
45
46#include "Teuchos_Comm.hpp"
47#include "Teuchos_CommHelpers.hpp"
48#include "Teuchos_dyn_cast.hpp"
49
50#include "PanzerDiscFE_config.hpp"
51#ifdef PANZER_HAVE_EPETRA_STACK
52#include "Epetra_LocalMap.h"
53#endif
54
55#include "Sacado_Traits.hpp"
56
57namespace panzer {
58
59template <typename EvalT>
61Response_Probe(const std::string & responseName, MPI_Comm comm,
62 const Teuchos::RCP<const panzer::LinearObjFactory<panzer::Traits> > & linObjFact)
63 : ResponseMESupport_Default<EvalT>(responseName,comm), value(0.0),
64 have_probe(false), linObjFactory_(linObjFact)
65{
66 if(linObjFactory_!=Teuchos::null) {
67 // requires thyra object factory
68 thyraObjFactory_ = Teuchos::rcp_dynamic_cast<const panzer::ThyraObjFactory<double> >(linObjFactory_,true);
69 setSolnVectorSpace(thyraObjFactory_->getThyraDomainSpace());
70
71 // build a ghosted container, with a solution vector
72 ghostedContainer_ = linObjFactory_->buildGhostedLinearObjContainer();
73
74 // set ghosted container (work space for assembly)
76 }
77}
78
79template <typename EvalT>
82{
83 value = 0.0;
84 have_probe = false;
85
86 if(ghostedContainer_!=Teuchos::null) ghostedContainer_->initialize();
87}
88
89template <typename EvalT>
92{
93 double glbValue = Sacado::scalarValue(value);
94
95 // find the minimum processor who has the probe value
96 int locProc = have_probe ? this->getComm()->getRank() : this->getComm()->getSize();
97 int glbProc = 0;
98 Teuchos::reduceAll(*this->getComm(), Teuchos::REDUCE_MIN, Thyra::Ordinal(1), &locProc, &glbProc);
99
100 TEUCHOS_ASSERT(glbProc < this->getComm()->getSize());
101
102 // now broadcast the value from proc glbProc
103 Teuchos::broadcast(*this->getComm(), glbProc, Thyra::Ordinal(1), &glbValue);
104
105 value = glbValue;
106
107 // built data in vectors
108#ifdef PANZER_HAVE_EPETRA_STACK
109 if(this->useEpetra()) {
110 // use epetra
111 this->getEpetraVector()[0] = glbValue;
112 }
113 else
114#endif
115 {
116 // use thyra
117 TEUCHOS_ASSERT(this->useThyra());
118
119 this->getThyraVector()[0] = glbValue;
120 }
121}
122
123template < >
126{
127 using Teuchos::rcp_dynamic_cast;
128
129 Teuchos::RCP<Thyra::MultiVectorBase<double> > dgdx_unique = getDerivative();
130
131 uniqueContainer_ = linObjFactory_->buildLinearObjContainer();
132 Teuchos::rcp_dynamic_cast<ThyraObjContainer<double> >(uniqueContainer_)->set_x_th(dgdx_unique->col(0));
133
134 linObjFactory_->ghostToGlobalContainer(*ghostedContainer_,*uniqueContainer_,LinearObjContainer::X);
135
136 uniqueContainer_ = Teuchos::null;
137}
138
139#ifdef Panzer_BUILD_HESSIAN_SUPPORT
140template < >
143{
144 using Teuchos::rcp_dynamic_cast;
145
146 Teuchos::RCP<Thyra::MultiVectorBase<double> > dgdx_unique = getDerivative();
147
148 uniqueContainer_ = linObjFactory_->buildLinearObjContainer();
149 Teuchos::rcp_dynamic_cast<ThyraObjContainer<double> >(uniqueContainer_)->set_x_th(dgdx_unique->col(0));
150
151 linObjFactory_->ghostToGlobalContainer(*ghostedContainer_,*uniqueContainer_,LinearObjContainer::X);
152
153 uniqueContainer_ = Teuchos::null;
154}
155#endif
156
157template < >
160{
161 const int n = value.size();
162 const int num_deriv = this->numDeriv();
163 TEUCHOS_ASSERT(n == 0 || n == num_deriv);
164 if (n == 0)
165 value.resize(num_deriv);
166
167 // find the minimum processor who has the probe value
168 if (num_deriv > 0) {
169 int locProc = have_probe ? this->getComm()->getRank() : this->getComm()->getSize();
170 int glbProc = 0;
171 Teuchos::reduceAll(*this->getComm(), Teuchos::REDUCE_MIN, Thyra::Ordinal(1), &locProc, &glbProc);
172
173 TEUCHOS_ASSERT(glbProc < this->getComm()->getSize());
174
175 // now broadcast the derivatives from proc glbProc
176 Teuchos::broadcast(*this->getComm(), glbProc, Thyra::Ordinal(num_deriv), &value.fastAccessDx(0));
177 }
178
179 // copy data in vectors
180#ifdef PANZER_HAVE_EPETRA_STACK
181 if(this->useEpetra()) {
182 // use epetra
183 Epetra_MultiVector& deriv = this->getEpetraMultiVector();
184 for (int i=0; i<num_deriv; ++i)
185 deriv[i][0] = value.dx(i);
186 }
187 else
188#endif
189 {
190 // use thyra
191 TEUCHOS_ASSERT(this->useThyra());
192 Thyra::ArrayRCP< Thyra::ArrayRCP<double> > deriv = this->getThyraMultiVector();
193 for (int i=0; i<num_deriv; ++i)
194 deriv[i][0] = value.dx(i);
195 }
196}
197
198// Do nothing unless derivatives are actually required
199template <typename EvalT>
201setSolnVectorSpace(const Teuchos::RCP<const Thyra::VectorSpaceBase<double> > & /* soln_vs */) { }
202
203// derivatives are required for
204template < >
206setSolnVectorSpace(const Teuchos::RCP<const Thyra::VectorSpaceBase<double> > & soln_vs)
207{
208 setDerivativeVectorSpace(soln_vs);
209}
210
211#ifdef Panzer_BUILD_HESSIAN_SUPPORT
212// derivatives are required for
213template < >
215setSolnVectorSpace(const Teuchos::RCP<const Thyra::VectorSpaceBase<double> > & soln_vs)
216{
217 setDerivativeVectorSpace(soln_vs);
218}
219#endif
220
221// Do nothing unless derivatives are required
222template <typename EvalT>
224adjustForDirichletConditions(const GlobalEvaluationData & /* localBCRows */, const GlobalEvaluationData & /* globalBCRows */) { }
225
226// Do nothing unless derivatives are required
227template < >
229adjustForDirichletConditions(const GlobalEvaluationData & localBCRows,const GlobalEvaluationData & globalBCRows)
230{
231 linObjFactory_->adjustForDirichletConditions(Teuchos::dyn_cast<const LinearObjContainer>(localBCRows),
232 Teuchos::dyn_cast<const LinearObjContainer>(globalBCRows),
233 *ghostedContainer_,true,true);
234}
235
236#ifdef Panzer_BUILD_HESSIAN_SUPPORT
237// Do nothing unless derivatives are required
238template < >
240adjustForDirichletConditions(const GlobalEvaluationData & localBCRows,const GlobalEvaluationData & globalBCRows)
241{
242 linObjFactory_->adjustForDirichletConditions(Teuchos::dyn_cast<const LinearObjContainer>(localBCRows),
243 Teuchos::dyn_cast<const LinearObjContainer>(globalBCRows),
244 *ghostedContainer_,true,true);
245}
246#endif
247
248}
249
250#endif
Teuchos::RCP< const panzer::ThyraObjFactory< double > > thyraObjFactory_
Teuchos::RCP< LinearObjContainer > ghostedContainer_
void setSolnVectorSpace(const Teuchos::RCP< const Thyra::VectorSpaceBase< double > > &soln_vs)
Set solution vector space.
Teuchos::RCP< const panzer::LinearObjFactory< panzer::Traits > > linObjFactory_
void adjustForDirichletConditions(const GlobalEvaluationData &localBCRows, const GlobalEvaluationData &globalBCRows)
virtual void scatterResponse()
This simply does global summation, then shoves the result into a vector.