OpenFPM_pdata  4.1.0
Project that contain the implementation of distributed structures
 
Loading...
Searching...
No Matches
Vandermonde.hpp
1//
2// Created by tommaso on 21/03/19.
3// Edited by Abhinav Singh on 24/01/2022
4
5#ifndef OPENFPM_PDATA_VANDERMONDE_HPP
6#define OPENFPM_PDATA_VANDERMONDE_HPP
7
8#include "MonomialBasis.hpp"
9#include "VandermondeRowBuilder.hpp"
10#include "Support.hpp"
11
12template<unsigned int dim, typename T, typename MatrixType>
14{
15private:
16 const Point<dim, T> point;
18 const MonomialBasis<dim> monomialBasis;
19 T eps,HOverEpsilon,minSpacing;
20
21public:
22/* Vandermonde(const Point<dim, T> &point, const std::vector<Point<dim, T>> &neighbours,
23 const MonomialBasis<dim> &monomialBasis);*/
24
25 template<typename vector_type,
26 typename vector_type2>
27 Vandermonde(const Support &support,
28 const MonomialBasis<dim> &monomialBasis,
29 const vector_type & particlesFrom,
30 const vector_type2 & particlesTo,T HOverEpsilon=0.5) //0.5 for the test
31 : point(particlesTo.getPosOrig(support.getReferencePointKey())),
32 monomialBasis(monomialBasis),HOverEpsilon(HOverEpsilon)
33 {
34 initialize(support,particlesFrom,particlesTo);
35 }
36
37
38 MatrixType &getMatrix(MatrixType &M)
39 {
40 // Build the Vandermonde matrix, row-by-row
41 VandermondeRowBuilder<dim, T> vrb(monomialBasis);
42 unsigned int row = 0;
43
44 size_t N = offsets.size();
45 for (size_t i = 0; i < N; ++i)
46 {
47 const auto& offset = offsets.get(i);
48 vrb.buildRow(M, row, offset, eps);
49 ++row;
50 }
51 return M;
52 }
53
54 T getEps()
55 {
56 return eps;
57 }
58 T getMinSpacing()
59 {
60 return minSpacing;
61 }
62
63private:
64
65
66 void computeEps(T factor)
67 {
68 T avgNeighbourSpacing = 0;
69 minSpacing=std::numeric_limits<T>::max();
70 size_t N = offsets.size();
71 for (size_t i = 0; i < N; ++i)
72 {
73 const auto& offset = offsets.get(i);
74 double dist=norm(offset);
75 avgNeighbourSpacing += computeAbsSum(offset);
76 if(minSpacing>dist)
77 {
78 minSpacing=dist;
79 }
80 }
81 avgNeighbourSpacing /= offsets.size();
82 eps = avgNeighbourSpacing/factor;
83 assert(eps != 0);
84 }
85
86 static T computeAbsSum(const Point<dim, T> &x)
87 {
88 T absSum = 0;
89 for (unsigned int i = 0; i < dim; ++i)
90 {
91 absSum += fabs(x.value(i));
92 }
93 return absSum;
94 }
95
96 template<typename vector_type, typename vector_type2>
97 void initialize(const Support &sup, const vector_type & particlesFrom, vector_type2 &particlesTo)
98 {
99 auto & keys = sup.getKeys();
100
101 for (int i = 0 ; i < keys.size() ; i++)
102 {
103 Point<dim,T> p = particlesTo.getPosOrig(sup.getReferencePointKey());
104 p -= particlesFrom.getPosOrig(keys.get(i));
105 offsets.add(p);
106 }
107
108 // First check that the number of points given is enough for building the Vandermonde matrix
109 if (offsets.size() < monomialBasis.size())
110 {
111 ACTION_ON_ERROR(std::length_error("Not enough neighbour points passed for Vandermonde matrix construction!"));
112 }
113 // Compute eps for this point
114 //factor here. This is C factor.
115 computeEps(HOverEpsilon);
116 }
117
118};
119
120
121
122#endif //OPENFPM_PDATA_VANDERMONDE_HPP
This class implement the point shape in an N-dimensional space.
Definition Point.hpp:28
Implementation of 1-D std::vector like structure.
size_t size()
Stub size.
auto getPosOrig(vect_dist_key_dx vec_key) -> decltype(vd.getPos(vec_key))
Get the position of an element.
Distributed vector.
auto getPosOrig(vect_dist_key_dx vec_key) const -> decltype(v_pos.template get< 0 >(vec_key.getKey()))
Get the position of an element.