Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_STK_PeriodicBC_MatchConditions.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_STK_PeriodicBC_MatchConditions_hpp__
44#define __Panzer_STK_PeriodicBC_MatchConditions_hpp__
45
46#include "Teuchos_Tuple.hpp"
47
48#include <vector>
49#include <string>
50
51namespace panzer_stk {
52
56 double error_;
57 int index_;
58 bool relative_; // compute error relative to length of domain
59 char labels_[3];
60
62 { labels_[0] = 'x'; labels_[1] = 'y'; labels_[2] = 'z'; }
63
64 void parseParams(const std::vector<std::string> & params)
65 {
66 std::string errStr = "CoordMatcher \"" + std::string(1,labels_[index_]) + "-coord\" takes at most two parameters <tol, relative>";
67 TEUCHOS_TEST_FOR_EXCEPTION(params.size()>2,std::logic_error,errStr);
68
69 // read in string, get double
70 if(params.size()>0) {
71 std::stringstream ss;
72 ss << params[0];
73 ss >> error_;
74 if(params.size()==2){
75 std::string errStr2 = params[1] + " is not a valid periodic option (try \"relative\")";
76 TEUCHOS_TEST_FOR_EXCEPTION(params[1]!="relative",std::logic_error,errStr2);
77 relative_ = true;
78 }
79 }
80 // else use default value for error
81 }
82
83public:
85 CoordMatcher(int index) : error_(1e-8),index_(index),relative_(false) { buildLabels(); }
86 CoordMatcher(int index,double error) : error_(error),index_(index),relative_(false) { buildLabels(); }
87 CoordMatcher(int index,const std::vector<std::string> & params) : error_(1e-8),index_(index),relative_(false)
88 { buildLabels(); parseParams(params); }
89
91
92 bool operator()(const Teuchos::Tuple<double,3> & a,
93 const Teuchos::Tuple<double,3> & b) const
94 {
95 double error = error_;
96 if(relative_) // scale error by length of domain
97 error*=std::fabs(a[1-index_]-b[1-index_]);
98 return std::fabs(a[index_]-b[index_])<error; /* I'm being lazy here! */
99 }
100
101 std::string getString() const
102 {
103 std::stringstream ss;
104 ss << labels_[index_] << "-coord <tol=" << error_ << ">";
105 return ss.str();
106 }
107
108 int getIndex() const {return index_;}
109
111 {
112 // This assumes a 2D x-y mesh even though the object supports the
113 // z direction. Once in 3D you need the PlaneMatcher.
114 TEUCHOS_ASSERT(index_ != 2);
115 if (index_ == 0)
116 return 1;
117 return 0;
118 }
119
120 double getAbsoluteTolerance() const {return error_;}
121
122 void transform(double * ptB, const std::vector<double> & centroidA) const
123 {
124 // Instead of matching directly, shift pt B given the centroid
125 // of side A
126 // For now, we assume at 2D x-y mesh as above so just need
127 // to overwrite the coordinate in the periodic direction
128
129 const int periodicIndex = this->getPeriodicDirection();
130 ptB[periodicIndex] = centroidA[periodicIndex];
131
132 return;
133 }
134};
135
139 double error_;
141 bool relative_; // compute error relative to length of domain
142 char labels_[3];
143
145 { labels_[0] = 'x'; labels_[1] = 'y'; labels_[2] = 'z'; }
146
147 void parseParams(const std::vector<std::string> & params)
148 {
149 std::string errStr = "PlaneMatcher \"" + std::string(1,labels_[index0_])+std::string(1,labels_[index1_])
150 + "-coord\" takes at most two parameter <tol, relative>";
151 TEUCHOS_TEST_FOR_EXCEPTION(params.size()>2,std::logic_error,errStr);
152
153 // read in string, get double
154 if(params.size()>0) {
155 std::stringstream ss;
156 ss << params[0];
157 ss >> error_;
158 if(params.size()==2){
159 if (params[1] == "3D") {
160 // Warn user but continue
161 std::cout << "WARNING : Keyword " << params[1] << " not needed for PlaneMatcher" << std::endl;
162 return;
163 }
164 std::string errStr2 = params[1] + " is not a valid periodic option (try \"relative\")";
165 TEUCHOS_TEST_FOR_EXCEPTION(params[1]!="relative",std::logic_error,errStr2);
166 relative_ = true;
167 }
168 }
169 // else use default value for error
170 return;
171 }
172
173public:
174 PlaneMatcher(int index0,int index1) : error_(1e-8),index0_(index0), index1_(index1), relative_(false)
175 { TEUCHOS_ASSERT(index0!=index1); buildLabels(); }
176
177 PlaneMatcher(int index0,int index1,double error) : error_(error),index0_(index0), index1_(index1), relative_(false)
178 { TEUCHOS_ASSERT(index0!=index1); buildLabels(); }
179
180 PlaneMatcher(int index0,int index1,const std::vector<std::string> & params)
181 : error_(1e-8), index0_(index0), index1_(index1), relative_(false)
182 { TEUCHOS_ASSERT(index0!=index1); buildLabels(); parseParams(params); }
183
185 { buildLabels(); }
186
187 bool operator()(const Teuchos::Tuple<double,3> & a,
188 const Teuchos::Tuple<double,3> & b) const
189 {
190 double error = error_;
191 if(relative_) // scale error by length of domain in normal direction
192 error*=std::fabs(a[3-index0_-index1_]-b[3-index0_-index1_]);
193 return (std::fabs(a[index0_]-b[index0_])<error_)
194 && (std::fabs(a[index1_]-b[index1_])<error_) ; /* I'm being lazy here! */
195 }
196
197 std::string getString() const
198 {
199 std::stringstream ss;
200 ss << labels_[index0_] << labels_[index1_] << "-coord <tol=" << error_ << ">";
201 return ss.str();
202 }
203
204 int getIndex0() const {return index0_;}
205 int getIndex1() const {return index1_;}
207 {
208 if (index0_ ==0) {
209 if (index1_ == 1)
210 return 2; // match x,y=periodic in z
211 else
212 return 1; // match x,z=periodic in y
213 }
214 else if (index0_ == 1) {
215 if (index1_ == 0)
216 return 2; // match y,x=periodic in z
217 else
218 return 0; // match y,z=periodic in x
219 }
220 else {
221 if (index1_ == 0)
222 return 1; // match z,x=periodic in y
223 else
224 return 0; // match z,y=periodic in x
225 }
226 }
227
228 double getAbsoluteTolerance() const {return error_;}
229
230 void transform(double * ptB, const std::vector<double> & centroidA) const
231 {
232 // Instead of matching directly, shift pt B given the centroid
233 // of side A
234 // For now, we assume the planes are aligned with one of the
235 // coordinate axes so we just need to overwrite the coordinate
236 // in the periodic direction
237
238 const int periodicIndex = this->getPeriodicDirection();
239 ptB[periodicIndex] = centroidA[periodicIndex];
240
241 return;
242 }
243};
244
248 double error_;
250 char labels_[3];
251
253 { labels_[0] = 'x'; labels_[1] = 'y'; labels_[2] = 'z'; }
254
255 void parseParams(const std::vector<std::string> & params)
256 {
257 std::string errStr = "QuarterPlaneMatcher \"(" + std::string(1,labels_[index0a_])+std::string(1,labels_[index0b_])+")"+std::string(1,labels_[index1_])
258 + "-quarter-coord\" takes only one parameter <tol>";
259 TEUCHOS_TEST_FOR_EXCEPTION(params.size()>1,std::logic_error,errStr);
260
261 // read in string, get double
262 if(params.size()==1) {
263 std::stringstream ss;
264 ss << params[0];
265 ss >> error_;
266 }
267 // else use default value for error
268 }
269
270public:
271 QuarterPlaneMatcher(int index0a,int index0b,int index1)
272 : error_(1e-8), index0a_(index0a), index0b_(index0b), index1_(index1)
273 { TEUCHOS_ASSERT(index0a!=index1); TEUCHOS_ASSERT(index0b!=index1); buildLabels(); }
274
275 QuarterPlaneMatcher(int index0a,int index0b,int index1,double error)
276 : error_(error), index0a_(index0a), index0b_(index0b), index1_(index1)
277 { TEUCHOS_ASSERT(index0a!=index1); TEUCHOS_ASSERT(index0b!=index1); buildLabels(); }
278
279 QuarterPlaneMatcher(int index0a,int index0b,int index1,const std::vector<std::string> & params)
280 : error_(1e-8), index0a_(index0a), index0b_(index0b), index1_(index1)
281 { TEUCHOS_ASSERT(index0a!=index1); TEUCHOS_ASSERT(index0b!=index1); buildLabels(); parseParams(params); }
282
285 { buildLabels(); }
286
287 bool operator()(const Teuchos::Tuple<double,3> & a,
288 const Teuchos::Tuple<double,3> & b) const
289 { return (std::fabs(a[index0a_]-b[index0b_])<error_)
290 && (std::fabs(a[index1_]-b[index1_])<error_) ; /* I'm being lazy here! */ }
291
292 std::string getString() const
293 {
294 std::stringstream ss;
295 ss << "(" << labels_[index0a_] << labels_[index0b_] << ")" << labels_[index1_] << "-quarter-coord <tol=" << error_ << ">";
296 return ss.str();
297 }
298
299 double getAbsoluteTolerance() const {return error_;}
300
301 void transform(double * ptB, const std::vector<double> & centroidA) const
302 {
303 // Instead of matching directly, shift pt B given the centroid
304 // of side A
305 // For now, we assume the planes are aligned with one of the
306 // coordinate axes
307 // We leave ptB[index1_] alone,
308 // put ptB[index0b_] in the index0a_ slot and replace it with
309 // centroidA[index0b_] which is the fixed value for plane B
310
311 ptB[index0a_] = ptB[index0b_];
312 ptB[index0b_] = centroidA[index0b_];
313
314 return;
315 }
316};
317
326 double error_;
331public:
332 enum class MirrorPlane : int {
333 XZ_PLANE=0,
334 YZ_PLANE=1
335 };
336 WedgeMatcher(MirrorPlane mp,const std::vector<std::string> & params )
337 : error_(1e-8),index0_(0),is_three_d_(true)
338 {
339 if (mp == MirrorPlane::XZ_PLANE)
340 index0_ = 1;
341 else // YZ_PLANE
342 index0_ = 0;
343
344 TEUCHOS_TEST_FOR_EXCEPTION(params.size() > 2,std::logic_error,"WedgeMatcher can only have one or two option parameters (tolerance and dimension)!");
345
346 // read in string, get double
347 if (params.size() > 0)
348 error_ = std::stod(params[0]);
349
350 if (params.size() > 1) {
351 if (params[1] == "2D")
352 is_three_d_ = false;
353 else if (params[1] == "3D")
354 is_three_d_ = true;
355 else {
356 TEUCHOS_TEST_FOR_EXCEPTION(true,std::runtime_error,"ERROR: WedgeMatcher::parsParams() - the second params must be iether \"2D\" or \"3D\", param=" << params[1] << "\n");
357 }
358 }
359 }
360 WedgeMatcher(const WedgeMatcher & cm) = default;
361
362 bool operator()(const Teuchos::Tuple<double,3> & a,
363 const Teuchos::Tuple<double,3> & b) const
364 {
365 if (is_three_d_) {
366 return ( (std::fabs(a[index0_]+b[index0_])<error_) &&
367 (std::fabs(a[1-index0_]-b[1-index0_])<error_) &&
368 (std::fabs(a[2]-b[2])<error_) );
369 }
370
371 // else 2D
372 return ( (std::fabs(a[index0_]+b[index0_])<error_) &&
373 (std::fabs(a[1-index0_]-b[1-index0_])<error_) );
374 }
375
376 std::string getString() const
377 {
378 std::stringstream ss;
379 if (index0_ == 0)
380 ss << "wy-coord <tol=" << error_ << ">";
381 else
382 ss << "wx-coord <tol=" << error_ << ">";
383 return ss.str();
384 }
385
386 int getIndex() const {return index0_;}
387
389 {
390 if (index0_ == 0)
393 }
394
395 bool isThreeD() const {return is_three_d_;}
396
397 double getAbsoluteTolerance() const {return error_;}
398
399 void transform(double * ptB, const std::vector<double> & centroidA) const
400 {
401 // Instead of matching directly, shift pt B given the centroid
402 // of side A
403 // For now, we assume the wedge is mirrored over the yz or xz plane
404 // Then we just need to mirror over the plane
405
406 ptB[index0_] = -ptB[index0_];
407
408 return;
409 }
410};
411
412} // end panzer_stk
413
414#endif
void parseParams(const std::vector< std::string > &params)
bool operator()(const Teuchos::Tuple< double, 3 > &a, const Teuchos::Tuple< double, 3 > &b) const
CoordMatcher(int index, const std::vector< std::string > &params)
void transform(double *ptB, const std::vector< double > &centroidA) const
CoordMatcher(int index)
index is the coordinate direction that will be compared to find matching nodes.
bool operator()(const Teuchos::Tuple< double, 3 > &a, const Teuchos::Tuple< double, 3 > &b) const
PlaneMatcher(int index0, int index1, double error)
PlaneMatcher(int index0, int index1, const std::vector< std::string > &params)
void transform(double *ptB, const std::vector< double > &centroidA) const
void parseParams(const std::vector< std::string > &params)
QuarterPlaneMatcher(int index0a, int index0b, int index1)
void transform(double *ptB, const std::vector< double > &centroidA) const
QuarterPlaneMatcher(int index0a, int index0b, int index1, const std::vector< std::string > &params)
QuarterPlaneMatcher(int index0a, int index0b, int index1, double error)
void parseParams(const std::vector< std::string > &params)
bool operator()(const Teuchos::Tuple< double, 3 > &a, const Teuchos::Tuple< double, 3 > &b) const
int index0_
index to compare - 0 for wy (mirrored over yz), 1 for wx (mirrored over xz)
bool operator()(const Teuchos::Tuple< double, 3 > &a, const Teuchos::Tuple< double, 3 > &b) const
WedgeMatcher::MirrorPlane getMirrorPlane() const
bool is_three_d_
Set to true if a 3D problem, set to false if 2D.
WedgeMatcher(const WedgeMatcher &cm)=default
void transform(double *ptB, const std::vector< double > &centroidA) const
WedgeMatcher(MirrorPlane mp, const std::vector< std::string > &params)