build_FEbuf.m
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This file is part of NS2DDV. %
% %
% Copyright(C) 2011-2018 C. Calgaro (caterina.calgaro@math.univ-lille1.fr) %
% E. Creusé (emmanuel.creuse@math.univ-lille1.fr) %
% T. Goudon (thierry.goudon@inria.fr) %
% A. Mouton (alexandre.mouton@math.univ-lille1.fr) %
% %
% NS2DDV is free software: you can redistribute it and/or modify it under the terms %
% of the GNU General Public License as published by the Free Software Foundation, %
% either version 3 of the License, or (at your option) any later version. %
% %
% NS2DDV is distributed in the hope that it will be useful, but WITHOUT ANY %
% WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A %
% PARTICULAR PURPOSE. See the GNU General Public License for more details. %
% %
% You should have received a copy of the GNU General Public License along with %
% NS2DDV. If not, see <http://www.gnu.org/licenses/>. %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function[usrbuf] = build_FEbuf(answers)
% This routine must be used after
% - build_PHYSICALbuf
% - build_DOMAINbuf
usrbuf = {};
l = 1;
switch answers.MODEL
case {'NS', 'NSDV'}
usrbuf{l} = '';
l = l+1;
usrbuf{l} = '% Finite Element parameters';
l = l+1;
usrbuf{l} = '%%%%%%%%%%%%%%%%%%%%%%%%%%%';
l = l+1;
% Time semi-discretization
%%%%%%%%%%%%%%%%%%%%%%%%%%
usrbuf{l} = '% Time semi-discretization of the Stokes equation';
l = l+1;
usrbuf{l} = '% ''BDF2'' Implicit 2-steps Backward Differentiation Formula';
l = l+1;
usrbuf{l} = '% ''BDF2_PROJ'' Implicit 2-steps Backward Differentiation Formula + Projection method';
l = l+1;
if strcmp(answers.MODEL, 'NSDV')
usrbuf{l} = '% ''BDF2_PROJ_VAR'' Implicit 2-steps Backward Differentiation Formula + Projection method with variable density';
l = l+1;
if (max(strcmp(answers.TESTCASE, answers.ANALYTIC_TESTCASES)) == 1)
usrbuf{l} = '% ''BDF2_PROJ_VAR_A'' Implicit 2-steps Backward Differentiation Formula + Projection method with variable analytical density';
l = l+1;
end
end
if max(strcmp(answers.TESTCASE, {'POIS', 'NONA', 'GTPSI'}))
usrbuf{l} = '% Warning: ''STRESS_BALANCE'' and ''NATURAL'' boundary conditions are compatible only with ''BDF2''';
l = l+1;
end
usrbuf{l} = 'PARAMETERS.FE.SCHEME = ''BDF2'';';
l = l+1;
usrbuf{l} = '% Amplitude of the rotational term in the pressure correction for ''BDF2_PROJ*'' methods';
l = l+1;
if strcmp(answers.CV_STUDY, 'IN_REYNOLDS')
usrbuf{l} = 'PARAMETERS.FE.ROT_CORRECTION = num2cell(0.5./cell2mat(PARAMETERS.PHYSICAL.RE));';
l = l+1;
else
usrbuf{l} = 'PARAMETERS.FE.ROT_CORRECTION = 0.5/PARAMETERS.PHYSICAL.RE;';
l = l+1;
end
% Finite element discretization in space
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
usrbuf{l} = '% Space discretization of velocity';
l = l+1;
usrbuf{l} = '% ''P1B'' P1-bubble finite elements approximation';
l = l+1;
usrbuf{l} = '% ''P2'' P2 finite elements approximation';
l = l+1;
usrbuf{l} = 'PARAMETERS.FE.TYPE = ''P2'';';
l = l+1;
usrbuf{l} = '% Space discretization of density in Stokes equation';
l = l+1;
usrbuf{l} = '% ''P1'' P1 finite elements approximation';
l = l+1;
usrbuf{l} = 'PARAMETERS.FE.TYPE_DENSITY = ''P1'';';
l = l+1;
% Closure constrain
%%%%%%%%%%%%%%%%%%%
usrbuf{l} = '% Constrain for closing the pressure problem';
l = l+1;
usrbuf{l} = '% ''ZERO_FIXED_POINT'' p = 0 is imposed on a specific point in the space domain';
l = l+1;
usrbuf{l} = '% ''ZERO_AVERAGE'' the average of p is set to 0 within the discrete system';
l = l+1;
usrbuf{l} = '% ''ZERO_AVERAGE_APOSTERIORI'' the average of p is set to 0 as an a posteriori constrain (only works with PARAMETERS.FE.SCHEME = ''BDF2'')';
l = l+1;
usrbuf{l} = 'PARAMETERS.FE.PRESSURE_CONSTRAIN = ''ZERO_FIXED_POINT'';';
l = l+1;
switch answers.TESTCASE
case {'CAEN'}
switch answers.GEOM
case {'RECTANGLE'}
xcoord = ' 0.5*(PARAMETERS.DOMAIN.XMIN+PARAMETERS.DOMAIN.XMAX)';
ycoord = ' 0.5*(PARAMETERS.DOMAIN.YMIN+PARAMETERS.DOMAIN.YMAX)';
otherwise
error('''CAEN'' test case is only compatible with ''RECTANGLE'' domain geometry');
end
case {'EXAC'}
switch answers.GEOM
case {'RECTANGLE'}
xcoord = ' 0.5*(PARAMETERS.DOMAIN.XMIN+PARAMETERS.DOMAIN.XMAX)';
ycoord = ' 0.5*(PARAMETERS.DOMAIN.YMIN+PARAMETERS.DOMAIN.YMAX)';
otherwise
error('''EXAC'' test case is only compatible with ''RECTANGLE'' domain geometry');
end
case {'RTIN'}
switch answers.GEOM
case {'RECTANGLE'}
xcoord = ' 0.5*(PARAMETERS.DOMAIN.XMIN+PARAMETERS.DOMAIN.XMAX)';
ycoord = ' 0.75*PARAMETERS.DOMAIN.YMAX+0.25*PARAMETERS.DOMAIN.YMIN';
otherwise
error('''RTIN'' test case is only compatible with ''RECTANGLE'' domain geometry');
end
case {'DROP'}
switch answers.GEOM
case {'RECTANGLE'}
xcoord = ' 0.5*(PARAMETERS.DOMAIN.XMIN+PARAMETERS.DOMAIN.XMAX)';
ycoord = ' 0.5*(PARAMETERS.DOMAIN.YMIN+PARAMETERS.PHYSICAL.DROP_INTERFACE_Y)';
otherwise
error('''DROP'' test case is only compatible with ''RECTANGLE'' domain geometry');
end
case {'NONA'}
switch answers.GEOM
case {'RECTANGLE', 'STEP', 'DIHEDRON'}
xcoord = ' 0.1*PARAMETERS.DOMAIN.XMIN+0.9*PARAMETERS.DOMAIN.XMAX';
ycoord = ' 0.1*PARAMETERS.DOMAIN.YMIN+0.9*PARAMETERS.DOMAIN.YMAX';
otherwise
error('''NONA'' test case is only compatible with ''RECTANGLE'', ''STEP'' and ''DIHEDRON'' domain geometries');
end
case {'POIS'}
switch answers.GEOM
case {'RECTANGLE', 'STEP'}
xcoord = ' 0.1*PARAMETERS.DOMAIN.XMIN+0.9*PARAMETERS.DOMAIN.XMAX';
ycoord = ' 0.1*PARAMETERS.DOMAIN.YMIN+0.9*PARAMETERS.DOMAIN.YMAX';
otherwise
error('''POIS'' test case is only compatible with ''RECTANGLE'' and ''STEP'' domain geometries');
end
case {'EXACNEU'}
switch answers.GEOM
case {'RECTANGLE'}
xcoord = ' 0.5*(PARAMETERS.DOMAIN.XMIN+PARAMETERS.DOMAIN.XMAX)';
ycoord = ' 0.5*(PARAMETERS.DOMAIN.YMIN+PARAMETERS.DOMAIN.YMAX)';
otherwise
error('''EXACNEU'' test case is only compatible with ''RECTANGLE'' domain geometry');
end
case {'GTPSI'}
switch answers.GEOM
case {'RECTANGLE'}
xcoord = ' 0.9*PARAMETERS.DOMAIN.XMIN+0.1*PARAMETERS.DOMAIN.XMAX';
ycoord = ' 0.5*(PARAMETERS.DOMAIN.YMIN+PARAMETERS.DOMAIN.YMAX)';
otherwise
error('''GTPSI'' test case is only compatible with ''RECTANGLE'' domain geometry');
end
otherwise
error('Unknown test case');
end
usrbuf{l} = '% Coordinates of the null pressure point';
l = l+1;
usrbuf{l} = '% WARNING: if the domain geometry parameters are obtained from an external file, it is not recommended to use';
l = l+1;
usrbuf{l} = '% the variables PARAMETERS.DOMAIN.XXX to define the coordinates of the null pressure point. Use numerical values instead';
l = l+1;
usrbuf{l} = strcat('PARAMETERS.FE.NFX_X = ', xcoord, ';');
l = l+1;
usrbuf{l} = strcat('PARAMETERS.FE.NFX_Y = ', ycoord, ';');
l = l+1;
otherwise
error('Unknown model');
end
end
