source: issm/trunk/src/m/classes/damage.m@ 24313

Last change on this file since 24313 was 24313, checked in by Mathieu Morlighem, 5 years ago

merged trunk-jpl and trunk for revision 24310

File size: 7.3 KB
RevLine 
[16160]1%DAMAGEICE class definition
2%
3% Usage:
4% damage=damage();
5
6classdef damage
7 properties (SetAccess=public)
8 %damage
[17989]9 isdamage = 0;
[16160]10 D = NaN;
[18301]11 law = 0;
[16181]12 spcdamage = NaN;
[18301]13 max_damage = 0;
[16188]14
15 %numerical
[18301]16 stabilization = 0;
17 maxiter = 0;
[17806]18 elementinterp = '';
[16160]19
[16198]20 %general parameters for evolution law:
[18301]21 stress_threshold = 0;
[24313]22 stress_ubound = 0;
[20500]23 kappa = 0;
[18301]24 c1 = 0;
25 c2 = 0;
26 c3 = 0;
27 c4 = 0;
28 healing = 0;
29 equiv_stress = 0;
[17806]30 requested_outputs = {};
[16218]31 end
[16160]32 methods
[19105]33 function self = damage(varargin) % {{{
[16160]34 switch nargin
35 case 0
[19105]36 self=setdefaultparameters(self);
[16160]37 case 1
38 inputstruct=varargin{1};
39 list1 = properties('damage');
40 list2 = fieldnames(inputstruct);
41 for i=1:length(list1)
42 fieldname = list1{i};
43 if ismember(fieldname,list2),
[19105]44 self.(fieldname) = inputstruct.(fieldname);
[16160]45 end
46 end
47 otherwise
48 error('constructor not supported');
49 end
50 end % }}}
[19105]51 function self = extrude(self,md) % {{{
52 self.D=project3d(md,'vector',self.D,'type','node');
53 self.spcdamage=project3d(md,'vector',self.spcdamage,'type','node');
54 end % }}}
55 function self = setdefaultparameters(self) % {{{
[16160]56
57 %damage parameters:
[19105]58 self.isdamage=0;
59 self.D=0;
60 self.law=0;
[16188]61
[19105]62 self.max_damage=1-1e-5; %if damage reaches 1, solve becomes singular, as viscosity becomes nil
[16181]63
64 %Type of stabilization used
[20500]65 self.stabilization=4;
[16188]66
67 %Maximum number of iterations
[19105]68 self.maxiter=100;
[16188]69
[17806]70 %finite element interpolation
[19105]71 self.elementinterp='P1';
[17806]72
[16198]73 %damage evolution parameters
[20500]74 self.stress_threshold=1.3e5;
75 self.kappa=2.8;
[19105]76 self.healing=0;
77 self.c1=0;
78 self.c2=0;
79 self.c3=0;
80 self.c4=0;
81 self.equiv_stress=0;
[16160]82
[17806]83 %output default:
[19105]84 self.requested_outputs={'default'};
[17806]85
[16160]86 end % }}}
[19105]87 function md = checkconsistency(self,md,solution,analyses) % {{{
[16160]88
[17806]89 md = checkfield(md,'fieldname','damage.isdamage','values',[1,0]);
[19105]90 if self.isdamage,
[24313]91 md = checkfield(md,'fieldname','damage.law','numel',[1],'values',[0,1,2,3]);
[19105]92 md = checkfield(md,'fieldname','damage.D','>=',0,'<=',self.max_damage,'size',[md.mesh.numberofvertices 1]);
[20500]93 md = checkfield(md,'fieldname','damage.spcdamage','Inf',1,'timeseries',1);
[17806]94 md = checkfield(md,'fieldname','damage.max_damage','<',1,'>=',0);
[20500]95 md = checkfield(md,'fieldname','damage.stabilization','numel',[1],'values',[0 1 2 4]);
[17806]96 md = checkfield(md,'fieldname','damage.maxiter','>=0',0);
97 md = checkfield(md,'fieldname','damage.elementinterp','values',{'P1','P2'});
[20500]98 md = checkfield(md,'fieldname','damage.stress_threshold','>=',0);
[24313]99 md = checkfield(md,'fieldname','damage.stress_ubound','>=',0);
[20500]100 md = checkfield(md,'fieldname','damage.kappa','>',1);
[17806]101 md = checkfield(md,'fieldname','damage.healing','>=',0);
102 md = checkfield(md,'fieldname','damage.c1','>=',0);
103 md = checkfield(md,'fieldname','damage.c2','>=',0);
104 md = checkfield(md,'fieldname','damage.c3','>=',0);
105 md = checkfield(md,'fieldname','damage.c4','>=',0);
106 md = checkfield(md,'fieldname','damage.equiv_stress','numel',[1],'values',[0 1]);
107 md = checkfield(md,'fieldname','damage.requested_outputs','stringrow',1);
[19105]108 elseif (self.law~=0),
[21341]109 if (strcmp(solution,'DamageEvolutionSolution')),
[16198]110 error('Invalid evolution law (md.damage.law) for a damage solution');
[16188]111 end
[16160]112 end
113 end % }}}
[17806]114 function list=defaultoutputs(self,md) % {{{
115
116 if strcmp(domaintype(md.mesh),'2Dhorizontal'),
117 list = {'DamageDbar'};
118 else
119 list = {'DamageD'};
120 end
121 end % }}}
[19105]122 function disp(self) % {{{
[16160]123 disp(sprintf(' Damage:\n'));
124
[19105]125 fielddisplay(self,'isdamage','is damage mechanics being used? {true,false}');
126 if self.isdamage,
[20500]127 fielddisplay(self,'law','damage law {''0: analytical'',''1: pralong''}');
[19105]128 fielddisplay(self,'D','damage tensor (scalar)');
129 fielddisplay(self,'spcdamage','damage constraints (NaN means no constraint)');
130 fielddisplay(self,'max_damage','maximum possible damage (0<=max_damage<1)');
[17806]131
[20500]132 fielddisplay(self,'stabilization','0: no, 1: artificial_diffusivity, 2: SUPG (not working), 4: flux corrected transport');
[19105]133 fielddisplay(self,'maxiter','maximum number of non linear iterations');
134 fielddisplay(self,'elementinterp','interpolation scheme for finite elements {''P1'',''P2''}');
[24313]135 fielddisplay(self,'stress_threshold','stress threshold for damage initiation (Pa)');
136 fielddisplay(self,'stress_ubound','stress upper bound for damage healing (Pa), arctan law');
[20500]137 fielddisplay(self,'kappa','ductility parameter for stress softening and damage');
[19105]138 fielddisplay(self,'c1','damage parameter 1');
139 fielddisplay(self,'c2','damage parameter 2');
140 fielddisplay(self,'c3','damage parameter 3');
141 fielddisplay(self,'c4','damage parameter 4');
142 fielddisplay(self,'healing','damage healing parameter');
143 fielddisplay(self,'equiv_stress','0: von Mises, 1: max principal');
144 fielddisplay(self,'requested_outputs','additional outputs requested');
[16160]145 end
146
147 end % }}}
[21341]148 function marshall(self,prefix,md,fid) % {{{
[16162]149
[21341]150 WriteData(fid,prefix,'object',self,'fieldname','isdamage','format','Boolean');
[19105]151 if self.isdamage,
[21341]152 WriteData(fid,prefix,'object',self,'fieldname','law','format','Integer');
153 WriteData(fid,prefix,'object',self,'fieldname','D','format','DoubleMat','mattype',1);
154 WriteData(fid,prefix,'object',self,'fieldname','spcdamage','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
155 WriteData(fid,prefix,'object',self,'fieldname','max_damage','format','Double');
[16188]156
[21341]157 WriteData(fid,prefix,'object',self,'fieldname','stabilization','format','Integer');
158 WriteData(fid,prefix,'object',self,'fieldname','maxiter','format','Integer');
159 WriteData(fid,prefix,'name','md.damage.elementinterp','data',self.elementinterp,'format','String');
160 WriteData(fid,prefix,'object',self,'fieldname','stress_threshold','format','Double');
[24313]161 WriteData(fid,prefix,'object',self,'fieldname','stress_ubound','format','Double');
[21341]162 WriteData(fid,prefix,'object',self,'fieldname','kappa','format','Double');
163 WriteData(fid,prefix,'object',self,'fieldname','c1','format','Double');
164 WriteData(fid,prefix,'object',self,'fieldname','c2','format','Double');
165 WriteData(fid,prefix,'object',self,'fieldname','c3','format','Double');
166 WriteData(fid,prefix,'object',self,'fieldname','c4','format','Double');
167 WriteData(fid,prefix,'object',self,'fieldname','healing','format','Double');
168 WriteData(fid,prefix,'object',self,'fieldname','equiv_stress','format','Integer');
[16160]169 end
170
[17806]171 %process requested outputs
[19105]172 outputs = self.requested_outputs;
[17806]173 pos = find(ismember(outputs,'default'));
174 if ~isempty(pos),
175 outputs(pos) = []; %remove 'default' from outputs
[19105]176 outputs = [outputs defaultoutputs(self,md)]; %add defaults
[17806]177 end
[19105]178 if self.isdamage,
[21341]179 WriteData(fid,prefix,'data',outputs,'name','md.damage.requested_outputs','format','StringArray');
[17806]180 end
181
[16160]182 end % }}}
[20500]183 function savemodeljs(self,fid,modelname) % {{{
184
185 writejsdouble(fid,[modelname '.damage.isdamage'],self.isdamage);
186 if self.isdamage,
187 error('savemodeljs error message: not implemented yet!');
188 end
189
190 end % }}}
[16160]191 end
192end
Note: See TracBrowser for help on using the repository browser.