49
50
51
57 USE format_mod , ONLY : lfield
58
59
60
61#include "implicit_f.inc"
62
63
64
65#include "analyse_name.inc"
66
67
68
69#include "scr17_c.inc"
70#include "com01_c.inc"
71#include "com04_c.inc"
72#include "com10_c.inc"
73#include "sphcom.inc"
74#include "units_c.inc"
75#include "titr_c.inc"
76#include "param_c.inc"
77
78
79
80 TYPE (UNIT_TYPE_),INTENT(IN) :: UNITAB
81 INTEGER ISKN(LISKN,*), ITAB(*), ITABM1(*), NSN(*)
82 my_real skew(lskew,*), x(3,*), rtrans(ntransf,*)
83 TYPE(SUBMODEL_DATA) LSUBMODEL(*)
84 INTEGER NOM_OPT(LNOPT1,*)
85
86
87
88 INTEGER I, N, IMOV, J, N1, N2, N3, K, NSK,
89 . IUN, SUB_ID,
90 . IDSUB,ITY,L,READPT,J1,J2,NUMSPH_TMP,SUB_LEVEL,CUR_SUBMOD,
91 . IDIR,IFLAGUNIT,ID,UID,CPT
92 my_real p(12), pnor1, pnor2, pnorm1, det1, det2, det3, det, pp,bid,
93 . x0(3),rot(9)
94 CHARACTER(LEN=NCHARTITLE)::NOMSG
95 CHARACTER(LEN=NCHARTITLE) :: TITR
96 CHARACTER MESS*40
97 CHARACTER(LEN=NCHARKEY) :: KEY
98 CHARACTER(LEN=NCHARFIELD) :: DIR
99 LOGICAL IS_AVAILABLE
100
101
102
103 INTEGER USR2SYS
104 DATA iun/1/
105 DATA mess/'MOVING SKEW SYSTEM DEFINITION '/
106 DATA nomsg/'global skew system '/
107
108 idsub = 0
109
110 DO i=1,liskn
111 iskn(i,1)=0
112 ENDDO
113
114 DO i=1,lskew
115 skew(i,1)=zero
116 ENDDO
117 skew(1,1)=one
118 skew(5,1)=one
119 skew(9,1)=one
120 nom_opt(1,1)=0
121 readpt=0
122 CALL fretitl(nomsg,nom_opt(lnopt1-ltitr+1,1),ltitr)
123
124 IF(numskw==0)GOTO 201
125
126
127
129 i = 0
130
131
132
133 DO 100 cpt=1,numskw
134 i = i + 1
135
136
137
140 . unit_id = uid,
141 . submodel_id = sub_id,
142 . option_titr = titr,
143 . keyword2 = key)
144
146 CALL fretitl(titr,nom_opt(lnopt1-ltitr+1,i+1),ltitr)
147 imov = 0
148
149 iflagunit = 0
150 DO j=1,unitab%NUNITS
151 IF (unitab%UNIT_ID(j) == uid) THEN
152 iflagunit = 1
153 EXIT
154 ENDIF
155 ENDDO
156 IF (uid/=0.AND.iflagunit==0) THEN
157 CALL ancmsg(msgid=659,anmode=aninfo,msgtype=msgerror,
158 . i2=uid,i1=
id,c1=
'SKEW SYSTEM',c2=
'SKEW SYSTEM',
159 . c3=titr)
160 ENDIF
161
162 IF(key(1:3)=='FIX')THEN
163
164
165
166 CALL hm_get_floatv(
'globaloriginx',p(10),is_available,lsubmodel,unitab)
167 CALL hm_get_floatv(
'globaloriginy',p(11),is_available,lsubmodel,unitab)
168 CALL hm_get_floatv(
'globaloriginz',p(12),is_available,lsubmodel,unitab)
169
170 CALL hm_get_floatv(
'globalyaxisx',p(4),is_available,lsubmodel,unitab)
171 CALL hm_get_floatv(
'globalyaxisy',p(5),is_available,lsubmodel,unitab)
172 CALL hm_get_floatv(
'globalyaxisz',p(6),is_available,lsubmodel,unitab)
173
174 CALL hm_get_floatv(
'globalzaxisx',p(7),is_available,lsubmodel,unitab)
175 CALL hm_get_floatv(
'globalzaxisy',p(8),is_available,lsubmodel,unitab)
176 CALL hm_get_floatv(
'globalzaxisz',p(9),is_available,lsubmodel,unitab)
177
178 ELSEIF (key(1:4)=='MOV2')THEN
179 imov=2
180
181
182
183 CALL hm_get_intv(
'originnodeid',n1,is_available,lsubmodel)
184 CALL hm_get_intv(
'axisnodeid',n2,is_available,lsubmodel)
185 CALL hm_get_intv(
'planenodeid',n3,is_available,lsubmodel)
186
187 ELSE
188 imov=1
189 idir = 1
190
191
192
193 CALL hm_get_intv(
'originnodeid',n1,is_available,lsubmodel)
194 CALL hm_get_intv(
'axisnodeid',n2,is_available,lsubmodel)
195 CALL hm_get_intv(
'planenodeid',n3,is_available,lsubmodel)
196
197
198
200 DO k = 1,lfield
201 IF(dir(k:k) == 'X'.OR.dir(k:k) == 'x')idir = 1
202 IF(dir(k:k) == 'Y'.OR.dir(k:k) == 'y')idir = 2
203 IF(dir(k:k) == 'Z'.OR.dir(k:k) == 'z')idir = 3
204 ENDDO
205 iskn(6,i+1)=idir
206
207 ENDIF
209
210
211
212 IF (imov == 2) THEN
219 iskn(1,i+1)=n1
220 iskn(2,i+1)=n2
221 iskn(3,i+1)=n3
222 iskn(5,i+1)=imov
223 p(7)=x(1,n2)-x(1,n1)
224 p(8)=x(2,n2)-x(2,n1)
225 p(9)=x(3,n2)-x(3,n1)
226 p(1)=x(1,n3)-x(1,n1)
227 p(2)=x(2,n3)-x(2,n1)
228 p(3)=x(3,n3)-x(3,n1)
229
230
231
232 p(4)=p(8)*p(3)-p(9)*p(2)
233 p(5)=p(9)*p(1)-p(7)*p(3)
234 p(6)=p(7)*p(2)-p(8)*p(1)
235
236
237
238 p(1)=p(5)*p(9)-p(6)*p(8)
239 p(2)=p(6)*p(7)-p(4)*p(9)
240 p(3)=p(4)*p(8)-p(5)*p(7)
241
242
243
244 p(10) = x(1,n1)
245 p(11) = x(2,n1)
246 p(12) = x(3,n1)
247
248
249
250 pnor1=sqrt(p(7)*p(7)+p(8)*p(8)+p(9)*p(9))
251 IF (pnor1 < em20) THEN
253 . msgtype=msgerror,
254 . anmode=aninfo_blind_1,
255 . i2=itab(n1),
256 . i1=n,c1=titr,
257 . i3=itab(n2))
258 ENDIF
259
260 pnor2=sqrt(p(1)*p(1)+p(2)*p(2)+p(3)*p(3))
261 IF (pnor2 > em20) THEN
262 pnorm1=one/(pnor1*pnor2)
263 det1=abs((p(8)*p(3)-p(9)*p(2))*pnorm1)
264 det2=abs((p(9)*p(1)-p(7)*p(3))*pnorm1)
265 det3=abs((p(7)*p(2)-p(8)*p(1))*pnorm1)
266 det=
max(det1,det2,det3)
267 ELSE
268 det=zero
269 ENDIF
270 IF (det < em5) THEN
272 . msgtype=msgwarning,
273 . anmode=aninfo_blind_1,
275 IF(abs(p(2)) < em5) THEN
276 p(4)=abs(p(1))+ten
277 ELSE
278 p(5)=ten
279 ENDIF
280 ENDIF
281
282
283
284 ELSEIF (imov==1) THEN
289 iskn(1,i+1)=n1
290 iskn(2,i+1)=n2
291 iskn(5,i+1)=imov
292
293
294
295 IF(n2d==0)THEN
296
297 IF (idir == 1) THEN
298 p(1)=x(1,n2)-x(1,n1)
299 p(2)=x(2,n2)-x(2,n1)
300 p(3)=x(3,n2)-x(3,n1)
301 ELSEIF (idir == 2) THEN
302 p(4)=x(1,n2)-x(1,n1)
303 p(5)=x(2,n2)-x(2,n1)
304 p(6)=x(3,n2)-x(3,n1)
305 ELSEIF (idir == 3) THEN
306 p(7)=x(1,n2)-x(1,n1)
307 p(8)=x(2,n2)-x(2,n1)
308 p(9)=x(3,n2)-x(3,n1)
309 ENDIF
310
313 iskn(3,i+1)=n3
314
315 IF (idir == 1) THEN
316 p(4)=x(1,n3)-x(1,n1)
317 p(5)=x(2,n3)-x(2,n1)
318 p(6)=x(3,n3)-x(3,n1)
319 ELSEIF (idir == 2) THEN
320 p(7)=x(1,n3)-x(1,n1)
321 p(8)=x(2,n3)-x(2,n1)
322 p(9)=x(3,n3)-x(3,n1)
323 ELSEIF (idir == 3) THEN
324 p(1)=x(1,n3)-x(1,n1)
325 p(2)=x(2,n3)-x(2,n1)
326 p(3)=x(3,n3)-x(3,n1)
327 ENDIF
328
329 ELSE
330 p(1)=one
331 p(2)=zero
332 p(3)=zero
333
334 p(4)=x(1,n2)-x(1,n1)
335 p(5)=x(2,n2)-x(2,n1)
336 p(6)=x(3,n2)-x(3,n1)
337 ENDIF
338
339 p(10) = x(1,n1)
340 p(11) = x(2,n1)
341 p(12) = x(3,n1)
342
343
344
345 pnor1 = zero
346 IF (idir == 1) pnor1=sqrt(p(1)*p(1)+p(2)*p(2)+p(3)*p(3))
347 IF (idir == 2) pnor1=sqrt(p(4)*p(4)+p(5)*p(5)+p(6)*p(6))
348 IF (idir == 3) pnor1=sqrt(p(7)*p(7)+p(8)*p(8)+p(9)*p(9))
349 IF(pnor1<1.e-20) THEN
351 . msgtype=msgerror,
352 . anmode=aninfo_blind_1,
353 . i2=itab(n1),
355 . i3=itab(n2))
356 ENDIF
357
358 IF (idir == 1) pnor2=sqrt(p(4)*p(4)+p(5)*p(5)+p(6)*p(6))
359 IF (idir == 2) pnor2=sqrt(p(7)*p(7)+p(8)*p(8)+p(9)*p(9))
360 IF (idir == 3) pnor2=sqrt(p(1)*p(1)+p(2)*p(2)+p(3)*p(3))
361 IF(pnor2>1.e-20) THEN
362 pnorm1=1./(pnor1*pnor2)
363 IF (idir == 1) THEN
364 det1=abs((p(1)*p(5)-p(2)*p(4))*pnorm1)
365 det2=abs((p(1)*p(6)-p(3)*p(4))*pnorm1)
366 det3=abs((p(2)*p(6)-p(3)*p(5))*pnorm1)
367 ELSEIF (idir == 2) THEN
368 det1=abs((p(4)*p(8)-p(5)*p(7))*pnorm1)
369 det2=abs((p(4)*p(9)-p(6)*p(7))*pnorm1)
370 det3=abs((p(5)*p(9)-p(6)*p(8))*pnorm1)
371 ELSEIF (idir == 3) THEN
372 det1=abs((p(7)*p(2)-p(8)*p(1))*pnorm1)
373 det2=abs((p(7)*p(3)-p(9)*p(1))*pnorm1)
374 det3=abs((p(8)*p(3)-p(9)*p(2))*pnorm1)
375 ENDIF
376 det=
max(det1,det2,det3)
377 ELSE
378 det=zero
379 ENDIF
380 IF(det<em5) THEN
382 . msgtype=msgwarning,
383 . anmode=aninfo_blind_1,
385 IF (idir == 1) THEN
386 IF(abs(p(2))>em5) THEN
387 p(4)=abs(p(1))+ten
388 ELSE
389 p(5)=ten
390 ENDIF
391 ELSEIF (idir == 2) THEN
392 IF(abs(p(5))>em5) THEN
393 p(7)=abs(p(4))+ten
394 ELSE
395 p(8)=ten
396 ENDIF
397 ELSEIF (idir == 3) THEN
398 IF(abs(p(8))>em5) THEN
399 p(1)=abs(p(7))+ten
400 ELSE
401 p(2)=ten
402 ENDIF
403 ENDIF
404 ENDIF
405
406
407
408 IF (idir == 1) THEN
409 p(7)=p(2)*p(6)-p(3)*p(5)
410 p(8)=p(3)*p(4)-p(1)*p(6)
411 p(9)=p(1)*p(5)-p(2)*p(4)
412 ELSEIF (idir == 2) THEN
413 p(1)=p(5)*p(9)-p(6)*p(8)
414 p(2)=p(6)*p(7)-p(4)*p(9)
415 p(3)=p(4)*p(8)-p(5)*p(7)
416 ELSEIF (idir == 3) THEN
417 p(4)=p(8)*p(3)-p(9)*p(2)
418 p(5)=p(9)*p(1)-p(7)*p(3)
419 p(6)=p(7)*p(2)-p(8)*p(1)
420 ENDIF
421
422
423
424 IF (idir == 1) THEN
425 p(4)=p(8)*p(3)-p(9)*p(2)
426 p(5)=p(9)*p(1)-p(7)*p(3)
427 p(6)=p(7)*p(2)-p(8)*p(1)
428 ELSEIF (idir == 2) THEN
429 p(7)=p(2)*p(6)-p(3)*p(5)
430 p(8)=p(3)*p(4)-p(1)*p(6)
431 p(9)=p(1)*p(5)-p(2)*p(4)
432 ELSEIF (idir == 3) THEN
433 p(1)=p(5)*p(9)-p(6)*p(8)
434 p(2)=p(6)*p(7)-p(4)*p(9)
435 p(3)=p(4)*p(8)-p(5)*p(7)
436 ENDIF
437 ELSE
438
439
440
441
442 iskn(1,i+1)=0
443 iskn(2,i+1)=0
444 iskn(3,i+1)=0
445 iskn(5,i+1)=0
446
447 IF(p(4)==zero.AND.p(6)==zero) p(5)=sign(one,p(5))
448 IF(p(7)==zero.AND.p(8)==zero) p(9)=sign(one,p(9))
449
450
451
452 p(1)=p(5)*p(9)-p(6)*p(8)
453 p(2)=p(6)*p(7)-p(4)*p(9)
454 p(3)=p(4)*p(8)-p(5)*p(7)
455
456
457
458 p(4)=p(8)*p(3)-p(9)*p(2)
459 p(5)=p(9)*p(1)-p(7)*p(3)
460 p(6)=p(7)*p(2)-p(8)*p(1)
461
462 IF (sub_id /= 0)THEN
463 x0(1:3) = zero
465 IF (lsubmodel(j)%NOSUBMOD == sub_id) idsub = j
466 ENDDO
467 cur_submod = idsub
468 sub_level = lsubmodel(idsub)%LEVEL
469 DO WHILE (sub_level /= 0)
470 DO j=1,lsubmodel(cur_submod)%NBTRANS
471 ity = rtrans(lsubmodel(cur_submod)%IDTRANS(j),2)
472 IF( ity == 2 .OR. ity == 3 ) THEN
473 DO k=1,9
474 rot(k) = rtrans(lsubmodel(cur_submod)%IDTRANS(j),k+2)
475 ENDDO
479 ENDIF
480 ENDDO
481 sub_level = sub_level - 1
482 cur_submod = lsubmodel(cur_submod)%IFATHER
483 ENDDO
484 IF(lsubmodel(idsub)%NBTRANS /=0)
485 .
CALL subrotpoint(p(10),p(11),p(12),rtrans,sub_id,lsubmodel)
486 ENDIF
487 ENDIF
488
489
490
491 pp=sqrt(p(1)*p(1)+p(2)*p(2)+p(3)*p(3))
492 p(1)=p(1)/pp
493 p(2)=p(2)/pp
494 p(3)=p(3)/pp
495 pp=sqrt(p(4)*p(4)+p(5)*p(5)+p(6)*p(6))
496 p(4)=p(4)/pp
497 p(5)=p(5)/pp
498 p(6)=p(6)/pp
499 pp=sqrt(p(7)*p(7)+p(8)*p(8)+p(9)*p(9))
500 p(7)=p(7)/pp
501 p(8)=p(8)/pp
502 p(9)=p(9)/pp
503
504 DO k=1,12
505 skew(k,i+1)=p(k)
506 ENDDO
507
508 100 CONTINUE
509
510 WRITE (iout,'(A)')titre(85)
511 WRITE (iout,'(A)')titre(81)
512 DO 110 i=1,numskw
513 j=i+1
514 nsk = iskn(4,j)
515 n1=iskn(1,j)
516 n2=iskn(2,j)
517 n3=iskn(3,j)
518 IF(n1/=0)THEN
519 n1=itab(n1)
520 n2=itab(n2)
521 n3=itab(n3)
522 ENDIF
523 WRITE(iout,1000)
524 WRITE(iout,'(1X,4I10,1X,3F16.7,3F16.7)')nsk,n1,n2,n3,
525 & (skew(k,j),k=1,3),(skew(k,j),k=10,12)
526 WRITE(iout,'(3(42X,3F16.7/))') (skew(k,j),k=4,9)
527 110 CONTINUE
528
529 DO 140 k = 1,nsnod
530 nsn(k) = iabs(nsn(k))
531 140 CONTINUE
532
533 201 CONTINUE
534
535
536
537 IF(nspcond/=0)THEN
538 DO j=(numskw+1)+1,(numskw+1)+numsph
539 DO k=1,lskew
540 skew(k,j)=zero
541 ENDDO
542 skew(1,j)=one
543 skew(5,j)=one
544 skew(9,j)=one
545 iskn(1,j)=0
546 iskn(2,j)=0
547 iskn(3,j)=0
548
549 iskn(4,j)=-(j-numskw)
550 ENDDO
551 ENDIF
552
553
554
556 WRITE (iout,'(A)')titre(118)
557 WRITE (iout,'(A)')titre(119)
558 WRITE (iout,1001)
559 IF (nspcond/=0) THEN
560 j1 = (numskw+1)+numsph+1
562 numsph_tmp = numsph
563 ELSE
564 numsph_tmp = 0
565 j1 = (numskw+1)+1
567 ENDIF
568 DO j=j1,j2
569 DO k=1,lskew
570 skew(k,j)=zero
571 ENDDO
572 skew(1,j)=one
573 skew(5,j)=one
574 skew(9,j)=one
575 iskn(1,j)=0
576 iskn(2,j)=0
577 iskn(3,j)=0
578 iskn(5,j)=0
579 idsub = j-(numskw+numsph_tmp+1)
580 cur_submod = idsub
581 sub_level = lsubmodel(idsub)%LEVEL
582 x0(1:3) = zero
583
584 DO WHILE (sub_level /= 0)
585 DO k=1,lsubmodel(cur_submod)%NBTRANS
586 ity = rtrans(lsubmodel(cur_submod)%IDTRANS(k),2)
587 IF( ity == 2 .OR. ity == 3 )THEN
588 DO l=1,9
589 rot(l) = rtrans(lsubmodel(cur_submod)%IDTRANS(k),l+2)
590 ENDDO
594 ENDIF
595 ENDDO
596 sub_level = sub_level - 1
597 cur_submod = lsubmodel(cur_submod)%IFATHER
598 ENDDO
599 IF(lsubmodel(idsub)%NBTRANS /=0)
600 .
CALL subrotpoint(skew(10,j),skew(11,j),skew(12,j),
601 . rtrans,lsubmodel(idsub)%NOSUBMOD,lsubmodel)
602
603 iskn(4,j)=1000000001 + (j-numskw-numsph_tmp-2)
604 WRITE(iout,'(1X,I10,1X,3F16.7,3F16.7)')iskn(4,j),
605 . (skew(k,j),k=1,3),(skew(k,j),k=10,12)
606 WRITE(iout,'(3(12X,3F16.7/))') (skew(k,j),k=4,9)
607 ENDDO
608 ENDIF
609
610
611
612
613 IF (numskw/=0)
614 .
CALL udouble(iskn(4,1),liskn,
616 . mess,0,bid)
617
618 RETURN
619
620 1000 FORMAT(5x,'NUMBER',8x,'N1',8x,'N2',8x,'N3',10x,'VECTORS',42x,
621 . 'ORIGIN')
622 1001 FORMAT(5x,'NUMBER',10x,'VECTORS',42x,'ORIGIN')
623
624 RETURN
void anodset(int *id, int *type)
subroutine euler_vrot(x0, x, rot)
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 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)
integer function usr2sys(iu, itabm1, mess, id)
subroutine udouble(list, ilist, nlist, mess, ir, rlist)
subroutine subrotpoint(x, y, z, rtrans, sub_id, lsubmodel)