KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
density_function_polynomial.h
Go to the documentation of this file.
1 /*
2 ==============================================================================
3 Kratos
4 A General Purpose Software for Multi-Physics Finite Element Analysis
5 Version 1.0 (Released on march 05, 2007).
6 
7 Copyright 2007
8 Pooyan Dadvand, Riccardo Rossi
9 pooyan@cimne.upc.edu
10 rrossi@cimne.upc.edu
11 CIMNE (International Center for Numerical Methods in Engineering),
12 Gran Capita' s/n, 08034 Barcelona, Spain
13 
14 Permission is hereby granted, free of charge, to any person obtaining
15 a copy of this software and associated documentation files (the
16 "Software"), to deal in the Software without restriction, including
17 without limitation the rights to use, copy, modify, merge, publish,
18 distribute, sublicense and/or sell copies of the Software, and to
19 permit persons to whom the Software is furnished to do so, subject to
20 the following condition:
21 
22 Distribution of this code for any commercial purpose is permissible
23 ONLY BY DIRECT ARRANGEMENT WITH THE COPYRIGHT OWNER.
24 
25 The above copyright notice and this permission notice shall be
26 included in all copies or substantial portions of the Software.
27 
28 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
29 EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
30 MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
31 IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
32 CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
33 TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
34 SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
35 
36 ==============================================================================
37  */
38 
39 //
40 // Project Name: Kratos
41 // Last Modified by: $Author: gcasas $
42 // Date: $Date: 2014-10-09 10:34:00 $
43 // Revision: $Revision: 0.1 $
44 //
45 //
46 
47 #if !defined(KRATOS_DENSITY_FUNCTION_POLYNOMIAL)
48 #define KRATOS_DENSITY_FUNCTION_POLYNOMIAL
49 // System includes
50 #include <string>
51 #include <iostream>
52 #include <cstdlib>
53 #include "density_function.h"
54 #include <cmath>
55 
56 namespace Kratos
57 {
58 
59 template< unsigned int TDim>
61 {
62 public:
64 
65 DensityFunctionPolynomial(const double range, const double shape_factor)
66  : mR(range),
67  mHeightOverR(shape_factor)
68 {
69  m6 = fm6(mR);
70  m4 = fm4(mR);
71  m2 = fm2(mR);
72  m0 = fm0(mR);
73 }
74 
76 
77 void ComputeWeights(std::vector<double> & distances, std::vector<double> & nodal_areas, std::vector<double> & weights)
78 {
79  double total_nodal_area_inv = 0.0;
80 
81  for (unsigned int i = 0; i != distances.size(); ++i){
82  total_nodal_area_inv += nodal_areas[i];
83  }
84 
85  total_nodal_area_inv = 1.0 / total_nodal_area_inv;
86 
87  double sum_of_weights_inv = 0.0;
88 
89  for (unsigned int i = 0; i != distances.size(); ++i){
90  double radius_2 = distances[i] * distances[i];
91  double weight = nodal_areas[i] * (m6 * std::pow(radius_2, 3) + m6 * m2 * radius_2 + m0);
92  weights[i] = weight;
93  sum_of_weights_inv += weight;
94  }
95 
96  sum_of_weights_inv = 1.0 / sum_of_weights_inv;
97 
98  // normalizing weights
99  for (unsigned int i = 0; i != distances.size(); ++i){
100  weights[i] *= sum_of_weights_inv;
101  }
102 }
103 
104 void ComputeWeights(std::vector<double> & distances, std::vector<double> & nodal_areas, const double max_nodal_area_inv, std::vector<double> & weights)
105 {
106  double total_nodal_area_inv = 0.0;
107 
108  for (unsigned int i = 0; i != distances.size(); ++i){
109  total_nodal_area_inv += nodal_areas[i];
110  }
111 
112  total_nodal_area_inv = 1.0 / total_nodal_area_inv;
113 
114  double sum_of_weights_inv = 0.0;
115 
116  for (unsigned int i = 0; i != distances.size(); ++i){
117  double radius_2 = distances[i] * distances[i];
118  double r = mR;// * pow(nodal_areas[i] * max_nodal_area_inv, 0.3333333333333333333333333333333);
119  double weight;
120 
121  weight = nodal_areas[i] * ((radius_2 > r * r) ? 0.0 : m6 * std::pow(radius_2, 3) + m6 * m2 * radius_2 + m0);
122  weights[i] = weight;
123  sum_of_weights_inv += weight;
124  }
125 
126  if(std::abs(sum_of_weights_inv) < std::numeric_limits<double>::epsilon()) sum_of_weights_inv = 0.0;
127  else sum_of_weights_inv = 1.0 / sum_of_weights_inv;
128 
129  // normalizing weights
130 
131  for (unsigned int i = 0; i != distances.size(); ++i){
132  weights[i] *= sum_of_weights_inv;
133  }
134 }
135 
136 private:
137 
138 double mR; // the radius of the support of the PDF;
139 double mHeightOverR;
140 double m6;
141 double m4;
142 double m2;
143 double m0;
144 
145 double fm6(double r)
146 {
147  return 7/(16*std::pow(mR, 7));
148 }
149 
150 double fm4(double r)
151 {
152  return -15*(std::pow(mR, 2) + 2)*std::exp(std::pow(mR, 2)/2)/(2*std::pow(mR, 4)*(2*std::pow(mR, 3)*std::exp(std::pow(mR, 2)/2) + 14*mR*std::exp(std::pow(mR, 2)/2) - 15*std::sqrt(2)*std::sqrt(Globals::Pi)*std::exp(std::pow(mR, 2))*erf(std::sqrt(2)*mR/2)));
153 }
154 
155 double fm2(double r)
156 {
157  return -3*std::pow(mR, 4);
158 }
159 
160 double fm0(double r)
161 {
162  return 7/(8*mR);
163 }
164 
165 }; // class DensityFunctionPolynomial
166 
167 } // namespace Kratos.
168 
169 #endif // KRATOS_DENSITY_FUNCTION_POLYNOMIAL
Definition: density_function.h:59
Definition: density_function_polynomial.h:61
DensityFunctionPolynomial(const double range, const double shape_factor)
Definition: density_function_polynomial.h:65
virtual ~DensityFunctionPolynomial()
Definition: density_function_polynomial.h:75
void ComputeWeights(std::vector< double > &distances, std::vector< double > &nodal_areas, const double max_nodal_area_inv, std::vector< double > &weights)
Definition: density_function_polynomial.h:104
void ComputeWeights(std::vector< double > &distances, std::vector< double > &nodal_areas, std::vector< double > &weights)
Definition: density_function_polynomial.h:77
KRATOS_CLASS_POINTER_DEFINITION(DensityFunctionPolynomial)
constexpr double Pi
Definition: global_variables.h:25
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
integer i
Definition: TensorModule.f:17