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,minval,maxval
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.AND. IF (MINVAL == ZERO 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_option_count(entity_type, hm_option_number)
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)