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

Last change on this file since 17806 was 17806, checked in by Mathieu Morlighem, 11 years ago

merged trunk-jpl and trunk for revision 17804
File size: 13.6 KB
Line 
1%DAMAGEICE class definition
2%
3% Usage:
4% damage=damage();
5
6classdef damage
7 properties (SetAccess=public)
8 %damage
9 isdamage = NaN;
10 D = NaN;
11 law = '';
12 spcdamage = NaN;
13 max_damage = NaN;
14
15 %numerical
16 stabilization = NaN;
17 maxiter = NaN;
18 elementinterp = '';
19 penalty_threshold = NaN;
20 penalty_lock = NaN;
21 penalty_factor = NaN;
22
23 %general parameters for evolution law:
24 stress_threshold = NaN;
25 c1 = NaN;
26 c2 = NaN;
27 c3 = NaN;
28 c4 = NaN;
29 healing = NaN;
30 equiv_stress = NaN;
31 requested_outputs = {};
32 end
33 methods
34 function createxml(obj,fid) % {{{
35 fprintf(fid, '\n\n');
36 fprintf(fid, '%s\n', '<!-- damage -->');
37 fprintf(fid, '%s\n', '<!-- Note: this class depends on different input of law -->');
38
39 %fprintf(fid,'%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="law" type="logical"', '" default="', num2str(obj.law),'">', ' <section name="damage" />',' <help> damage law (string) from {"undamaged","pralong"} </help>','</parameter>');
40 % drop-down
41 fprintf(fid,'%s%s%s%s%s\n\t%s\n','<parameter key ="law" type="','alternative','" optional="','false','">','<section name="damage" />');
42
43 % law = 'undamage'
44 fprintf(fid,'\t%s%s%s%s%s\n\t\t%s\n','<option value="undamage" type="','string','" default="','true','">','<help> law = undamage </help>');
45 % footer for option
46 fprintf(fid,'\t%s\n%s\n','</option>');
47
48 % law = 'pralong'
49 fprintf(fid,'\t%s%s%s%s%s\n\t\t%s\n','<option value="pralong" type="','string','" default="','false','">','<help> law = pralong </help>');
50
51 fprintf(fid,'\t\t%s%s%s%s%s\n\t\t\t%s\n\t\t%s\n', '<parameter key ="stress_threshold" type="',class(obj.stress_threshold),'" default="',num2str(obj.stress_threshold),'">','<help> damage stress threshold [Pa] </help>','</parameter>');
52 fprintf(fid,'\t\t%s%s%s%s%s\n\t\t\t%s\n\t\t%s\n', '<parameter key ="c1" type="', class(obj.c1),'" default="', num2str(obj.c1),'">', '<help> damage parameter 1 </help>','</parameter>');
53 fprintf(fid,'\t\t%s%s%s%s%s\n\t\t\t%s\n\t\t%s\n', '<parameter key ="c2" type="', class(obj.c2),'" default="', num2str(obj.c2),'">','<help> damage parameter 2 </help>','</parameter>');
54 fprintf(fid,'\t\t%s%s%s%s%s\n\t\t\t%s\n\t\t%s\n', '<parameter key ="c3" type="', class(obj.c3),'" default="', num2str(obj.c3),'">','<help> damage parameter 3 [W/m^2] </help>','</parameter>');
55 fprintf(fid,'\t\t%s%s%s%s%s\n\t\t\t%s\n\t\t%s\n', '<parameter key ="c4" type="', class(obj.c4),'" default="', num2str(obj.c4),'">','<help> damage parameter 4 </help>','</parameter>');
56 fprintf(fid,'\t\t%s%s%s%s%s\n\t\t\t%s\n\t\t%s\n', '<parameter key ="healing" type="', class(obj.healing),'" default="', num2str(obj.healing),'">','<help> damage healing parameter 1 </help>','</parameter>');
57 fprintf(fid,'\t\t%s%s%s%s%s\n\t\t\t%s\n\t\t%s\n', '<parameter key ="equiv_stress" type="', class(obj.equiv_stress),'" default="',convert2str(obj.equiv_stress),'">','<help> 0: von Mises </help>','</parameter>');
58 fprintf(fid,'\t\t%s%s%s%s%s\n\t\t\t%s\n\t\t%s\n', '<parameter key ="requested_outputs" type="', class(obj.requested_outputs),'" default="',convert2str(obj.requested_outputs),'">','<help> additional outputs requested </help>','</parameter>');
59
60
61 % footer for option
62 fprintf(fid,'\t%s\n%s\n','</option>');
63
64
65 % footer for drop-down
66 fprintf(fid,'\t%s\n%s\n%s','<help> damage law (string) from {"undamaged","pralong"} </help>','</parameter>');
67
68
69 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="D" type="', class(obj.D),'" default="', num2str(obj.D),'">', ' <section name="damage" />',' <help> damage tensor (scalar) </help>','</parameter>');
70 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="law" type="', class(obj.law),'" default="', num2str(obj.law),'">', ' <section name="damage" />',' <help> damage law (string) from {"undamaged","pralong"} </help>','</parameter>');
71 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="spcdamage" type="', class(obj.spcdamage),'" default="', num2str(obj.spcdamage),'">', ' <section name="damage" />',' <help> damage constraints (NaN means no constraint) </help>','</parameter>');
72 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="max_damage" type="', class(obj.max_damage),'" default="', num2str(obj.max_damage),'">', ' <section name="damage" />',' <help> maximum possible damage (0&amp;lt;=max_damage&amp;lt;1) </help>','</parameter>');
73
74 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="stabilization" type="', class(obj.stabilization),'" default="', num2str(obj.stabilization),'">', ' <section name="damage" />',' <help> 0: no, 1: artificial_diffusivity, 2: SUPG </help>','</parameter>');
75 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="maxiter" type="', class(obj.maxiter),'" default="', num2str(obj.maxiter),'">', ' <section name="damage" />',' <help> maximum number of non linear iterations </help>','</parameter>');
76 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="penalty_lock" type="', class(obj.penalty_lock),'" default="', num2str(obj.penalty_lock),'">', ' <section name="damage" />',' <help> stabilize unstable damage constraints that keep zigzagging after n iteration (default is 0, no stabilization) </help>','</parameter>');
77 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="penalty_threshold" type="', class(obj.penalty_threshold),'" default="', num2str(obj.penalty_threshold),'">', ' <section name="damage" />',' <help> threshold to declare convergence of damage evolution solution (default is 0) </help>','</parameter>');
78 fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n', '<parameter key ="penalty_factor" type="', class(obj.penalty_factor),'" default="', num2str(obj.penalty_factor),'">', ' <section name="damage" />',' <help> scaling exponent (default is 3) </help>','</parameter>');
79
80 end % }}}
81
82 function obj = damage(varargin) % {{{
83 switch nargin
84 case 0
85 obj=setdefaultparameters(obj);
86 case 1
87 inputstruct=varargin{1};
88 list1 = properties('damage');
89 list2 = fieldnames(inputstruct);
90 for i=1:length(list1)
91 fieldname = list1{i};
92 if ismember(fieldname,list2),
93 obj.(fieldname) = inputstruct.(fieldname);
94 end
95 end
96 otherwise
97 error('constructor not supported');
98 end
99 end % }}}
100 function obj = setdefaultparameters(obj) % {{{
101
102 %damage parameters:
103 obj.isdamage=0;
104 obj.D=0;
105 obj.law='undamaged';
106
107 obj.max_damage=1-1e-5; %if damage reaches 1, solve becomes singular, as viscosity becomes nil
108
109 %Type of stabilization used
110 obj.stabilization=2;
111
112 %Maximum number of iterations
113 obj.maxiter=100;
114
115 %finite element interpolation
116 obj.elementinterp='P1';
117
118 %factor used to compute the values of the penalties: kappa=max(stiffness matrix)*10^penalty_factor
119 obj.penalty_factor=3;
120
121 %stabilize unstable damage constraints that keep zigzagging after n iteration (default is 0, no stabilization)
122 obj.penalty_lock=0;
123
124 %threshold to declare convergence of thermal solution (default is 0)
125 obj.penalty_threshold=0;
126
127 %damage evolution parameters
128 obj.stress_threshold=0;
129 obj.healing=0;
130 obj.c1=0;
131 obj.c2=0;
132 obj.c3=0;
133 obj.c4=0;
134 obj.equiv_stress=0;
135
136 %output default:
137 obj.requested_outputs={'default'};
138
139 end % }}}
140 function md = checkconsistency(obj,md,solution,analyses) % {{{
141
142 md = checkfield(md,'fieldname','damage.isdamage','values',[1,0]);
143 if obj.isdamage,
144 md = checkfield(md,'fieldname','damage.law','values',{'undamaged','pralong'});
145 md = checkfield(md,'fieldname','damage.D','>=',0,'<=',obj.max_damage,'size',[md.mesh.numberofvertices 1]);
146 md = checkfield(md,'fieldname','damage.spcdamage','forcing',1);
147 md = checkfield(md,'fieldname','damage.max_damage','<',1,'>=',0);
148
149 md = checkfield(md,'fieldname','damage.stabilization','numel',[1],'values',[0 1 2]);
150 md = checkfield(md,'fieldname','damage.maxiter','>=0',0);
151 md = checkfield(md,'fieldname','damage.elementinterp','values',{'P1','P2'});
152 md = checkfield(md,'fieldname','damage.penalty_factor','>=0',0);
153 md = checkfield(md,'fieldname','damage.penalty_lock','>=0',0);
154 md = checkfield(md,'fieldname','damage.penalty_threshold','>=0',0);
155 end
156
157 if strcmpi(obj.law,'pralong'),
158 md = checkfield(md,'fieldname','damage.healing','>=',0);
159 md = checkfield(md,'fieldname','damage.c1','>=',0);
160 md = checkfield(md,'fieldname','damage.c2','>=',0);
161 md = checkfield(md,'fieldname','damage.c3','>=',0);
162 md = checkfield(md,'fieldname','damage.c4','>=',0);
163 md = checkfield(md,'fieldname','damage.stress_threshold','>=',0);
164 md = checkfield(md,'fieldname','damage.equiv_stress','numel',[1],'values',[0 1]);
165 md = checkfield(md,'fieldname','damage.requested_outputs','stringrow',1);
166 elseif strcmpi(obj.law,'undamaged'),
167 if (solution==DamageEvolutionSolutionEnum),
168 error('Invalid evolution law (md.damage.law) for a damage solution');
169 end
170 else
171 error('invalid damage evolution law');
172 end
173
174 end % }}}
175 function list=defaultoutputs(self,md) % {{{
176
177 if strcmp(domaintype(md.mesh),'2Dhorizontal'),
178 list = {'DamageDbar'};
179 else
180 list = {'DamageD'};
181 end
182 end % }}}
183 function disp(obj) % {{{
184 disp(sprintf(' Damage:\n'));
185
186 fielddisplay(obj,'isdamage','is damage mechanics being used? {true,false}');
187 if obj.isdamage,
188 fielddisplay(obj,'law','damage law (string) from {''undamaged'',''pralong''}');
189 fielddisplay(obj,'D','damage tensor (scalar)');
190 fielddisplay(obj,'spcdamage','damage constraints (NaN means no constraint)');
191 fielddisplay(obj,'max_damage','maximum possible damage (0<=max_damage<1)');
192
193 fielddisplay(obj,'stabilization','0: no, 1: artificial_diffusivity, 2: SUPG');
194 fielddisplay(obj,'maxiter','maximum number of non linear iterations');
195 fielddisplay(obj,'elementinterp','interpolation scheme for finite elements {''P1'',''P2''}');
196 fielddisplay(obj,'penalty_lock','stabilize unstable damage constraints that keep zigzagging after n iteration (default is 0, no stabilization)');
197 fielddisplay(obj,'penalty_threshold','threshold to declare convergence of damage evolution solution (default is 0)');
198 fielddisplay(obj,'penalty_factor','scaling exponent (default is 3)');
199 end
200
201 if strcmpi(obj.law,'pralong'),
202 fielddisplay(obj,'c1','damage parameter 1');
203 fielddisplay(obj,'c2','damage parameter 2');
204 fielddisplay(obj,'c3','damage parameter 3');
205 fielddisplay(obj,'c4','damage parameter 4');
206 fielddisplay(obj,'healing','damage healing parameter');
207 fielddisplay(obj,'stress_threshold','damage stress threshold [Pa]');
208 fielddisplay(obj,'equiv_stress','0: von Mises, 1: max principal');
209 fielddisplay(obj,'requested_outputs','additional outputs requested');
210 end
211
212 end % }}}
213 function marshall(obj,md,fid) % {{{
214
215 WriteData(fid,'object',obj,'fieldname','isdamage','format','Boolean');
216 if obj.isdamage,
217 WriteData(fid,'object',obj,'fieldname','law','format','String');
218 WriteData(fid,'object',obj,'fieldname','D','format','DoubleMat','mattype',1);
219 WriteData(fid,'object',obj,'fieldname','spcdamage','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1);
220 WriteData(fid,'object',obj,'fieldname','max_damage','format','Double');
221
222 WriteData(fid,'object',obj,'fieldname','stabilization','format','Integer');
223 WriteData(fid,'object',obj,'fieldname','maxiter','format','Integer');
224 WriteData(fid,'object',obj,'fieldname','elementinterp','format','String');
225 WriteData(fid,'object',obj,'fieldname','penalty_threshold','format','Integer');
226 WriteData(fid,'object',obj,'fieldname','penalty_lock','format','Integer');
227 WriteData(fid,'object',obj,'fieldname','penalty_factor','format','Double');
228 end
229
230 if strcmpi(obj.law,'pralong'),
231 WriteData(fid,'object',obj,'fieldname','c1','format','Double');
232 WriteData(fid,'object',obj,'fieldname','c2','format','Double');
233 WriteData(fid,'object',obj,'fieldname','c3','format','Double');
234 WriteData(fid,'object',obj,'fieldname','c4','format','Double');
235 WriteData(fid,'object',obj,'fieldname','stress_threshold','format','Double');
236 WriteData(fid,'object',obj,'fieldname','healing','format','Double');
237 WriteData(fid,'object',obj,'fieldname','equiv_stress','format','Integer');
238 end
239
240 %process requested outputs
241 outputs = obj.requested_outputs;
242 pos = find(ismember(outputs,'default'));
243 if ~isempty(pos),
244 outputs(pos) = []; %remove 'default' from outputs
245 outputs = [outputs defaultoutputs(obj,md)]; %add defaults
246 end
247 if obj.isdamage,
248 WriteData(fid,'data',outputs,'enum',DamageEvolutionRequestedOutputsEnum,'format','StringArray');
249 end
250
251 end % }}}
252 end
253end
Note: See TracBrowser for help on using the repository browser.