% Copyright (C) 2017-2017 Yves Renard.
%
% This file is a part of GetFEM++
%
% GetFEM++  is  free software;  you  can  redistribute  it  and/or modify it
% under  the  terms  of the  GNU  Lesser General Public License as published
% by  the  Free Software Foundation;  either version 3 of the License,  or
% (at your option) any later version along with the GCC Runtime Library
% Exception either version 3.1 or (at your option) any later version.
% This program  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 Lesser General Public
% License and GCC Runtime Library Exception for more details.
% You  should  have received a copy of the GNU Lesser General Public License
% along  with  this program;  if not, write to the Free Software Foundation,
% Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301, USA.
% -*- matlab -*- (enables emacs matlab mode)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% parameters for program bilaplacian                      %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%% pde parameters :	                              %%%%%
LX = 1.0; LY = LX; LZ = LX;	% sizes of the domain.
N = 2;
QUAD = 0;
MESH_NOISED = 0;        % Set to one if you want to "shake" the mesh
NX = 40;                % space step.
DIRICHLET_VERSION = 1;  % 0=multipliers, 1=penalization
EPS_DIRICHLET_PENAL = 1e10; % Parametre de penalisation
FT = 10.0;              % parameter for the exact solution
D = 1.0;    % "flexion modulus"
KL = 1;     % 0 for pure bilaplacian problem or 1 for Kirchhoff-Love
            % plate model (N=2 only).
NU = 0.3;   % Poisson ratio (0 <= NU <= 1), for Kirchhoff-Love plate model only

BOUNDARY_REF = 1;       % 0 :  non-homogeneous Dirichlet 
                        % 1 :  mixed boundary conditions 

% MESH_FILE='../../tests/meshes/disc_P2_h2.mesh'

if (N == 1)
  MESH_TYPE = 'GT_PK(1,1)';       % segments
  FEM_TYPE = 'FEM_HERMITE(1)';
  DATA_FEM_TYPE = 'FEM_PK(1,3)';
  INTEGRATION = 'IM_GAUSS1D(6)';  % Gauss-Legendre integration on segments
end

if (N == 2 && ~QUAD)
  MESH_TYPE = 'GT_PK(2,1)';        % triangles
  % INTEGRATION = 'IM_TRIANGLE(13)'; % quadrature of order 10 on triangles
    INTEGRATION = 'IM_STRUCTURED_COMPOSITE(IM_TRIANGLE(10), 3)';
  DATA_FEM_TYPE = 'FEM_PK(2, 5)';
  if 0   % hermite element
    FEM_TYPE = 'FEM_HERMITE(2)';   % not conformal -> problem
    DIRICHLET_FEM_TYPE = 'FEM_PK(2,3)';
  elseif 1  % argyris element
    FEM_TYPE = 'FEM_ARGYRIS';
    DIRICHLET_FEM_TYPE = 'FEM_PK(2,3)';
  else
%     FEM_TYPE = 'FEM_HCT_TRIANGLE';
    FEM_TYPE = 'FEM_REDUCED_HCT_TRIANGLE';
    DIRICHLET_FEM_TYPE = 'FEM_PK(2,3)';
    % INTEGRATION = 'IM_STRUCTURED_COMPOSITE(IM_TRIANGLE(5), 3)';
    INTEGRATION = 'IM_HCT_COMPOSITE(IM_TRIANGLE(6))';
    DATA_FEM_TYPE = 'FEM_PK(2, 3)';
  end;
end

if (N == 2 && QUAD)
  MESH_TYPE = 'GT_QK(2,1)';
  % FEM_TYPE = 'FEM_REDUCED_QUADC1_COMPOSITE';  
  FEM_TYPE = 'FEM_QUADC1_COMPOSITE';
  DATA_FEM_TYPE = 'FEM_QK(2,2)';
  INTEGRATION = 'IM_QUADC1_COMPOSITE(IM_TRIANGLE(6))';
  % DIRICHLET_FEM_TYPE = 'FEM_QK(2,1)';
end

if (N == 3 && ~QUAD)
  MESH_TYPE = 'GT_PK(3,1)';        % tetrahedrons
  INTEGRATION = 'IM_TETRAHEDRON(8)';
  DATA_FEM_TYPE = 'FEM_PK(3,3)';
  FEM_TYPE = 'FEM_HERMITE(3)';     % not conformal -> problem
  DIRICHLET_FEM_TYPE = 'FEM_PK(3,2)';
end

if (N == 3 && QUAD)
  error('Sorry, no hermite elements on hexahedrons for the moment');
end

if (N > 3)
  error('Correct the parameter file for N > 3');
end

RESIDUAL = 1E-9;     	% residual for conjugate gradient.

%%%%%   saving parameters                             %%%%%
ROOTFILENAME = 'bilaplacian';     % Root of data files.
VTK_EXPORT = 1; % export solution to a .vtk file ?