Panzer  Version of the Day
Panzer_Normals_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_NORMALS_IMPL_HPP
44 #define PANZER_NORMALS_IMPL_HPP
45 
46 #include <algorithm>
49 #include "Intrepid2_FunctionSpaceTools.hpp"
50 #include "Intrepid2_CellTools.hpp"
51 
52 namespace panzer {
53 
54 //**********************************************************************
56  : normalize(true)
57 {
58  // Read from parameters
59  const std::string name = p.get<std::string>("Name");
60  side_id = p.get<int>("Side ID");
63  if(p.isParameter("Normalize")) // set default
64  normalize = p.get<bool>("Normalize");
65 
66  // grab information from quadrature rule
67  Teuchos::RCP<PHX::DataLayout> vector_dl = quadRule->dl_vector;
68  quad_order = quadRule->cubature_degree;
69 
70  // build field, set as evaluated type
71  normals = PHX::MDField<ScalarT,Cell,Point,Dim>(name, vector_dl);
72  this->addEvaluatedField(normals);
73 
74  std::string n = "Normals: " + name;
75  this->setName(n);
76 }
77 
78 //**********************************************************************
80 {
81  this->utils.setFieldData(normals,fm);
82 
83  num_qp = normals.dimension(1);
84  num_dim = normals.dimension(2);
85 
86  quad_index = panzer::getIntegrationRuleIndex(quad_order,(*sd.worksets_)[0], this->wda);
87 }
88 
89 //**********************************************************************
90 PHX_EVALUATE_FIELDS(Normals,workset)
91 {
92  // ECC Fix: Get Physical Side Normals
93 
94  if(workset.num_cells>0) {
95  Intrepid2::CellTools<ScalarT>::getPhysicalSideNormals(normals,
96  this->wda(workset).int_rules[quad_index]->jac,
97  side_id, *this->wda(workset).int_rules[quad_index]->int_rule->topology);
98 
99  if(normalize) {
100  // normalize vector: getPhysicalSideNormals does not
101  // return normalized vectors
102  for(index_t c=0;c<workset.num_cells;c++) {
103  for(std::size_t q=0;q<num_qp;q++) {
104  ScalarT norm = 0.0;
105 
106  // compute squared norm
107  for(std::size_t d=0;d<num_dim;d++)
108  norm += normals(c,q,d)*normals(c,q,d);
109 
110  // adjust for length of vector, now unit vectors
111  norm = sqrt(norm);
112  for(std::size_t d=0;d<num_dim;d++)
113  normals(c,q,d) /= norm;
114  }
115  }
116  }
117  // else normals correspond to differential
118  }
119 }
120 
121 //**********************************************************************
122 
123 }
124 
125 #endif
std::size_t num_qp
T * get() const
std::size_t quad_index
PHX::MDField< ScalarT, Cell, Point, Dim > normals
PHX_EVALUATOR_CTOR(BasisValues_Evaluator, p)
PHX_EVALUATE_FIELDS(BasisValues_Evaluator, workset)
std::vector< std::string >::size_type getIntegrationRuleIndex(int ir_degree, panzer::Workset &workset, WorksetDetailsAccessor &wda)
PHX_POST_REGISTRATION_SETUP(BasisValues_Evaluator, sd, fm)