1 /*
2 XLiFE++ is an extended library of finite elements written in C++
3 Copyright (C) 2014 Lunéville, Eric; Kielbasiewicz, Nicolas; Lafranche, Yvon; Nguyen, Manh-Ha; Chambeyron, Colin
4
5 This program is free software: you can redistribute it and/or modify
6 it under the terms of the GNU General Public License as published by
7 the Free Software Foundation, either version 3 of the License, or
8 (at your option) any later version.
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
13 You should have received a copy of the GNU General Public License
14 along with this program. If not, see <http://www.gnu.org/licenses/>.
15 */
16 /*!
17 \file Laplace3dKernel.cpp
18 \author E. Lunéville
19 \since 17 jul 2014
20 \date 17 jul 2014
21
22 \brief Implementation of Laplace 3D kernels
23
24 Laplace3D kernel K(k; x, y)=(1/4 pi*|x-y|)
25 */
26
27 #include "Laplace3dKernel.hpp"
28 #include "utils.h"
29
30 namespace xlifepp
31 {
32
33
34 //--------------------------------------------------------------------------------------------
35 //construct Laplace3d Kernel : G(k; x, y)=(1/4 pi*|x-y|)
36 //--------------------------------------------------------------------------------------------
37
38 //construct a Laplace3d Kernel
Laplace3dKernel(Parameters & pars)39 Kernel Laplace3dKernel(Parameters& pars)
40 {
41 Kernel K;
42 K.name="Laplace 3D kernel";
43 K.shortname="Lap3D";
44 K.singularType =_r;
45 K.singularOrder = -1;
46 K.singularCoefficient = over4pi_;
47 K.symmetry=_symmetric;
48 K.userData.push(pars);
49 K.kernel = Function(Laplace3d, K.userData);
50 K.gradx = Function(Laplace3dGradx, K.userData);
51 K.grady = Function(Laplace3dGrady, K.userData);
52 K.ndotgradx = Function(Laplace3dNxdotGradx, K.userData);
53 K.ndotgrady = Function(Laplace3dNydotGrady, K.userData);
54 return K;
55 }
56
Laplace3d(const Point & x,const Point & y,Parameters & pars)57 real_t Laplace3d(const Point& x, const Point& y, Parameters& pars)
58 {
59 real_t r = x.distance(y);
60 return over4pi_ / r;
61 }
62
63 //derivatives
Laplace3dGradx(const Point & x,const Point & y,Parameters & pars)64 Vector<real_t> Laplace3dGradx(const Point& x, const Point& y,Parameters& pars)
65 {
66 real_t r2 = x.squareDistance(y);
67 real_t r = std::sqrt(r2);
68 Vector<real_t> g1(3);
69 scaledVectorTpl(-over4pi_ / (r*r2), x.begin(), x.end(), y.begin(), g1.begin());
70 return g1;
71 }
72
Laplace3dGrady(const Point & x,const Point & y,Parameters & pars)73 Vector<real_t> Laplace3dGrady(const Point& x, const Point& y,Parameters& pars)
74 {
75 real_t r2 = x.squareDistance(y);
76 real_t r = std::sqrt(r2);
77 Vector<real_t> g1(3);
78 scaledVectorTpl(over4pi_ / (r*r2), x.begin(), x.end(), y.begin(), g1.begin());
79 return g1;
80 }
81
Laplace3dNxdotGradx(const Point & x,const Point & y,Parameters & pars)82 real_t Laplace3dNxdotGradx(const Point& x, const Point& y, Parameters& pars)
83 {
84 Vector<real_t>& nxp = getNx();
85 std::vector<real_t>::const_iterator itx=x.begin(),ity=y.begin(), itn=nxp.begin();
86 real_t d1=(*itx++ - *ity++), d2=(*itx++ - *ity++), d3=(*itx - *ity);
87 real_t r = d1* *itn++; r+= d2* *itn++; r+=d3* *itn;
88 real_t r3= d1*d1+d2*d2+d3*d3; r3*=std::sqrt(r3);
89 return -r*over4pi_/r3;
90 }
91
Laplace3dNydotGrady(const Point & x,const Point & y,Parameters & pars)92 real_t Laplace3dNydotGrady(const Point& x, const Point& y, Parameters& pars)
93 {
94 Vector<real_t>& nyp = getNy();
95 std::vector<real_t>::const_iterator itx=x.begin(),ity=y.begin(), itn=nyp.begin();
96 real_t d1=(*itx++ - *ity++), d2=(*itx++ - *ity++), d3=(*itx - *ity);
97 real_t r = d1* *itn++; r+= d2* *itn++; r+=d3* *itn;
98 real_t r3= d1*d1+d2*d2+d3*d3; r3*=std::sqrt(r3);
99 return r*over4pi_/r3;
100 }
101
102 }
103