function [discret_data,...
bound_data] = bht_boundaries_init(gene_data,...
discret_data,...
bound_data)
% BHT_BOUNDARIES_INIT
%
% Build boundaries' meshes, compute the tangent and normal vectors,
% compute right hand side boundary terms of the NBVP.
%
% Syntax
%
% [discret_data,bound_data] =
% BHT_BOUNDARIES_INIT(gene_data,discret_data,bound_data)
%
% Description
%
% The function is ran by BHT_DATA_COMPILE. It calls the M-Files the
% boundaries are described in and computes the discretization points
% (using the parameter mesh_size of the DAT-File) as well as the
% tangent and normal vectors. Right hand side boundary terms arising in
% the NBVP (Neumann boundary value problem for the fluid potential) are
% also computed.
%
% Output variables
%
% The structure array variables 'discret_data' and 'bound_data' are
% preallocated in the M-File BHT_TRANSLATION. The following fields are
% set up by BHT_BOUNDARIES_INIt (we refer to the list of output
% variables for a description).
%
% o Field(s) of discret_data: last_line
%
% o Field(s) of bound_data: nbr_modes, nbr_points, t, dt, X, dl, Nn, N,
% K, flux, fflux
%
% See also BHT_DATA_COMPILE, BHT_TRANSLATION
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %
% BIOHYDRODYNAMICS MATLAB TOOLBOX %
% %
% A. MUNNIER and B. PINCON %
% %
% alexandre.munnier@iecn.u-nancy.fr bruno.pincon@iecn.u-nancy.fr %
% http://www.iecn.u-nancy.fr/~munnier http://www.iecn.u-nancy.fr/~pincon %
% %
% %
% INSTITUT ELIE CARTAN, NANCY 1 %
% http://www.iecn.u-nancy.fr/ %
% %
% INRIA Lorraine, Projet CORIDA %
% http://www.iecn.u-nancy.fr/~corida/ %
% %
% %
% %
% %
% August 15th 2008 %
% %
% GNU GPL v3 licence %
% %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%##########################################################################
% Preallocation of local variables
%---------------------------------
nb_bords = gene_data.nb_boundaries; % total number of boundaries
nb_solids = gene_data.nb_solids; % total number of solids
h = discret_data.mesh_size; % mesh size for boundaries'
% discretization
total_nb_points = 0; % Total number of points'
% preallocation
discret_data.last_line = []; % this variable is used for
% solving BOUNDED Neumann
% boundary value problems
% with Nystrom's method
%##########################################################################
% Main loop (on the boundaries)
%------------------------------
for k = 1:nb_bords;
param = bound_data(k).settings;
DL = @(t) sqrt(sum(feval(bound_data(k).mfilename,t,1,param).^2)); % compute the length
L = quad(DL,0,2*pi); % of boundary k
%======================================================================
n = floor(0.5*L/h); % number of modes
% (used for solving
% Neumann boundary
% value problem with
%++++++++++++++++++++++++++ % Nystrom's method)
bound_data(k).nbr_modes = n;
%++++++++++++++++++++++++++
%======================================================================
m = 2*n+1 ; % number of disretization
% points for
%+++++++++++++++++++++++++++ % boundary k
bound_data(k).nbr_points = m;
%+++++++++++++++++++++++++++
total_nb_points = total_nb_points + m;
%======================================================================
t = linspace(0,2*pi,m+1); t(m+1) = []; % points' locations
% on the boundary k
%++++++++++++++++++
bound_data(k).t = t;
%++++++++++++++++++
%======================================================================
bound_data(k).dt = 2*pi/m; % discretization size on
% boundary k
%======================================================================
bound_data(k).X = feval(bound_data(k).mfilename,t,0,param); % points' coordinates
%======================================================================
Xp = feval(bound_data(k).mfilename,t,1,param); % tangent vector
Xpp = feval(bound_data(k).mfilename,t,2,param); % second derivative
bound_data(k).dl = sqrt(sum((Xp).^2)); % area element
discret_data.last_line = [discret_data.last_line,...
bound_data(k).dt*bound_data(k).dl]; % used in Nystrom's method
bound_data(k).Nn = [-Xp(2,:);Xp(1,:)]./...
[bound_data(k).dl;bound_data(k).dl]; % Unitary normal vector
bound_data(k).N = [-Xp(2,:);Xp(1,:)]; % Normal vector (not unitary)
bound_data(k).K = sum(Xpp.*bound_data(k).Nn)./bound_data(k).dl.^2; % curvature
bound_data(k).flux = sum(bound_data(k).X.*Xp)./bound_data(k).dl; % flux
bound_data(k).fflux = sum(bound_data(k).X.*bound_data(k).Nn); % flux for the shape derivative
%======================================================================
end;
%##########################################################################
% Compute diameters of exterior boundaries and
% change order of exterior boundaries : from largest to smallest diameter
%------------------------------------------------------------------------
diam = Inf.*ones(1,nb_bords);
for k = nb_solids+1:nb_bords
m = bound_data(k).nbr_points;
dx1 = repmat(bound_data(k).X(1,:),m,1) - repmat((bound_data(k).X(1,:))',1,m);
dx2 = repmat(bound_data(k).X(2,:),m,1) - repmat((bound_data(k).X(2,:))',1,m);
diam(k) = max(max(dx1.^2 + dx2.^2));
end
[bb,XI] = sort(diam,'descend');
bound_data(XI) = bound_data;
%--------------------------------------------------------------------------
disp(['Total number of discretization poins: ',num2str(total_nb_points)]);
%==========================================================================
end
