46
47
48
55 USE format_mod , ONLY : lfield
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 TYPE (UNIT_TYPE_),INTENT(IN) ::UNITAB
73 . rnoise(nperturb,numelc+numeltg+numels+numsph),
74 . qp_rperturb(nperturb,4)
75 INTEGER IPART(LIPART1,*),IPARTC(*),
76 . IPARTG(*),IPM(NPROPMI,*),
77 . PERTURB(NPERTURB),OFFS,
78 . IDPERTURB(NPERTURB),INDEX(NUMELC+NUMELTG+NUMELS+NUMSPH),
79 . INDEX_ITYP(NUMELC+NUMELTG+NUMELS+NUMSPH),NPART_SHELL,
80 . QP_IPERTURB(NPERTURB,6)
81 TYPE(SUBMODEL_DATA) LSUBMODEL(*)
82
83 TYPE (GROUP_) , DIMENSION(NGRPART) :: IGRPART
84
85
86
87 INTEGER I,J,K,NUMA,I_METHOD,MAX_PART,
88 . CPT_PART,NB_RANDOM,I_SEED,DISTRIB(50),
89 . II,NB_INTERV,IGRPRT,N,IOK,,SEED_RANDOM,
90 . ITYP,L,I_PERTURB_VAR,SIZEY,
91 CHARACTER(LEN=NCHARTITLE) :: TITR
92 CHARACTER(LEN=NCHARKEY) :: KEY
93 CHARACTER MES*40
94 INTEGER, DIMENSION(:,:), ALLOCATABLE :: TAB_PART
95 INTEGER, DIMENSION(:), ALLOCATABLE :: A_SEED
96 INTEGER, DIMENSION(1:8) :: DT_SEED
98 . mean,sd,mean_input,sd_input,max_distrib,temp,min_value,
99 . max_value,interv,VALUE,max_value1,minval,maxval,bid
100 my_real,
DIMENSION(:),
ALLOCATABLE :: array
101 CHARACTER*100 CHAR(100)
102 CHARACTER*100 CHAR1(100)
103 CHARACTER*100 CHAR2
104 CHARACTER MESS*40
105 CHARACTER(LEN=NCHARFIELD)::CHVAR
106 LOGICAL IS_AVAILABLE
107
108
109
110 DATA mess/'PERTURBATION DEFINITION '/
111
112
113
114
115
116 max_part = 0
117 ityp = 0
118 bid = 0
119 is_available = .false.
121
122
123
124 DO i=1+offs,npart_shell+offs
125
126
127 titr = ''
129 . option_id = idperturb(i),
130 . option_titr = titr)
131
132
133 ityp = 1
134
135 i_perturb_var = 0
136 cpt_part = 0
137 igrprt = 0
138 iok = 0
139 chvar = ''
140
141
142 CALL hm_get_intv(
'grpart_ID' ,igrprt ,is_available,lsubmodel)
144 IF (chvar(1:5) == 'thick' .OR. chvar(1:5) == 'THICK') i_perturb_var = 1
145
146 empty = 1
147 DO k = 1,lfield
148 IF(chvar(k:k) /= ' ') empty = 0
149 ENDDO
150
151
152 IF (i_perturb_var /= 1 .AND. empty == 0)
CALL ancmsg(msgid=1194,
153 . msgtype=msgerror,
154 . anmode=aninfo,
155 . i1=idperturb(i),
156 . c1=titr,
157 . c2=chvar)
158
159
160 IF (igrprt /= 0) THEN
161 DO n=1,ngrpart
162 IF (igrpart(n)%ID == igrprt) THEN
163 igrprt=n
164 iok = 1
165 ityp = 1
166 EXIT
167 END IF
168 END DO
169 ENDIF
170
171
172 perturb(i) = ityp
173
174
175 IF (iok == 0) THEN
177 . msgtype=msgerror,
178 . anmode=aninfo,
179 . i1=idperturb(i),
180 . c1=titr,
181 . i2=igrprt,
182 . c2='GROUP OF PART')
183 ELSEIF(iok == 1)THEN
184 cpt_part = igrpart(igrprt)%NENTITY
185 ENDIF
186 max_part =
max(max_part,cpt_part)
187 ENDDO
188
189 ALLOCATE(tab_part(nperturb,max_part))
190
191
192
193
195 DO i = 1+offs,npart_shell+offs
196
197
198 index(1:(numelc+numeltg+numels+numsph)) = 0
199 index_ityp(1:numelc+numeltg+numels+numsph) = 0
200
201
202 titr = ''
204 . option_id = idperturb(i),
205 . option_titr = titr)
206
207
208 ityp = 1
209
210
211 CALL hm_get_floatv(
'F_Mean' ,mean ,is_available, lsubmodel, unitab)
212 CALL hm_get_floatv(
'Deviation' ,sd ,is_available, lsubmodel, unitab)
213 CALL hm_get_floatv(
'Min_cut' ,minval ,is_available, lsubmodel, unitab)
214 CALL hm_get_floatv(
'Max_cut' ,maxval ,is_available, lsubmodel, unitab)
216 CALL hm_get_intv(
'Idistri' ,i_method,is_available, lsubmodel)
217
218
219 IF(i_method == 0) i_method = 2
220 IF(minval == zero .AND. maxval == zero) THEN
221 IF(i_method == 1) THEN
222 ELSEIF(i_method == 2)THEN
223 minval = -ep30
224 maxval = ep30
225 ENDIF
226 ENDIF
227 sd_input = sd
228 mean_input = mean
229
230
231 qp_iperturb(i,1) = idperturb(i
232 qp_iperturb(i,2) = ityp
233 qp_iperturb(i,3) =
seed
234 qp_iperturb(i,4) = i_method
235 qp_rperturb(i,1) = mean
236 qp_rperturb(i,2) = sd
237 qp_rperturb(i,3) = minval
238 qp_rperturb(i,4) = maxval
239
240
241 cpt_part = 0
242 igrprt = 0
243 iok = 0
244
245
246 CALL hm_get_intv(
'grpart_ID',igrprt,is_available,lsubmodel)
247 qp_iperturb(i,5) = igrprt
249 IF (chvar(1:5) == 'thick' .OR. chvar(1:5) == 'THICK') qp_iperturb(i,6) = 1
250
251
252 IF (igrprt /= 0) THEN
253 DO n=1,ngrpart
254 IF (igrpart(n)%ID == igrprt) THEN
255 igrprt = n
256 iok = 1
257 ityp = 1
258 EXIT
259 END IF
260 END DO
261 ENDIF
262
263
264 IF (iok == 1) THEN
265 cpt_part = 0
266 DO j=1,igrpart(igrprt)%NENTITY
267 cpt_part = cpt_part + 1
268 numa = igrpart(igrprt)%ENTITY(j)
269 tab_part(i,cpt_part) = numa
270 END DO
271 ENDIF
272
273
274 IF(i_method == 2) THEN
275 WRITE (iout,1000)
276 . idperturb(i),
'GAUSSIAN',mean_input,sd_input,
seed
277 WRITE (iout,'(10I10)') ipart(4,tab_part(i,1:cpt_part))
278 WRITE(iout,*) ' '
279 WRITE(iout,*) ' '
280 ELSEIF(i_method == 1) THEN
281 WRITE (iout,1100)
282 . idperturb(i),
'RANDOM',
seed
283 WRITE (iout,'(10I10)') ipart(4,tab_part(i,1:cpt_part))
284 WRITE(iout,*) ' '
285 WRITE(iout,*) ' '
286 ENDIF
287
288
289 nb_random = 0
290 DO ii=1,numelc
291 DO k=1,cpt_part
292 IF (ipartc(ii) == tab_part(i,k)) THEN
293 nb_random = nb_random + 1
294 index(nb_random) = ii
295 index_ityp(nb_random) = 3
296 ENDIF
297 ENDDO
298 ENDDO
299 DO ii=1,numeltg
300 DO k=1,cpt_part
301 IF(ipartg(ii) == tab_part(i,k)) THEN
302 nb_random = nb_random + 1
303 index(nb_random) = ii
304 index_ityp(nb_random) = 7
305 ENDIF
306 ENDDO
307 ENDDO
308
309
311 CALL random_seed(size=i_seed)
312 ALLOCATE(a_seed(1:i_seed))
313 CALL random_seed(get=a_seed)
314 CALL date_and_time(values=dt_seed)
315 a_seed=dt_seed(8)*dt_seed(7)*dt_seed(6)
316 seed=dt_seed(8)*dt_seed(7)*dt_seed(6)
317 CALL random_seed(put=a_seed)
318 seed_random = 1
319 DEALLOCATE(a_seed)
320 ELSE
321 CALL random_seed(size=i_seed)
322 ALLOCATE(a_seed(1:i_seed))
324 CALL random_seed(put=a_seed)
325 seed_random = 0
326 DEALLOCATE(a_seed)
327 ENDIF
328
329
330 char=''
331 char1=''
332 char2=''
333 distrib(1:50) = 0
334 ALLOCATE(array(nb_random+2))
335 CALL random_number(array)
336
337
338 max_value = -ep30
339 min_value = ep30
340 IF ( i_method == 2) THEN
341 DO ii = 1, nb_random+1, 2
342 temp = sd * sqrt(-2.0*log(array(ii))) * cos(2*pi*array(ii+1)) +
343 . mean
344 array(ii+1) =
345 . sd * sqrt(-2.0*log(array(ii))) * sin(2*pi*array(ii+1)) + mean
346 array(ii) = temp
347 END DO
348 DO ii = 1, nb_random
349 array(ii) =
max(
min(maxval,array(ii)),minval)
350 max_value =
max(array(ii),max_value)
351 min_value =
min(array(ii),min_value)
352 END DO
353 ELSEIF(i_method == 1)THEN
354 DO ii = 1, nb_random
355 array(ii) = array(ii)*(maxval-minval)+minval
356 max_value =
max(array(ii),max_value)
357 min_value =
min(array(ii),min_value)
358 END DO
359 ENDIF
360
361
362 DO ii = 1, nb_random
363 IF (index_ityp(ii) == 3) THEN
364 rnoise(i,index(ii)) = array(ii)
365 ELSEIF (index_ityp(ii) == 7) THEN
366 rnoise(i,index(ii)+numelc) = array(ii)
367 ENDIF
368 ENDDO
369
370
371 mean = sum(array)/nb_random
372 sd = sqrt(sum((array - mean)**2)/nb_random)
373
374
375 IF(i_method == 2) THEN
376 max_distrib = one /(sd*sqrt(two * pi))
377 ELSEIF(i_method == 1) THEN
378 max_distrib = one /(max_value-min_value)
379 ENDIF
380 WRITE (iout,3000)
381 WRITE(iout,*) ' '
382 nb_interv = 50
383 sizey = 20
384 IF (minval /= -ep30 .AND. maxval /= ep30)THEN
385 min_value = minval
386 max_value = maxval
387 ENDIF
388 CALL plot_distrib( array,nb_random, nb_interv,sizey,min_value,
389 . max_value,max_distrib,'#')
390 IF(i_method == 2) THEN
391 WRITE (iout,2000) mean,sd
392 ELSEIF(i_method == 1) THEN
393 WRITE (iout,2050) mean
394 ENDIF
395 IF(seed_random == 1)
WRITE (iout,2100)
seed
396 WRITE(iout,*) ' '
397 WRITE(iout,*) ' '
398 IF (ALLOCATED(array)) DEALLOCATE(array)
399 ENDDO
400 DEALLOCATE(tab_part)
401
402 1000 FORMAT(/' PERTURBATION ID',i10/
403 + ' ---------------'/
404 + ' TYPE . . . . . . . . . . . . . . .',a/
405 + ' INPUT MEAN VALUE . . . . . . . . .',1pg20.13/
406 + ' INPUT STANDARD DEVIATION . . . . .',1pg20.13/
407 + ' INPUT SEED VALUE . . . . . . . . .',i10/
408 + ' SHELL THICKNESSES, PARTS:')
409 1100 FORMAT(/' PERTURBATION ID',i10/
410 + ' ---------------'/
411 + ' TYPE . . . . . . . . . . . . . . .',a/
412 + ' INPUT SEED VALUE . . . . . . . . .',i10/
413 + ' SHELL THICKNESSES, PARTS:')
414
415 2000 FORMAT(/
416 + ' GENERATED MEAN VALUE . . . . . . .',1pg20.13/
417 + ' GENERATED STANDARD DEVIATION . . .',1pg20.13)
418 2050 FORMAT(/
419 + ' GENERATED MEAN VALUE . . . . . . .',1pg20.13)
420 2100 FORMAT(/
421 + ' GENERATED SEED VALUE . . . . . . .',i10/)
422
423 3000 FORMAT(/
424 + ' DISTRIBUTION OF SCALE FACTORS APPLIED TO THICKNESSES OF SHELLS')
425
426
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_start(entity_type)
integer, parameter nchartitle
integer, parameter ncharkey
integer, parameter ncharfield
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)