%--------------------------------------------------------------------------
% PURPOSE
% Uniaxial test- ZAFATI's model
%--------------------------------------------------------------------------
% REFERENCES
% Eliass ZAFATI
% 05-12-2016
%--------------------------------------------------------------------------
% COMMENTS
%
%
%
%--------------------------------------------------------------------------
% MIT License
%
% Copyright (c) 2018 Benjamin Richard
%
% Permission is hereby granted, free of charge, to any person obtaining a
% copy of this software and associated documentation files (the "Software
% "), to deal in the Software without restriction, including without
% limitation the rights to use, copy, modify, merge, publish, distribute,
% sublicense, and/or sell copies of the Software, and to permit persons
% to whom the Software is furnished to do so, subject to the following
% conditions:
%
% The above copyright notice and this permission notice shall be included
% in all copies or substantial portions of the Software.
%
% THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
% OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
% MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
% IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
% CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
% TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
% SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
%--------------------------------------------------------------------------
%% Clearing off
fclose all;
clear
close all
%% Declaration de variables global
global options ME TP;
%% Definition des options
options.mode = 'TRID';
%% Loading of the input datafile
FILE = '69.mail';
ME = INPUT.ACQU(FILE,'MAIL');
%% Definition of the model
MO1 = MODEL('VTOT','MECHANICS','ELASTICITY','ISOTROPIC','DAMAGE','ZAFATI');
MOT = MO1;
%% Topology
TP = TOPOLOGY(MOT);
%% Definition of the material
MA1 = CHAMELEM.MATE(MOT,'youn',36000e6,'nu',0.2,'rho',2200,'k0',1e-4,'S',1e-4,'s',5,'B',2,'eta',0,'c2',1e4);
MAT = MA1 ;
%% Boundary conditions
% Line L1 fixed
CL1 = MATRICE('DIRI','P1',1);
CL2 = MATRICE('DIRI','P4',1);
CL3 = MATRICE('DIRI','P5',1);
CL4 = MATRICE('DIRI','P8',1);
CLT1 = [CL1 [[CL2 CL3] CL4]];
CL5 = MATRICE('DIRI','P2',1);
CL6 = MATRICE('DIRI','P3',1);
CL7 = MATRICE('DIRI','P6',1);
CL8 = MATRICE('DIRI','P7',1);
CLT2 = [CL5 [[CL6 CL7] CL8]];
CL9 = MATRICE('DIRI','P1',2,3);
CLT3 = CL9;
CLT = [CLT1 [CLT2 CLT3 ]];
% Definition of a prescribed displacement
FO1 = CHPOINT.DEPI(CLT2,1);
%% Loading
EV1 = EVOL([0 1 2],[0 3e-4 -3e-4],'Time','Displacement (m)');
CHT = TIMELOAD(FO1,EV1,'DIRI');
%% Static analysis
PB1 = PROBLEM('model',MOT,'mater',MAT,'diric',CLT,'loadt',CHT,'comp_time',0:0.01:2,...
'solve_type','QUASI_NEWTON','precision',1.0e-4,'prec_glob',1.0e-4);
SOL = SOLVE(PB1);
%% Post-treatment
% Reaction curve
EV_OUT = EVOL.REAC(SOL,CLT2,EV1,1);
plot(EV_OUT);
close all
%% Non regression test
if abs(EV_OUT.ordo(end) + 1.086979361716937e+07) > 1.0e-4
error('TEST IS NOT SUCCESSFUL')
else
disp('---------------------------------')
disp('069_INCR_NLINE_2D_TRI3_ZAFATI_UNIAXIAL_TEST')
disp('TEST IS SUCCESSFUL')
disp('---------------------------------')
end