-
Notifications
You must be signed in to change notification settings - Fork 1
/
lapfparam_adap.m
121 lines (85 loc) · 3.34 KB
/
lapfparam_adap.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
function val = lapfparam_adap(kernel,xy,coord,norml,rho,legecoeff,lpanel,...
nlege,center_coord,funcurv,funcurv_par,nana_w,naeps,adapeps)
%LAPFPARAM_ADAP Adaptively evaluate potential generated by charges/dipoles on the boundary.
%
% Only adaptively evaluate near boundary points.
%
% Input parameters:
%
% kernel - The kernel for evaluating the potential.
% xy - An n by 2 matrix that represents the coordinates of points to
% evaluate.
% coord - A (npanel*nlege) by 2 matrix that represents the coordinates of
% Gaussian nodes on the boundary.
% norml - A (npanel*nlege) by 2 matrix that represents the normal of
% Gaussian nodes to the boundary.
% rho - A (npanel*lege) by 1 array that represents the unweighted density
% on the boundary.
% lpanel - The arclength of each leaf panel.
% nlege - The number of Legendre nodes on each leaf panel.
% center_coord - An npanel by 2 matrix that represents the coordinates of
% centers of circles on each leaf panel.
% funcurv - A function that parametrizes the boundary of the domain.
% funcurv_par - Parameters used by funcurv.
% nana_w - a large array containing, among other things, nested Legendre
% expansions of the arc length of the curve as a function of the
% curve parameter. Can be used as an input to the function
% nana_fparam (see).
% naeps - The absolute error tolerance used by nanafast.
% adapeps - The absolute error tolerance used by adapgaus.
%
% Output parameters:
%
% val - An n by 1 vector that contains potential.
%... init
npanel=size(center_coord,1);
radius=lpanel*ones(npanel,1);
val=zeros(size(xy,1),1);
bdpar{2}=kernel;
bdpar{3}=nana_w;
bdpar{4}=naeps;
bdpar{5}=funcurv;
for j=1:size(xy,1)
cur_coord=xy(j,:);
bdpar{7}=cur_coord;
in_arr=is_in_circle(cur_coord,center_coord,radius);
for k=1:length(in_arr)
irange=(k-1)*nlege+1:k*nlege;
if in_arr(k)
%... process nearby point by adapgaus
bdpar{1}=legecoeff(irange);
lrange=[(k-1)*lpanel k*lpanel];
bdpar{6}=lrange;
val(j)=val(j)+...
adapgaus(lrange(1),lrange(2),@density_expansion,...
bdpar,funcurv_par,12,adapeps);
else
%... regular point. processed by quadrature
val(j)=val(j)+...
dot(kernel(cur_coord,coord(irange,:),norml(irange,:)),rho(irange));
end
end
end
end
function val=density_expansion(x,bdpar,funcurv_par)
cur_legecoeff=bdpar{1};
kernel=bdpar{2};
nana_w=bdpar{3};
naeps=bdpar{4};
funcurv=bdpar{5};
lrange=bdpar{6};
cur_coord=bdpar{7};
xx=-1+2*(x-lrange(1))/(lrange(2)-lrange(1));
legeval=legeexev(xx,cur_legecoeff);
tout=nana_fparam(nana_w,naeps,x,'silent');
[coord_x,coord_y,dxdt,dydt]=funcurv(tout,funcurv_par);
coord=[coord_x coord_y];
norml=[-dydt dxdt];
val=kernel(cur_coord,coord,norml)*legeval;
end
function in = is_in_circle(pt,center_coord,radius)
n=size(center_coord,1);
assert(n==size(radius,1),'size(car_center,1)~=size(radius,1)');
d=sqrt(sum((pt-center_coord).^2,2));
in=d<radius-(1e-15);
end