OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
hm_read_perturb_fail.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| hm_read_perturb_fail ../starter/source/general_controls/computation/hm_read_perturb_fail.F
25!||--- called by ------------------------------------------------------
26!|| hm_read_perturb ../starter/source/general_controls/computation/hm_read_perturb.F
27!||--- calls -----------------------------------------------------
28!|| ancmsg ../starter/source/output/message/message.F
29!|| hm_get_floatv ../starter/source/devtools/hm_reader/hm_get_floatv.F
30!|| hm_get_intv ../starter/source/devtools/hm_reader/hm_get_intv.F
31!|| hm_get_string ../starter/source/devtools/hm_reader/hm_get_string.F
32!|| hm_option_count ../starter/source/devtools/hm_reader/hm_option_count.F
33!|| hm_option_read_key ../starter/source/devtools/hm_reader/hm_option_read_key.F
34!|| hm_option_start ../starter/source/devtools/hm_reader/hm_option_start.F
35!|| plot_distrib ../starter/source/general_controls/computation/plot_distrib.F
36!||--- uses -----------------------------------------------------
37!|| hm_option_read_mod ../starter/share/modules1/hm_option_read_mod.F
38!|| message_mod ../starter/share/message_module/message_mod.F
39!|| submodel_mod ../starter/share/modules1/submodel_mod.f
40!||====================================================================
41 SUBROUTINE hm_read_perturb_fail(MAT_PARAM,
42 . IPART ,RNOISE ,IPARTC ,IPARTG ,IPARTSP ,
43 . IGRPART ,IPARTS ,PERTURB ,IDPERTURB,
44 . INDEX ,INDEX_ITYP,NPART_SHELL,OFFS ,QP_IPERTURB,
45 . QP_RPERTURB,LSUBMODEL,UNITAB)
46C-----------------------------------------------
47C M o d u l e s
48C-----------------------------------------------
49 USE message_mod
50 USE groupdef_mod
51 USE unitab_mod
52 USE submodel_mod
54 USE mat_elem_mod
56C-----------------------------------------------
57C I m p l i c i t T y p e s
58C-----------------------------------------------
59#include "implicit_f.inc"
60C-----------------------------------------------
61C C o m m o n B l o c k s
62C-----------------------------------------------
63#include "com04_c.inc"
64#include "scr17_c.inc"
65#include "units_c.inc"
66#include "param_c.inc"
67#include "sphcom.inc"
68C-----------------------------------------------
69C D u m m y A r g u m e n t s
70C-----------------------------------------------
71 INTEGER OFFS
72 my_real :: RNOISE(NPERTURB,NUMELC+NUMELTG+NUMELS+NUMSPH),
73 . QP_RPERTURB(NPERTURB,4)
74 INTEGER IPART(LIPART1,*),IPARTC(*),IPARTSP(*),IPARTG(*),IPARTS(*),
75 . perturb(nperturb),
76 . idperturb(nperturb),index(numelc+numeltg+numels+numsph),
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 (GROUP_) ,DIMENSION(NGRPART):: IGRPART
82 TYPE (MATPARAM_STRUCT_) ,DIMENSION(NUMMAT) ,INTENT(INOUT) :: MAT_PARAM
83C-----------------------------------------------
84C L o c a l V a r i a b l e s
85C-----------------------------------------------
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
93 my_real
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 '/
103C=======================================================================
104 CALL hm_option_count('/PERTURB/FAIL',nperturb_fail)
105c-----------------------------------------------------------------------
106c 1st loop over /PERTURB/FAIL/BIQUAD for computing table dimension
107c-----------------------------------------------------------------------
108 CALL hm_option_start('/PERTURB/FAIL')
109 max_part = 0
110 DO icount = 1+offs,nperturb_fail+offs
111 ifailtype = 0
112 cpt_part = 0
113 ityp = 2
114 CALL hm_option_read_key(lsubmodel,
115 . option_id = opt_id,
116 . unit_id = unit_id ,
117 . keyword2 = key1,
118 . keyword3 = key2,
119 . option_titr = titr)
120c
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
132c
133 CALL hm_get_intv ('fail_ID' ,fail_id ,is_available,lsubmodel)
134c
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
154c
155 perturb(icount) = ityp
156c
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
172c-------- end first loop over /PERTURB/FAIL
173c
174 ALLOCATE(tab_part(nperturb,max_part))
175c
176c-----------------------------------------------------------------------
177c 2nd loop over /PERTURB/FAIL/BIQUAD for reading and computing perturbation
178c-----------------------------------------------------------------------
179 CALL hm_option_start('/PERTURB/FAIL')
180 DO icount = 1+offs,nperturb_fail+offs
181 i_perturb_var = 0
182 cpt_part = 0
183 ityp = 2
184 CALL hm_option_read_key(lsubmodel,
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
191c
192 klen = len_trim(key2)
193 IF (key2(1:klen) == 'BIQUAD') THEN
194 ifailtype = 30
195card1
196 CALL hm_get_floatv('F_Mean' ,mean ,is_available,lsubmodel,unitab)
197 CALL hm_get_floatv('Deviation' ,stdev ,is_available,lsubmodel,unitab)
198 CALL hm_get_floatv('Min_cut' ,minval ,is_available,lsubmodel,unitab)
199 CALL hm_get_floatv('Max_cut' ,maxval ,is_available,lsubmodel,unitab)
200 CALL hm_get_intv ('Seed' ,seed ,is_available,lsubmodel)
201 CALL hm_get_intv ('Idistri' ,i_method ,is_available,lsubmodel)
202card2
203 CALL hm_get_intv ('fail_ID' ,fail_id ,is_available,lsubmodel)
204 CALL hm_get_string('parameter' ,param ,ncharkey,is_available)
205c--------------------------------------------
206c Input check
207c--------------------------------------------
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
217c
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
227c
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
237c
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
259c-------------------------------
260c PRINT OUT
261c-------------------------------
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
274c-------------------------------
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
312c-------------------------------------------------------
313c Set up random seed
314c-------------------------------------------------------
315 IF( seed == 0 )THEN
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))
328 a_seed=seed
329 CALL random_seed(put=a_seed)
330 seed_random = 0
331 DEALLOCATE(a_seed)
332 ENDIF
333c-------------------------------------------------------
334c build uniform distributions
335c-------------------------------------------------------
336 distrib(1:50) = 0
337 ALLOCATE(array(nb_random+2))
338 CALL random_number(array) ! Uniform distribution
339c-------------------------------------------------------
340c build normal distribution
341c-------------------------------------------------------
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
364c
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
376c-------------------------------------------------------
377c check mean and standard deviation
378c-------------------------------------------------------
379 mean = sum(array)/nb_random
380 stdev = sqrt(sum((array - mean)**2)/nb_random)
381c-------------------------------------------------------
382c plot normal distribution
383c-------------------------------------------------------
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,'#')
401c
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)
412c
413 END IF ! biquad
414c------------------------
415 ENDDO ! NPERTURB_FAIL
416c
417 DEALLOCATE(tab_part)
418C-------------------------------------------------------------
419C-------------------------------------------------------------
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:')
438C-------------------------------------------------------------
439 2000 FORMAT(/
440 . ' generated mean VALUE . . . . . . .',1PG20.13/
441 . ' generated standard deviation . . .',1PG20.13)
442 2050 FORMAT(/
443 . ' generated mean VALUE . . . . . . .',1PG20.13)
444 2100 FORMAT(/
445 . ' generated seed VALUE . . . . . . .',I10/)
446C-------------------------------------------------------------
447 5000 FORMAT(/
448 . ' distribution of scale factors applied to ',A,' VALUE'/
449 . ' of failure criteria id= . . . . . .',I10)
450C------------------------------
451 RETURN
452 END
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)
subroutine hm_read_perturb_fail(mat_param, ipart, rnoise, ipartc, ipartg, ipartsp, igrpart, iparts, perturb, idperturb, index, index_ityp, npart_shell, offs, qp_iperturb, qp_rperturb, lsubmodel, unitab)
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
initmumps id
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)
Definition message.F:889
program starter
Definition starter.F:39