46
47
48
54 USE mat_elem_mod
56
57
58
59#include "implicit_f.inc"
60
61
62
63#include "com04_c.inc"
64#include "scr17_c.inc"
65#include "units_c.inc"
66#include "param_c.inc"
67#include "sphcom.inc"
68
69
70
71 INTEGER OFFS
72 my_real :: rnoise(nperturb,numelc+numeltg+numels+numsph),
73 . qp_rperturb(nperturb,4)
74 INTEGER IPART(,*),IPARTC(*),IPARTSP(*),IPARTG(*),IPARTS(*),
75 . PERTURB(NPERTURB),
76 . IDPERTURB(NPERTURB),INDEX(NUMELC+NUMELTG+NUMELS+),
77 . INDEX_ITYP(NUMELC+NUMELTG+NUMELS+NUMSPH),NPART_SHELL,
78 . QP_IPERTURB(NPERTURB,6)
79 TYPE (UNIT_TYPE_) ,INTENT(IN) :: UNITAB
80 TYPE (SUBMODEL_DATA) ,INTENT(IN) :: LSUBMODEL(*)
81 TYPE () ,DIMENSION(NGRPART):: IGRPART
82 TYPE (MATPARAM_STRUCT_) ,DIMENSION(NUMMAT) ,INTENT(INOUT) :: MAT_PARAM
83
84
85
86 INTEGER ICOUNT,II,J,K,N,I_METHOD,MAX_PART,OPT_ID,FAIL_ID,UNIT_ID,KLEN,
87 . CPT_PART,NB_RANDOM,I_SEED,NPERTURB_FAIL,
88 . NB_INTERV,SEED,SEED_RANDOM,IFAILMAT,IFAILTYPE,ITYP,I_PERTURB_VAR,SIZEY
89 INTEGER, DIMENSION(50) :: DISTRIB
90 INTEGER, DIMENSION(8) :: DT_SEED
91 INTEGER, DIMENSION(:,:), ALLOCATABLE :: TAB_PART
92 INTEGER, DIMENSION(:), ALLOCATABLE :: A_SEED
94 . mean,stdev,mean_input,sd_input,max_distrib,temp,min_value,
95 . max_value,interv,VALUE,max_value1,minval,maxval,bid
96 my_real,
DIMENSION(:),
ALLOCATABLE :: array
97 CHARACTER(LEN=NCHARTITLE)::TITR
98 CHARACTER*100 KEY1,KEY2
99 CHARACTER(LEN=NCHARKEY) :: PARAM
100 CHARACTER MESS*40
101 LOGICAL IS_AVAILABLE
102 DATA mess/'PERTURBATION DEFINITION '/
103
105
106
107
109 max_part = 0
110 DO icount = 1+offs,nperturb_fail+offs
111 ifailtype = 0
112 cpt_part = 0
113 ityp = 2
115 . option_id = opt_id,
116 . unit_id = unit_id ,
117 . keyword2 = key1,
118 . keyword3 = key2,
119 . option_titr = titr)
120
121 idperturb(icount) = opt_id
122 klen = len_trim(key2)
123 IF (key2(1:klen) == 'BIQUAD') THEN
124 ifailtype = 30
125 ELSE
126 CALL ancmsg(msgid=1192, msgtype=msgerror, anmode=aninfo,
127 . i1=opt_id,
128 . c1=titr,
129 . c2=key2)
130 cycle
131 ENDIF
132
133 CALL hm_get_intv (
'fail_ID' ,fail_id ,is_available,lsubmodel)
134
135 ifailmat = 0
136 IF (fail_id > 0) THEN
137 DO n=1,nummat
138 DO j=1,mat_param(n)%NFAIL
139 IF (mat_param(n)%FAIL(j)%FAIL_ID == fail_id)THEN
140 IF (ifailtype /= mat_param(n)%FAIL(j)%IRUPT) THEN
141 CALL ancmsg(msgid=1193, msgtype=msgerror, anmode=aninfo,
142 . i1=opt_id,
143 . c1=titr,
144 . i2=fail_id,
145 . c2=key2)
146 ENDIF
147 ifailmat = n
148 EXIT
149 END IF
150 END DO
151 IF (ifailmat > 0) EXIT
152 END DO
153 ENDIF
154
155 perturb(icount) = ityp
156
157 IF (ifailmat > 0) THEN
158 DO n=1,npart
159 IF(ipart(1,n) == ifailmat) THEN
160 cpt_part = cpt_part + 1
161 ENDIF
162 ENDDO
163 ELSE
164 CALL ancmsg(msgid=1137, msgtype=msgerror, anmode=aninfo,
165 . i1=opt_id,
166 . c1=titr,
167 . i2=fail_id,
168 . c2='FAILURE CRITERIA')
169 ENDIF
170 max_part =
max(max_part,cpt_part)
171 ENDDO
172
173
174 ALLOCATE(tab_part(nperturb,max_part))
175
176
177
178
180 DO icount = 1+offs,nperturb_fail+offs
181 i_perturb_var = 0
182 cpt_part = 0
183 ityp = 2
185 . option_id = opt_id,
186 . unit_id = unit_id ,
187 . keyword2 = key1,
188 . keyword3 = key2,
189 . option_titr = titr)
190 idperturb(icount) = opt_id
191
192 klen = len_trim(key2)
193 IF (key2(1:klen) == 'BIQUAD') THEN
194 ifailtype = 30
195
196 CALL hm_get_floatv(
'F_Mean' ,mean ,is_available,lsubmodel,unitab)
197 CALL hm_get_floatv(
'Deviation' ,stdev ,is_available,lsubmodel,unitab)
199 CALL hm_get_floatv(
'Max_cut' ,maxval ,is_available,lsubmodel,unitab)
201 CALL hm_get_intv (
'Idistri' ,i_method ,is_available,lsubmodel)
202
203 CALL hm_get_intv (
'fail_ID' ,fail_id ,is_available,lsubmodel)
205
206
207
208 IF (i_method == 0) i_method = 2
209 IF (i_method == 2) THEN
210 IF (minval == zero .AND. maxval == zero) THEN
211 minval =-ep20
212 maxval = ep20
213 ENDIF
214 ENDIF
215 sd_input = stdev
216 mean_input = mean
217
218 qp_iperturb(icount,1) = opt_id
219 qp_iperturb(icount,2) = ityp
220 qp_iperturb(icount,3) =
seed
221 qp_iperturb(icount,4) = i_method
222 qp_iperturb(icount,5) = fail_id
223 qp_rperturb(icount,1) = mean
224 qp_rperturb(icount,2) = stdev
225 qp_rperturb(icount,3) = minval
226 qp_rperturb(icount,4) = maxval
227
228 IF (param(1:2) == 'c3' .or. param(1:2) == 'C3') THEN
229 i_perturb_var = 1
230 qp_iperturb(icount,6) = i_perturb_var
231 ELSE
232 CALL ancmsg(msgid=1194,msgtype=msgerror,anmode=aninfo,
233 . i1=opt_id,
234 . c1=titr,
235 . c2=param)
236 END IF
237
238 ifailmat = 0
239 IF (fail_id > 0)THEN
240 DO n=1,nummat
241 DO j=1,mat_param(n)%NFAIL
242 IF (mat_param(n)%FAIL(j)%FAIL_ID == fail_id) THEN
243 ifailmat = n
244 EXIT
245 END IF
246 END DO
247 IF (ifailmat > 0) EXIT
248 END DO
249 ENDIF
250 IF (ifailmat > 0) THEN
251 cpt_part = 0
252 DO n=1,npart
253 IF(ipart(1,n) == ifailmat) THEN
254 cpt_part = cpt_part + 1
255 tab_part(icount,cpt_part) = n
256 ENDIF
257 ENDDO
258 END IF
259
260
261
262 IF(i_method == 2) THEN
263 WRITE (iout,4000)
264 . opt_id,
'GAUSSIAN',mean_input,sd_input,
seed,key2,fail_id,param
265 WRITE (iout,'(10I10)') ipart(4,tab_part(icount,1:cpt_part))
266 WRITE(iout,*) ' '
267 WRITE(iout,*) ' '
268 ELSEIF(i_method == 1) THEN
269 WRITE (iout,4100) opt_id,
'RANDOM',
seed,key2,fail_id,param
270 WRITE (iout,'(10I10)') ipart(4,tab_part(icount,1:cpt_part))
271 WRITE(iout,*) ' '
272 WRITE(iout,*) ' '
273 ENDIF
274
275 nb_random = 0
276 DO ii=1,numelc
277 DO k=1,cpt_part
278 IF (ipartc(ii) == tab_part(icount,k)) THEN
279 nb_random = nb_random + 1
280 index(nb_random) = ii
281 index_ityp(nb_random) = 3
282 ENDIF
283 ENDDO
284 ENDDO
285 DO ii=1,numeltg
286 DO k=1,cpt_part
287 IF(ipartg(ii) == tab_part(icount,k)) THEN
288 nb_random = nb_random + 1
289 index(nb_random) = ii
290 index_ityp(nb_random) = 7
291 ENDIF
292 ENDDO
293 ENDDO
294 DO ii=1,numels
295 DO k=1,cpt_part
296 IF (iparts(ii) == tab_part(icount,k)) THEN
297 nb_random = nb_random + 1
298 index(nb_random) = ii
299 index_ityp(nb_random) = 1
300 ENDIF
301 ENDDO
302 ENDDO
303 DO ii=1,numsph
304 DO k=1,cpt_part
305 IF (ipartsp(ii) == tab_part(icount,k)) THEN
306 nb_random = nb_random + 1
307 index(nb_random) = ii
308 index_ityp(nb_random) = 51
309 ENDIF
310 ENDDO
311 ENDDO
312
313
314
316 CALL random_seed(size=i_seed)
317 ALLOCATE(a_seed(1:i_seed))
318 CALL random_seed(get=a_seed)
319 CALL date_and_time(values=dt_seed)
320 a_seed=dt_seed(8)*dt_seed(7)*dt_seed(6)
321 seed=dt_seed(8)*dt_seed(7)*dt_seed(6)
322 CALL random_seed(put=a_seed)
323 seed_random = 1
324 DEALLOCATE(a_seed)
325 ELSE
326 CALL random_seed(size=i_seed)
327 ALLOCATE(a_seed(1:i_seed))
329 CALL random_seed(put=a_seed)
330 seed_random = 0
331 DEALLOCATE(a_seed)
332 ENDIF
333
334
335
336 distrib(1:50) = 0
337 ALLOCATE(array(nb_random+2))
338 CALL random_number(array)
339
340
341
342 max_value = -ep30
343 min_value = ep30
344 IF ( i_method == 2) THEN
345 DO ii = 1, nb_random+1, 2
346 temp = stdev * sqrt(-2.0*log(array(ii))) * cos(2*pi*array(ii+1)) +
347 . mean
348 array(ii+1) =
349 . stdev * sqrt(-2.0*log(array(ii))) * sin(2*pi*array(ii+1)) + mean
350 array(ii) = temp
351 END DO
352 DO ii = 1, nb_random
353 array(ii) =
max(
min(maxval,array(ii)),minval)
354 max_value =
max(array(ii),max_value)
355 min_value =
min(array(ii),min_value)
356 END DO
357 ELSEIF(i_method == 1)THEN
358 DO ii = 1, nb_random
359 array(ii) = array(ii)*(maxval-minval)+minval
360 max_value =
max(array(ii),max_value)
361 min_value =
min(array(ii),min_value)
362 END DO
363 ENDIF
364
365 DO ii = 1, nb_random
366 IF (index_ityp(ii) == 3) THEN
367 rnoise(icount,index(ii)) = array(ii)
368 ELSEIF (index_ityp(ii) == 7) THEN
369 rnoise(icount,index(ii)+numelc) = array(ii)
370 ELSEIF (index_ityp(ii) == 1) THEN
371 rnoise(icount,index(ii)+numelc+numeltg) = array(ii)
372 ELSEIF (index_ityp(ii) == 51) THEN
373 rnoise(icount,index(ii)+numelc+numeltg+numels) = array(ii)
374 ENDIF
375 ENDDO
376
377
378
379 mean = sum(array)/nb_random
380 stdev = sqrt(sum((array - mean)**2)/nb_random)
381
382
383
384 IF (i_method == 2) THEN
385 max_distrib = one /(stdev*sqrt(two * pi))
386 ELSEIF(i_method == 1) THEN
387 max_distrib = one /(max_value-min_value)
388 ENDIF
389 IF(ityp == 2)THEN
390 WRITE (iout,5000)'C3',fail_id
391 ENDIF
392 WRITE(iout,*) ' '
393 nb_interv = 50
394 sizey = 20
395 IF (minval /= -ep30 .AND. maxval /= ep30) THEN
396 min_value = minval
397 max_value = maxval
398 ENDIF
399 CALL plot_distrib( array,nb_random, nb_interv,sizey,min_value,
400 . max_value,max_distrib,'#')
401
402 IF (i_method == 2) THEN
403 WRITE (iout,2000) mean,stdev
404 ELSEIF (i_method == 1) THEN
405 WRITE (iout,2050) mean
406 ENDIF
407
408 IF (seed_random == 1)
WRITE (iout,2100)
seed
409 WRITE(iout,*) ' '
410 WRITE(iout,*) ' '
411 IF (ALLOCATED(array)) DEALLOCATE(array)
412
413 END IF
414
415 ENDDO
416
417 DEALLOCATE(tab_part)
418
419
420 4000 FORMAT(/' PERTURBATION ID',i10/
421 . ' ---------------'/
422 . ' TYPE . . . . . . . . . . . . . . .',a/
423 . ' INPUT MEAN VALUE . . . . . . . . .',1pg20.13/
424 . ' INPUT STANDARD DEVIATION . . . . .',1pg20.13/
425 . ' INPUT SEED VALUE . . . . . . . . .',i10/
426 . ' FAILURE CRITERIA . . . . . . . . .',a/
427 . ' FAILURE CRITERIA ID. . . . . . . .',i10/
428 . ' APPLIED ON PARAMETER . . . . . . .',a/
429 . ' PARTS:')
430 4100 FORMAT(/' PERTURBATION ID',i10/
431 . ' ---------------'/
432 . ' TYPE . . . . . . . . . . . . . . .',a/
433 . ' INPUT SEED VALUE . . . . . . . . .',i10/
434 . ' FAILURE CRITERIA . . . . . . . . .',a/
435 . ' FAILURE CRITERIA ID. . . . . . . .',i10/
436 . ' APPLIED ON PARAMETER . . . . . . .',a/
437 . ' PARTS:')
438
439 2000 FORMAT(/
440 . ' GENERATED MEAN VALUE . . . . . . .',1pg20.13/
441 . ' GENERATED STANDARD DEVIATION . . .',1pg20.13)
442 2050 FORMAT(/
443 . ' GENERATED MEAN VALUE . . . . . . .',1pg20
444 2100 FORMAT(/
445 . ' GENERATED SEED VALUE . . . . . . .',i10/)
446
447 5000 FORMAT(/
448 . ' DISTRIBUTION OF SCALE FACTORS APPLIED TO ',a,' VALUE'/
449 . ' OF FAILURE CRITERIA ID= . . . . . .',i10)
450
451 RETURN
subroutine hm_get_floatv(name, rval, is_available, lsubmodel, unitab)
subroutine hm_get_intv(name, ival, is_available, lsubmodel)
subroutine hm_get_string(name, sval, size, is_available)
subroutine hm_option_count(entity_type, hm_option_number)
subroutine hm_option_start(entity_type)
integer, parameter nchartitle
integer, parameter ncharkey
subroutine plot_distrib(array, s_array, nb_interv, sizey, x_minvalue, x_maxvalue, y_maxvalue, ecrit)
subroutine ancmsg(msgid, msgtype, anmode, i1, i2, i3, i4, i5, i6, i7, i8, i9, i10, i11, i12, i13, i14, i15, i16, i17, i18, i19, i20, r1, r2, r3, r4, r5, r6, r7, r8, r9, c1, c2, c3, c4, c5, c6, c7, c8, c9, prmode)