1/*
2 * Laplacian.hpp
3 *
4 * Created on: Oct 5, 2015
5 * Author: i-bird
6 */
7
8#ifndef OPENFPM_NUMERICS_SRC_FINITEDIFFERENCE_LAPLACIAN_HPP_
9#define OPENFPM_NUMERICS_SRC_FINITEDIFFERENCE_LAPLACIAN_HPP_
10
11#include "FD_util_include.hpp"
12#include "util/util_num.hpp"
13#include "FiniteDifference/eq.hpp"
14
15/*! \brief Laplacian second order on h (spacing)
16 *
17 * \tparam Field which field derive
18 * \tparam impl which implementation
19 *
20 */
21template<typename arg, typename Sys_eqs, unsigned int impl=CENTRAL>
22class Lap
23{
24 /*! \brief Calculate which colums of the Matrix are non zero
25 *
26 * stub_or_real it is just for change the argument type when testing, in normal
27 * conditions it is a distributed map
28 *
29 * \param g_map map grid
30 * \param kmap position in the grid
31 * \param spacing of the grid
32 * \param gs Grid info
33 * \param cols non-zero colums calculated by the function
34 * \param coeff coefficent (constant in front of the derivative)
35 *
36 * ### Example
37 *
38 * \snippet FDScheme_unit_tests.hpp Laplacian usage
39 *
40 */
41 inline static void value(const typename stub_or_real<Sys_eqs,Sys_eqs::dims,typename Sys_eqs::stype,typename Sys_eqs::b_grid::decomposition::extended_type>::type & g_map,
42 grid_dist_key_dx<Sys_eqs::dims> & kmap ,
43 const grid_sm<Sys_eqs::dims,void> & gs,
44 std::unordered_map<long int,typename Sys_eqs::stype > & cols,
45 typename Sys_eqs::stype coeff)
46 {
47 std::cerr << "Error " << __FILE__ << ":" << __LINE__ << " only CENTRAL, FORWARD, BACKWARD derivative are defined";
48 }
49};
50
51/*! \brief Laplacian second order approximation CENTRAL Scheme
52 *
53 * \verbatim
54
55 1.0
56 *
57 | -4.0
58 1.0 *---+---* 1.0
59 |
60 *
61 1.0
62
63 * \endverbatim
64 *
65 *
66 */
67template<typename arg, typename Sys_eqs>
68class Lap<arg,Sys_eqs,CENTRAL>
69{
70public:
71
72
73 /*! \brief Calculate which colums of the Matrix are non zero
74 *
75 * stub_or_real it is just for change the argument type when testing, in normal
76 * conditions it is a distributed map
77 *
78 * \param g_map map grid
79 * \param kmap position in the grid
80 * \param spacing of the grid
81 * \param gs Grid info
82 * \param cols non-zero colums calculated by the function
83 * \param coeff coefficent (constant in front of the derivative)
84 *
85 * ### Example
86 *
87 * \snippet FDScheme_unit_tests.hpp Laplacian usage
88 *
89 */
90 inline static void value(const typename stub_or_real<Sys_eqs,Sys_eqs::dims,typename Sys_eqs::stype,typename Sys_eqs::b_grid::decomposition::extended_type>::type & g_map,
91 grid_dist_key_dx<Sys_eqs::dims> & kmap ,
92 const grid_sm<Sys_eqs::dims,void> & gs,
93 typename Sys_eqs::stype (& spacing )[Sys_eqs::dims] ,
94 std::unordered_map<long int,typename Sys_eqs::stype > & cols,
95 typename Sys_eqs::stype coeff)
96 {
97 // for each dimension
98 for (size_t i = 0 ; i < Sys_eqs::dims ; i++)
99 {
100 long int old_val = kmap.getKeyRef().get(i);
101 kmap.getKeyRef().set_d(i, kmap.getKeyRef().get(i) + 1);
102 arg::value(g_map,kmap,gs,spacing,cols,coeff/spacing[i]/spacing[i]);
103 kmap.getKeyRef().set_d(i,old_val);
104
105 old_val = kmap.getKeyRef().get(i);
106 kmap.getKeyRef().set_d(i, kmap.getKeyRef().get(i) - 1);
107 arg::value(g_map,kmap,gs,spacing,cols,coeff/spacing[i]/spacing[i]);
108 kmap.getKeyRef().set_d(i,old_val);
109
110 arg::value(g_map,kmap,gs,spacing,cols, - 2.0 * coeff/spacing[i]/spacing[i]);
111 }
112 }
113};
114
115
116/*! \brief Laplacian second order approximation CENTRAL Scheme (with central derivative in the single)
117 *
118 * \verbatim
119
120 1.0
121 *
122 | -4.0
123 1.0 *---+---* 1.0
124 |
125 *
126 1.0
127
128 * \endverbatim
129 *
130 *
131 */
132template<typename arg, typename Sys_eqs>
133class Lap<arg,Sys_eqs,CENTRAL_SYM>
134{
135public:
136
137
138 /*! \brief Calculate which colums of the Matrix are non zero
139 *
140 * stub_or_real it is just for change the argument type when testing, in normal
141 * conditions it is a distributed map
142 *
143 * \param g_map map grid
144 * \param kmap position in the grid
145 * \param spacing of the grid
146 * \param gs Grid info
147 * \param cols non-zero colums calculated by the function
148 * \param coeff coefficent (constant in front of the derivative)
149 *
150 * ### Example
151 *
152 * \snippet FDScheme_unit_tests.hpp Laplacian usage
153 *
154 */
155 inline static void value(const typename stub_or_real<Sys_eqs,Sys_eqs::dims,typename Sys_eqs::stype,typename Sys_eqs::b_grid::decomposition::extended_type>::type & g_map, grid_dist_key_dx<Sys_eqs::dims> & kmap , const grid_sm<Sys_eqs::dims,void> & gs, typename Sys_eqs::stype (& spacing )[Sys_eqs::dims] , std::unordered_map<long int,typename Sys_eqs::stype > & cols, typename Sys_eqs::stype coeff)
156 {
157 // for each dimension
158 for (size_t i = 0 ; i < Sys_eqs::dims ; i++)
159 {
160 long int old_val = kmap.getKeyRef().get(i);
161 kmap.getKeyRef().set_d(i, kmap.getKeyRef().get(i) + 2);
162 arg::value(g_map,kmap,gs,spacing,cols,coeff/spacing[i]/spacing[i]/4.0);
163 kmap.getKeyRef().set_d(i,old_val);
164
165 old_val = kmap.getKeyRef().get(i);
166 kmap.getKeyRef().set_d(i, kmap.getKeyRef().get(i) - 2);
167 arg::value(g_map,kmap,gs,spacing,cols,coeff/spacing[i]/spacing[i]/4.0);
168 kmap.getKeyRef().set_d(i,old_val);
169
170 arg::value(g_map,kmap,gs,spacing,cols, - 2.0 * coeff/spacing[i]/spacing[i]/4.0);
171 }
172 }
173};
174
175#endif /* OPENFPM_NUMERICS_SRC_FINITEDIFFERENCE_LAPLACIAN_HPP_ */
176