Télécharger la plaquette Cast3M


fr en

Télécharger hhopre.eso

Retour à la liste

Numérotation des lignes : 

hhopre

  1. C HHOPRE    SOURCE    OF166741  24/06/19    21:15:09     11942          
  2.  
  3.       SUBROUTINE HHOPRE (charHHO, mailHHO, lentHHO, iret)
  4.  
  5.       IMPLICIT INTEGER(I-N)
  6.       IMPLICIT REAL*8(A-H,O-Z)
  7.  
  8. -INC PPARAM
  9. -INC CCOPTIO
  10. -INC CCGEOME
  11.  
  12. -INC CCHHOPA
  13. -INC CCHHOPR
  14.  
  15. -INC SMCOORD
  16. -INC SMELEME
  17. -INC SMLENTI
  18.       POINTEUR mleCEL.mlenti,mleSQE.mlenti
  19.  
  20.       EXTERNAL LONG
  21.  
  22.       CHARACTER*(*) charHHO
  23.       CHARACTER*(3) hyp_z
  24.  
  25.       iret = 0
  26.  
  27. C=----------------------------------------------------------------------
  28. C= ORDRE DES CELLULES ET DES FACES :
  29. C=----------------------------------------------------------------------
  30. C= Decodage de la chaine pour avoir ordre des cellules et des faces :
  31.       n_z = LONG(charHHO)
  32.       l_z = LEN(charHHO)
  33. c-AV      IF (l_z.LT.13) THEN
  34.       IF (l_z.LT.11) THEN
  35.         write(ioimp,*) 'LENGTH charHHO incorrect (PBM 10)'
  36.         iret = 5
  37.         RETURN
  38.       END IF
  39. C= Passage en minuscules :
  40.       CALL CHCASS(charHHO(1:n_z),0,charHHO(1:n_z))
  41. C= Transformation des " " en "_" --> "hho_c_f" ou "hho_c_f_hyp" ?
  42.       DO i = 1, n_z
  43.         IF (charHHO(i:i).EQ.' ') charHHO(i:i) = '_'
  44.       END DO
  45. c-dbg      write(ioimp,*) 'charHHO =>'//charHHO(1:n_z)//'<=',n_z,l_z
  46. C= La chaine charHHO doit etre de la forme :
  47. C= "hho_c_f" ou "hho_c_f_hyp" avec c et f entiers positifs, et hyp=sp/ft
  48. C= Petits tests :
  49.       IF (n_z.LT.7) THEN
  50.         write(ioimp,*) 'String HHO too short (PBM 0)'
  51.         iret = 5
  52.         RETURN
  53.       END IF
  54.       IF (charHHO(1:3).NE.'hho') THEN
  55.         write(ioimp,*) 'String HHO incorrect (PBM 1)'
  56.         iret = 21
  57.         RETURN
  58.       END IF
  59.       i_c = INDEX(charHHO(1:n_z),'_')
  60.       IF (i_c.NE.4) THEN
  61.         write(ioimp,*) 'String HHO incorrect (PBM 2)',i_c
  62.         iret = 21
  63.         RETURN
  64.       END IF
  65.       i_f = INDEX(charHHO(i_c+1:n_z),'_')
  66.       IF (i_f.LT.1) THEN
  67.         write(ioimp,*) 'Sring HHO incorrect (PBM 3)',i_f
  68.         iret = 21
  69.         RETURN
  70.       END IF
  71.       i_f = i_f + i_c
  72.       i_h = INDEX(charHHO(i_f+1:n_z),'_')
  73.       IF (i_h.LT.1) THEN
  74.         i_h = n_z + 1
  75.         i_3 = 0
  76.       ELSE
  77.         i_h = i_f + i_h
  78.         i_3 = n_z
  79.       ENDIF
  80. c-dbg      write(ioimp,*) 'HHOPRE 1 ',i_c,i_f,i_h,i_3
  81.  
  82.       i_z = (i_f - 1) - (i_c + 1) + 1
  83.       IF (i_z.LT.1) THEN
  84.         write(ioimp,*) 'CELL ORDER undefined (PBM 4)'
  85.         iret = 21
  86.         RETURN
  87.       END IF
  88.       i_z = (i_h - 1) - (i_f + 1) + 1
  89.       IF (i_z.LT.1) THEN
  90.         write(ioimp,*) 'FACE ORDER undefined (PBM 5)'
  91.         iret = 21
  92.         RETURN
  93.       END IF
  94.  
  95.       n_o_cell = -999
  96.       n_o_face = -999
  97.       READ(charHHO(i_c+1:i_f-1),*,ERR=901) n_o_cell
  98.  901  CONTINUE
  99.       READ(charHHO(i_f+1:i_h-1),*,ERR=902) n_o_face
  100.  902  CONTINUE
  101.       IF (n_o_cell.LT.0 .OR. n_o_face.LT.0) THEN
  102.         write(ioimp,*) 'CELL/FACE ORDER incorrect (PBM 6)'
  103.         write(ioimp,*) '   =>',charHHO(i_c+1:i_f-1),'<=',
  104.      &                 '   =>',charHHO(i_f+1:i_h-1),'<='
  105.         iret = 21
  106.         RETURN
  107.       END IF
  108. C= Test de coherence :
  109.       IF (IDIM.EQ.1) THEN
  110.         IF (n_o_face.NE.0) THEN
  111.           write(ioimp,*) 'IDIM = 1 : FACE ORDER must be 0 (PBM 7)'
  112.           iret = 21
  113.         END IF
  114.       END IF
  115.       IF ( (IDIM.EQ.2) .OR. (IDIM.EQ.3) ) THEN
  116.         IF ( n_o_cell.LT.(n_o_face-1) .OR.
  117.      &       n_o_cell.GT.(n_o_face+1) ) THEN
  118.           write(ioimp,*) 'IDIM = 2/3 : CELL ORDER uncorrect (PBM 8)'
  119.           write(ioimp,*) '  l cell in [ k face - 1 ; k face + 1]'
  120.           iret = 21
  121.         END IF
  122.       END IF
  123. C= Validation de l'hypothese de calcul ("sp" ou "ft")
  124.       IF (i_3.GT.0) THEN
  125.         hyp_z = '   '
  126.         hyp_z = charHHO(i_h:n_z)
  127.         IF (hyp_z.NE.'_ft' .AND. hyp_z.NE.'_sp') THEN
  128.           write(ioimp,*) 'HYPOTHESIS SP/FT incorrect (PBM 9)'
  129.           iret = 21
  130.           RETURN
  131.         END IF
  132.       ELSE
  133.         hyp_z = '_sp'
  134.       END IF
  135. *AV      charHHO(1:13) = 'hho_00_00_** '
  136. *AV      WRITE(charHHO(5:6),'(I2.2)') n_o_cell
  137. *AV      WRITE(charHHO(8:9),'(I2.2)') n_o_face
  138. *AV      charHHO(10:12) = hyp_z
  139. *AVC= La chaine charHHO est de la forme : "hho_cc_ff_hy ".
  140.       charHHO(1:11) = 'hho_0_0_** '
  141.       WRITE(charHHO(5:5),'(I1.1)') n_o_cell
  142.       WRITE(charHHO(7:7),'(I1.1)') n_o_face
  143.       charHHO(8:10) = hyp_z
  144. C= La chaine charHHO est de la forme : "hho_c_f_hy ".
  145. C= Elle sera utile pour appeler les fonctions HHO adequates.
  146.  
  147. C= Nombre de ddl par face et par cellule selon ordre et dime
  148. C= Attention pour les faces : dime = idim - 1 !
  149. C= nddl_dir = Produit(i = 1 a dime) [ (ordre + dime + 1 - i) / i ]
  150.       IF (IDIM.EQ.1) THEN
  151.         n_d_face = 1
  152.         n_d_cell = n_o_cell + 1
  153.       ELSE IF (IDIM.EQ.2) THEN
  154.         n_d_face = n_o_face + 1
  155.         n_d_cell = (n_o_cell + 2) * (n_o_cell + 1) / 2
  156.       ELSE IF (IDIM.EQ.3) THEN
  157.         n_d_face = (n_o_face + 2) * (n_o_face + 1) / 2
  158.         n_d_cell = (n_o_face + 3) * (n_o_face + 2) * (n_o_face + 1) / 6
  159.       END IF
  160. C= nddl tot = nddl_dir * idfo (idfo = idim pour MECA, = 1 pour THER ou PFM)
  161.  
  162. C= Quelques restrictions qui pourront etre levees par la suite :
  163.       IF (n_o_face.NE.1) THEN
  164.         write(ioimp,*) 'FACE ORDER must be equal to 1 for the moment'
  165.         iret = 21
  166.         RETURN
  167.       END IF
  168.       IF (n_d_face.GT.10) THEN
  169.         write(ioimp,*) 'FACE ORDER is too big (PBM 12F)'
  170.         iret = 21
  171.         RETURN
  172.       END IF
  173.       IF (n_d_cell.GT.10) THEN
  174.         write(ioimp,*) 'CELL ORDER is too big (PBM 12C)'
  175.         iret = 21
  176.         RETURN
  177.       END IF
  178.  
  179. c-dbg      write(ioimp,*)
  180. c-dbg      write(ioimp,*) 'HHO : IDIM, IFOUR =', IDIM,IFOUR
  181. c-dbg      write(ioimp,*) '      =>'//charHHO(1:LONG(charHHO))//'<='
  182. c-dbg      write(ioimp,*) '      FACE ORDER / DOF DIR = ',n_o_face,n_d_face
  183. c-dbg      write(ioimp,*) '      CELL ORDER / DOF DIR = ',n_o_cell,n_d_cell
  184. c-dbg      write(ioimp,*) '      HYPOTHESIS = ', hyp_z(2:3)
  185.  
  186. C=----------------------------------------------------------------------
  187. C= Traitement du maillage total du modele (fourni en entree) :
  188. C=----------------------------------------------------------------------
  189.       meleme  = mailHHO
  190. c*      segact,meleme*nomod (segment actif en entree)
  191.  
  192. C=----------------------------------------------------------------------
  193. C= Construction du maillage des faces :
  194. C=----------------------------------------------------------------------
  195.       CALL ECROBJ('MAILLAGE',mailHHO)
  196.       IF (IDIM.EQ.2) THEN
  197.         CALL CHANLG
  198.       ELSE IF (IDIM.EQ.3) THEN
  199.         CALL ECRCHA('NOID')
  200.         CALL ENVVO2(1)
  201.       ELSE
  202.         IERR = 5
  203.       END IF
  204.       IF (IERR.NE.0) THEN
  205.         iret = 21
  206.         RETURN
  207.       END IF
  208.       CALL LIROBJ('MAILLAGE',mailSQE,1,iretc)
  209.  
  210. C= On reactive les maillages :
  211.       CALL ACTOBJ('MAILLAGE',mailHHO,1)
  212.       CALL ACTOBJ('MAILLAGE',mailSQE,1)
  213.  
  214. C= ---------------------------
  215. C= Compatibilite des maillages avec la formulation HHO
  216. C= ---------------------------
  217. C= Les elements geometriques (cellules/faces) autorises :
  218.       IF (IDIM.EQ.3) THEN
  219.         indc = IC3MAX+1
  220.         lonc = NC3MAX
  221.         indf = IF3MAX+1
  222.         lonf = NF3MAX
  223.       ELSE IF (IDIM.EQ.2) THEN
  224.         indc = IC2MAX+1
  225.         lonc = NC2MAX
  226.         indf = IF2MAX+1
  227.         lonf = NF2MAX
  228. c-      ELSE IF (IDIM.EQ.1) THEN
  229.       ELSE
  230.         indc = IC1MAX+1
  231.         lonc = NC1MAX
  232.         indf = IF1MAX+1
  233.         lonf = NF1MAX
  234.       END IF
  235.  
  236. C= Verifications du maillage :
  237. C= ---------------------------
  238.       meleme  = mailHHO
  239.       nbsCEL = meleme.LISOUS(/1)
  240.       ipt1 = meleme
  241.       nbs1 = MAX(1,nbsCEL)
  242.       n_z = 0
  243.       DO i = 1, nbs1
  244.         IF (nbsCEL.NE.0) ipt1 = meleme.LISOUS(i)
  245.         ity1 = ipt1.itypel
  246.         if (ity1.eq.32) ity1 = ity1 * 100 + ipt1.num(/1)
  247.         i_z = 0
  248.         CALL PLACE2(LICHHO(indc),lonc,i_z,ity1)
  249.         IF (i_z.EQ.0) THEN
  250.           write(ioimp,*) 'HHOPRE: cell not defined',ipt1,NOMS(ity1),ity1
  251.           n_z = n_z + 1
  252.         END IF
  253.       END DO
  254.       IF (n_z.GT.0) THEN
  255.         iret = 251
  256.         RETURN
  257.       END IF
  258.  
  259.       meleme = mailSQE
  260.       nbsSQE = meleme.lisous(/1)
  261.       if (nbsSQE.ne.0) then
  262.         IF (IDIM.EQ.1 .OR. IDIM.EQ.2) THEN
  263.           write(ioimp,*) 'HHO 1D/2D : skeleton not simple'
  264.           iret = 5
  265.           return
  266.         END IF
  267.       end if
  268.       ipt1 = meleme
  269.       nbs1 = MAX(1,nbsSQE)
  270.       n_z = 0
  271.       DO i = 1, nbs1
  272.         IF (nbsSQE.NE.0) ipt1 = meleme.lisous(i)
  273.         ity1 = ipt1.itypel
  274. c*face3Dpoly : inconnu a ce jour       if (ity1.eq.  ) ity1 = formule a ecrire si necessaire
  275.         i_z = 0
  276.         CALL PLACE2(LIFHHO(indf),lonf,i_z,ity1)
  277.         IF (i_z.EQ.0) THEN
  278.           write(ioimp,*) 'HHOPRE: face not defined',ipt1,NOMS(ity1),ity1
  279.           n_z = n_z + 1
  280.         END IF
  281.       END DO
  282.       IF (n_z.GT.0) THEN
  283.         iret = 5
  284.         RETURN
  285.       END IF
  286.  
  287. C=----------------------------------------------------------------------
  288. C= Initialisations/Verifications DIMENSION et MODE DE CALCUL HHO
  289. C=----------------------------------------------------------------------
  290.       i1_HHO = 0
  291.       IF (IDIHHO.LT.0) THEN
  292.         i1_HHO = 1
  293.         IDIHHO = IDIM
  294.         IFOHHO = IFOUR
  295.       END IF
  296.       IF (IDIHHO .NE. IDIM) THEN
  297.         write(ioimp,*) 'HHOPRE: IDIM cannot be changed'
  298.         iret = 5
  299.         RETURN
  300.       END IF
  301.       IF (IFOHHO .NE. IFOUR) THEN
  302.         write(ioimp,*) 'HHOPRE: IFOUR cannot be changed'
  303.         iret = 5
  304.         RETURN
  305.       END IF
  306.  
  307. C- Juste pour le debogage :
  308.       if (nbchho(0).ne.0) then
  309.         write(ioimp,*) 'HHOPRE: nbchho(0) not 0'
  310.         iret = 5
  311.       end if
  312.       if (nbfhho(0).ne.0) then
  313.         write(ioimp,*) 'HHOPRE: nbfhho(0) not 0'
  314.         iret = 5
  315.       end if
  316.       if (iret.ne.0) return
  317.  
  318. c-dbg c Si on souhaite surveiller un maillage
  319. c-dbg      CALL OOOSUR(M..HHO)
  320. c-dbg      msurve = M...HO
  321.  
  322. C= Construction du maillage HHO global MCEHHO :
  323. C= --------------------------------------------
  324.       JG = 2 * NCEMAX
  325.       SEGINI,mleCEL
  326.       DO i = 1, NCEMAX
  327.         mleCEL.lect(i) = NBCHHO(i)
  328.         mleCEL.lect(i+NCEMAX) = 0
  329.       END DO
  330.  
  331. C= Remplissage initial de MCEHHO :
  332.       IF (MCEHHO.LT.0) THEN
  333. c-dbg        write(ioimp,*) 'HHOPRE: Initialisation MCEHHO'
  334. C= Petites verifications juste pour le debogage :
  335.         if (i1_HHO.NE.1) then
  336.           write(ioimp,*) 'HHOPRE-MCEHHO: Bad initialization (1)'
  337.           iret = 5
  338.           return
  339.         end if
  340.         n_z = 0
  341.         DO i = 1, NCEMAX
  342.           IF (MACHHO(i).GT.0 .OR. NBCHHO(i).NE.0) THEN
  343.             n_z = n_z + 1
  344.             write(ioimp,*) 'HHOPRE-MCEHHO: Init.',i,MACHHO(i),NBCHHO(i)
  345.           END IF
  346.         END DO
  347.         if (n_z.gt.0) then
  348.           iret = 5
  349.           return
  350.         end if
  351.  
  352. C- Si plusieurs zones, on duplique entete du maillage
  353.         IF (nbsCEL.GT.1) THEN
  354.           ipt1 = mailHHO
  355.           SEGINI,meleme=ipt1
  356.           MCEHHO = meleme
  357.         ELSE
  358.           MCEHHO = mailHHO
  359.         END IF
  360.  
  361.         meleme = MCEHHO
  362.         ipt1 = meleme
  363.         nbs1 = MAX(1,nbsCEL)
  364.         DO i = 1, nbs1
  365.           IF (nbsCEL.NE.0) ipt1 = meleme.lisous(i)
  366.           ity1 = ipt1.itypel
  367.           if (ity1.eq.32) ity1 = ity1 * 100 + ipt1.num(/1)
  368.           i_z = 0
  369.           CALL PLACE2(LICHHO(indc),lonc,i_z,ity1)
  370.           nbe1 = ipt1.num(/2)
  371. c* a revoir quand il y aura cumul
  372.           NBCHHO(i_z) = nbe1
  373.           MACHHO(i_z) = ipt1
  374.           if (mleCEL.lect(i_z+NCEMAX).ne.0)
  375.      &      write(ioimp,*) 'HHOPRE: bizarre',i,i_z,ipt1,ity1
  376.           mleCEL.lect(i_z+NCEMAX) = nbe1
  377.         END DO
  378.  
  379. C= Mise a jour de MCEHHO :
  380.       ELSE
  381.         if (i1_HHO.NE.0) then
  382.           write(ioimp,*) 'HHOPRE-MCEHHO: Bad initialization (2)'
  383.           iret = 5
  384.           return
  385.         end if
  386.         write(ioimp,*) 'MCEHHO deja defini --> A completer'
  387.         write(ioimp,*) 'MAIS CAS EN COURS D IMPLEMENTATION'
  388.         iret = 5
  389.         RETURN
  390.       END IF
  391.  
  392.       nelCEL = 0
  393.       nbsCEL = 0
  394.       DO i = 1, NCEMAX
  395.         nelCEL = nelCEL + NBCHHO(i)
  396.         IF (MACHHO(i).GT.0) nbsCEL = nbsCEL + 1
  397.       END DO
  398.       NCEHHO = nelCEL
  399.       NUCHHO = nbsCEL
  400. c-dbg      write(ioimp,*) 'HHOPRE-MCEHHO:',NCEHHO,NUCHHO
  401.  
  402. C* On reordonne MCEHHO selon liste type :
  403.       ipt2 = MCEHHO
  404.       IF (nbsCEL.GT.1) THEN
  405.         segact,ipt2*MOD
  406.         isou = 0
  407.         DO i = 1, NCEMAX
  408.           ipt1 = MACHHO(i)
  409.           IF (ipt1.GT.0) THEN
  410.             isou = isou + 1
  411.             ipt2.lisous(isou) = ipt1
  412.           END IF
  413.         END DO
  414.         segact,ipt2*NOMOD
  415.       END IF
  416. C= Doit-on faire un savseg de MCEHHO et des sous-zones eventuelles ?
  417.       DO i = 1, NCEMAX
  418.         NBCHHO(i) = NBCHHO(i) + NBCHHO(i-1)
  419.       END DO
  420.       if (nbchho(ncemax).ne.NCEHHO)
  421.      &  write(ioimp,*) 'Bizarre nbchho(ncemax) != NCEHHO'
  422.  
  423. C= Construction du squelette HHO global MSQHHO :
  424. C= ---------------------------------------------
  425. c== Remplissage de ...HHO : l'indice i correspond aux donnees des faces
  426. C== qui sont des polygones a i cotes (indice =1 fixe a 0, =2 pour 2D, 
  427. C== =3 et superieurs pour 3D, limite = faces a moins de hho_max_edge cotes)
  428.       IF (IDIM.EQ.3) THEN
  429.         idebf = 3
  430.         ifinf = NFAMAX
  431.       ELSE IF (IDIM.EQ.2) THEN
  432.         idebf = 2
  433.         ifinf = 2
  434.       ELSE IF (IDIM.EQ.1) THEN
  435.         idebf = 1
  436.         ifinf = 1
  437.       END IF
  438.  
  439.       JG = 2 * NFAMAX
  440.       SEGINI,mleSQE
  441.       DO i = 1, NFAMAX
  442.         mleSQE.lect(i) = NBFHHO(i)
  443.         mleSQE.lect(i+NFAMAX) = 0
  444.       END DO
  445.  
  446. C= Remplissage initial de MSQHHO :
  447.       IF (MSQHHO.LT.0) THEN
  448. c-dbg        write(ioimp,*) 'HHOPRE: Initialisation MSQHHO'
  449. C= Petites verifications juste pour le debogage :
  450.         if (i1_HHO.NE.1) then
  451.           write(ioimp,*) 'HHOPRE-MSQHHO: Bad initialization (1)'
  452.           iret = 5
  453.           return
  454.         end if
  455.         n_z = 0
  456.         DO i = 1, NFAMAX
  457.           IF (MAFHHO(i).GT.0 .OR. NBFHHO(i).NE.0) THEN
  458.             n_z = n_z + 1
  459.             write(ioimp,*) 'HHOPRE-MSQHHO: Init.',i,MAFHHO(i),NBFHHO(i)
  460.           END IF
  461.         END DO
  462.         if (n_z.gt.0) then
  463.           iret = 5
  464.           return
  465.         end if
  466.  
  467. C- Si plusieurs zones, on duplique entete du maillage
  468.         IF (nbsSQE.GT.1) THEN
  469.           ipt1 = mailSQE
  470.           SEGINI,meleme=ipt1
  471.           MSQHHO = meleme
  472.         ELSE
  473.           MSQHHO = mailSQE
  474.         END IF
  475.  
  476.         meleme = mailSQE
  477.         ipt1 = meleme
  478.         nbs1 = MAX(1,nbsSQE)
  479.         DO i = 1, nbs1
  480.           IF (nbsSQE.NE.0) ipt1 = meleme.lisous(i)
  481.           nbn1 = ipt1.num(/1)
  482.           nbe1 = ipt1.num(/2)
  483.           if (mleSQE.lect(NFAMAX+nbn1).ne.0)
  484.      &      write(ioimp,*) 'HHOPRE: bizarre SQE',i,nb1,ipt1,nbe1
  485.           mleSQE.lect(NFAMAX+nbn1) = nbe1
  486.           JG = 2
  487.           SEGINI,mlent1
  488.           mlent1.lect(1) = -nbe1
  489.           mlent1.lect(2) = n_o_face
  490.           NBFHHO(nbn1) = nbe1
  491.           MAFHHO(nbn1) = ipt1
  492.           LOFHHO(nbn1) = mlent1
  493.         END DO
  494.  
  495. C= Mise a jour de MSQHHO :
  496.       ELSE
  497.         if (i1_HHO.NE.0) then
  498.           write(ioimp,*) 'HHOPRE-MSQHHO: Bad initialization (2)'
  499.           iret = 5
  500.           return
  501.         end if
  502.         write(ioimp,*) 'MSQHHO deja definie --> A completer'
  503.         write(ioimp,*) 'MAIS CAS EN COURS D IMPLEMENTATION'
  504.         iret = 5
  505.         RETURN
  506.  
  507. C----        meleme = mailSQE
  508. C----        ipt1 = meleme
  509. C----        DO i = 1, MAX(1,nbsSQE)
  510. C----          IF (nbsSQE.NE.0) ipt1 = meleme.lisous(i)
  511. C----          nbn1 = ipt1.num(/1)
  512. C----          nbe1 = ipt1.num(/2)
  513. C----          IF (MAFHHO(nbn1).GT.0) THEN
  514. C----* il faut fusionner les maillages de maniere unique...
  515. C----* en n'ajoutant que les nouvelles faces a la suite des existantes !
  516. C----            ipt2   = MAFHHO(nbn1)
  517. C----            mlent2 = LOFHHO(nbn1)
  518. C----            segact,mlent2*MOD
  519. C----            iadj2 = 0
  520. C----            IF (mlent2.lect(1).LT.0) THEN
  521. C----            END IF
  522. C----c*            CALL INTERB(ipt2,ipt1,ipti,ivid)
  523. C----c*            CALL OUEXCL(ipt1,ipti,iptc,ivid)
  524. C----c*            nbec = iptc.num(/2)
  525. C----c*            CALL FUSMAIL(ipt2,iptc,ipt3,ivid)
  526. C----c*
  527. C----c*            nbe1 = nbe1 + nbec
  528. C----          ELSE
  529. C----            JG = 2
  530. C----            SEGINI,mlent1
  531. C----            mlent1.lect(1) = -nbe1
  532. C----            mlent1.lect(2) = n_o_face
  533. C----          END IF
  534. C----          NBFHHO(nbn1) = nbe1
  535. C----          MAFHHO(nbn1) = ipt1
  536. C----          LOFHHO(nbn1) = mlent1
  537. C----        END DO
  538. C----
  539. C----        nbs1 = 0
  540. C----        DO i = idebf, ifinf
  541. C----          IF (MAFHHO(i).GT.0) nbs1 = nbs1 + 1
  542. C----        END DO
  543. C----        if (nbs1.le.0) then
  544. C----          write(ioimp,*) 'HHOPRE : MSQHHO update incorrect (1)'
  545. C----          iret = 5
  546. C----          return
  547. C----        end if
  548. C----        if (nbs1.lt.MAX(1,nbsSQE)) then
  549. C----          write(ioimp,*) 'HHOPRE : MSQHHO update incorrect (2)'
  550. C----          iret = 5
  551. C----          return
  552. C----        end if
  553.       END IF
  554.  
  555. C= Verification du tableau NBFHHO
  556.       IF (IDIM.EQ.1) THEN
  557.         if (nbfhho(1).le.0) then
  558.           write(ioimp,*) 'HHOPRE: no POI1 in DIME 1?'
  559.           iret = 5
  560.         end if
  561.         do i = 2, NFAMAX
  562.           if (nbfhho(i).gt.0) then
  563.             write(ioimp,*) 'HHOPRE:',i,'-side face in DIME 1'
  564.             iret = 5
  565.           end if
  566.         end do
  567.       END IF
  568.       IF (IDIM.EQ.2) THEN
  569.         if (nbfhho(1).ne.0) then
  570.           write(ioimp,*) 'HHOPRE: POI1 in DIME 2'
  571.           iret = 5
  572.         end if
  573.         NBFHHO(1) = 0
  574.         if (nbfhho(2).le.0) then
  575.           write(ioimp,*) 'HHOPRE: no SEG2 in DIME 2 ?'
  576.           iret = 5
  577.         end if
  578.         do i = 3, NFAMAX
  579.           if (nbfhho(i).gt.0) then
  580.             write(ioimp,*) 'HHOPRE: polygonal face in DIME 2',i,' sides'
  581.             iret = 5
  582.           end if
  583.         end do
  584.       END IF
  585.       IF (IDIM.EQ.3) THEN
  586.         if (nbfhho(1).ne.0) then
  587.           write(ioimp,*) 'HHOPRE: POI1 in DIME 3 ?'
  588.           iret = 5
  589.         end if
  590.         if (nbfhho(2).ne.0) then
  591.           write(ioimp,*) 'HHOPRE: SEG2 in DIME 3 ?'
  592.           iret = 5
  593.         end if
  594.         NBFHHO(1) = 0
  595.         NBFHHO(2) = 0
  596.       END IF
  597.       if (iret.gt.0) return
  598.  
  599.       nfaSQE = 0
  600.       nbsSQE = 0
  601.       DO i = idebf, ifinf
  602.         nfaSQE = nfaSQE + NBFHHO(i)
  603.         IF (MAFHHO(i).GT.0) nbsSQE = nbsSQE + 1
  604.       END DO
  605.       NFAHHO = nfaSQE
  606.       NUFHHO = nbsSQE
  607. c-dbg      write(ioimp,*) 'HHOPRE-MSQHHO:',NFAHHO,NUFHHO
  608.  
  609. C* On reordonne MSQHHO par faces de sommets croissants :
  610.       ipt2 = MSQHHO
  611.       IF (nbsSQE.GT.1) THEN
  612.         segact,ipt2*MOD
  613.         isou = 0
  614.         DO i = idebf, ifinf
  615.           ipt1 = MAFHHO(i)
  616.           IF (ipt1.GT.0) THEN
  617.             isou = isou + 1
  618.             ipt2.lisous(isou) = ipt1
  619.           END IF
  620.         END DO
  621.         segact,ipt2*NOMOD
  622.       END IF
  623. C= Doit-on faire un savseg de MSQHHO et des sous-zones eventuelles ?
  624.  
  625.       NBNN   = 1
  626.       NBELEM = NCEHHO
  627.       NBSOUS = 0
  628.       NBREF  = 0
  629.       IF (i1_HHO.EQ.1) THEN
  630. c-dbg        write(ioimp,*) 'HHOPRE: Initialisation MPCHHO',NCEHHO
  631.         SEGINI,ipt2
  632.         ipt2.itypel = 1
  633.         MPCHHO = ipt2
  634.         NBNEWP = NBELEM
  635.       ELSE
  636. c-dbg        write(ioimp,*) 'HHOPRE: Ajustement MPFHHO'
  637.         ipt2 = MPCHHO
  638.         segact,ipt2*MOD
  639.         NBNEWP = NBELEM - ipt2.num(/2)
  640.         SEGADJ,ipt2
  641. c* Il faut tout decaler par type de face i !
  642.       END IF
  643.       ipoi1 = nbpts
  644.       nbpts = nbpts + NBNEWP
  645.       SEGADJ,MCOORD
  646.       iel2 = 0
  647.       DO i = 1, NCEMAX
  648.         ipt1 = MACHHO(i)
  649.         IF (ipt1.LE.0) GOTO 100
  650.         nel1 = ipt1.num(/2)
  651. c*        nel1 = mleCEL.lect(i+NCEMAX)
  652. c*        jel1 = 1 + mleCEL.lect(i)
  653.         jel1 = 1
  654.         nbn1 = ipt1.num(/1)
  655. c-dbg        write(ioimp,*) 'MACHHO',i,ipt1,nel1,nbn1,iel2,ipoi1
  656.         DO j = jel1, nel1
  657.           jpoi1 = (IDIM+1)*ipoi1
  658.           ipoi1 = ipoi1 + 1
  659.           DO k = 1, IDIM
  660.             r_z = 0.D0
  661.             DO l = 1, nbn1
  662.               lpoi1 = (IDIM+1)*(ipt1.num(l,j)-1)
  663.               r_z = r_z + XCOOR(lpoi1+k)
  664.             END DO
  665.             XCOOR(jpoi1+k) = r_z / nbn1
  666.           END DO
  667.           iel2 = iel2 + 1
  668.           ipt2.num(1,iel2) = ipoi1
  669.         END DO
  670. c-dbg          write(ioimp,*) 'HHOPRE: Verif.',i,iel2,nbchho(i),ipoi1
  671. c---        if (iel2.ne.nbfhho(i)) then
  672. c---          write(ioimp,*) 'HHOPRE(1): inconsistent iel2'
  673. c---        end if
  674.  100    CONTINUE
  675.       END DO
  676.       if (iel2.ne.NCEHHO) then
  677.         iret = 5
  678.         return
  679.       end if
  680.  
  681.       NBNN   = 1
  682.       NBELEM = NFAHHO
  683.       NBSOUS = 0
  684.       NBREF  = 0
  685.       IF (i1_HHO.EQ.1) THEN
  686. c-dbg        write(ioimp,*) 'HHOPRE: Initialisation MPFHHO',NFAHHO
  687.         SEGINI,ipt2
  688.         ipt2.itypel = 1
  689.         MPFHHO = ipt2
  690.         NBNEWP = NBELEM
  691.       ELSE
  692. c-dbg        write(ioimp,*) 'HHOPRE: Ajustement MPFHHO'
  693.         ipt2 = MPFHHO
  694.         segact,ipt2*MOD
  695.         NBNEWP = NBELEM - ipt2.num(/2)
  696.         SEGADJ,ipt2
  697. c* Il faut tout decaler par type de face i !
  698.       END IF
  699.       ipoi1 = nbpts
  700.       nbpts = nbpts + NBNEWP
  701.       SEGADJ,MCOORD
  702.       iel2 = 0
  703.       DO i = idebf, ifinf
  704.         ipt1 = MAFHHO(i)
  705.         IF (ipt1.LE.0) GOTO 200
  706. c*        nel1 = NBFHHO(i)
  707.         nel1 = ipt1.num(/2)
  708. c*        jel1 = 1 + mlent2.lect(i)
  709.         jel1 = 1
  710.         nbn1 = i
  711.         DO j = jel1, nel1
  712.           jpoi1 = (IDIM+1)*ipoi1
  713.           ipoi1 = ipoi1 + 1
  714.           DO k = 1, IDIM
  715.             r_z = 0.D0
  716.             DO l = 1, nbn1
  717.               lpoi1 = (IDIM+1)*(ipt1.num(l,j)-1)
  718.               r_z = r_z + XCOOR(lpoi1+k)
  719.             END DO
  720.             XCOOR(jpoi1+k) = r_z / nbn1
  721.           END DO
  722.           iel2 = iel2 + 1
  723.           ipt2.num(1,iel2) = ipoi1
  724.         END DO
  725. c-dbg          write(ioimp,*) 'HHOPRE: Verif.',i,iel2,nbfhho(i),ipoi1
  726. c---        if (iel2.ne.nbfhho(i)) then
  727. c---          write(ioimp,*) 'HHOPRE(1): inconsistent iel2'
  728. c---        end if
  729.  200    CONTINUE
  730.       END DO
  731.       if (iel2.ne.NFAHHO) then
  732.         iret = 5
  733.         return
  734.       end if
  735.  
  736.       SEGACT,MCOORD*NOMOD
  737.  
  738. C= Remplissage de segments :
  739.       IPOSR = 0
  740.  
  741.       i_z = 0
  742.       CALL HHOLI2('INIT_IPOS',i_z,IPOSR,i_z,iret)
  743.       if (iret.ne.0) return
  744.  
  745.       NMAXR = 0
  746.       DO i = idebf, ifinf
  747.         ipt1 = MAFHHO(i)
  748.         IF (ipt1.GT.0) THEN
  749.           CALL HHOLI2('REMP_IPOS',ipt1,IPOSR,i_z,iret)
  750.           if (iret.ne.0) return
  751.           NMAXR = MAX(NMAXR,i_z)
  752.         END IF
  753.       END DO
  754.       NISFHO = NMAXR
  755.  
  756.       NMAXR = 0
  757.       DO i = 1, NCEMAX
  758.         ipt1 = MACHHO(i)
  759.         IF (ipt1.GT.0) THEN
  760.           CALL HHOLI2('REMP_IPOS',ipt1,IPOSR,i_z,iret)
  761.           if (iret.ne.0) return
  762.           NMAXR = MAX(NMAXR,i_z)
  763.         END IF
  764.       END DO
  765.       NISCHO = NMAXR
  766.  
  767.       i_z = 0
  768.       INDSR = 0
  769.       NMAXR = MAX(NISCHO,NISFHO)
  770.       CALL HHOLI2('INIT_INDS',i_z,NMAXR,INDSR,iret)
  771.       if (iret.ne.0) return
  772.  
  773. c== Remplissage de lentHHO :
  774.       JG = 10
  775.       SEGINI,mlenti
  776.       DO i = 1, JG
  777.         mlenti.lect(i) = -999
  778.       END DO
  779. C= Dimension du probleme :
  780.       mlenti.lect(1) = IDIM
  781. C= Ordre et ddl par face  :
  782.       mlenti.lect(2) = n_o_face
  783.       mlenti.lect(3) = n_d_face
  784. C= Ordre et ddl par cellule :
  785.       mlenti.lect(4) = n_o_cell
  786.       mlenti.lect(5) = n_d_cell
  787. C= Segments de travail :
  788.       mlenti.lect(6) = IPOSR
  789.       mlenti.lect(7) = INDSR
  790. C= Indices utilises ulterieurement :
  791. c-dbg      mlenti.lect( 8) = ...
  792. c-dbg      mlenti.lect( 9) = ...
  793. c-dbg      mlenti.lect(10) = ...
  794.  
  795.       lentHHO = mlenti
  796.  
  797.       SEGSUP,mleCEL,mleSQE
  798.  
  799. c-dbgC= Test independant de la bibliotheque :
  800. c-dbg      write(ioimp,*)
  801. c-dbg      write(ioimp,*) 'HHOPRE - DEBUT TEST_HHO'
  802. c-dbg      CALL TEST_HHO
  803. c-dbg      write(ioimp,*)
  804. c-dbg      write(ioimp,*) 'HHOPRE - FIN TEST_HHO'
  805. c-dbg      write(ioimp,*)
  806.  
  807.       RETURN
  808.       END
  809.  
  810.  
  811.  
  812.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales