Sources Esope | Cast3M

biokam
C BIOKAM    SOURCE    CB215821  22/07/20    15:39:38     11411          
        SUBROUTINE BIOKAM(MCHPO1,IPT1,ENT)
      IMPLICIT INTEGER(I-N)
        IMPLICIT REAL*8(A-H,O-Z)
 
-INC PPARAM
-INC CCOPTIO
-INC CCREEL
-INC SMELEME
-INC SMCOORD
-INC SMCHPOI
-INC TMTRAV
        SEGMENT ICPR(nbpts)
        SEGMENT STAB
           REAL*8 ICOS(INMAX), INTELL(INM,6)
        END SEGMENT
        INTEGER ENTIER,ENT,I,IVAL,NMAX,NM
        character*4 ichaine
        DIMENSION TAB1(101),TAB2(101)
 
        NMAX=50
        NM=99
        NE=NM+1
 
C============================================
C     ON COMMENCE PAR EXPLORER LE MAILLAGE  =
C============================================
      SEGACT IPT1
      IF (IPT1.ITYPEL.NE.1.OR.IPT1.LISOUS(/1).NE.0) CALL CHANGE(IPT1,1)
 
C==============================================================
C    ON VERIFIE QUE CHAQUE POINT N'APPARAIT QU'UNE SEULE FOIS =
C==============================================================
 
      NNNOE=IPT1.NUM(/2)
      SEGINI ICPR
      ICON=0
      DO 1 I=1,NNNOE
         IKI=IPT1.NUM(1,I)
         IF(ICPR(IKI).NE.0) GO TO 1
         ICON=ICON+1
         ICPR(IKI)=ICON
 1    CONTINUE
      SEGSUP ICPR
 
      IF(ICON.NE.NNNOE) THEN
         CALL ERREUR(75)
         SEGDES IPT1
         RETURN
      ENDIF
 
 
 
C===========================================
C INITIALISATION DES NOMS DES COMPOSANTES  =
C===========================================
 
        NNIN=IDIM+1
        SEGINI MTRAV
        NHRM=NIFOUR
        IF (NIFOUR.EQ.1) NHRM=IFOUR
        INCO(1)='BX  '
        INCO(2)='BY  '
        INCO(3)='BZ  '
        INCO(4)='FLUX'
        NHAR(1)=NHRM
        NHAR(2)=NHRM
        NHAR(3)=NHRM
        NHAR(4)=NHRM
 
 
C================================================
C REMPLISSAGE DE TAB1(T_ELIP) ET TAB2(T_DELIP)  =
C================================================
 
 
        CALL ELPTI7 (NMAX,NM,ENTIER)
        STAB=ENTIER
        SEGACT STAB
        TAB1(1)=0.
        TAB2(1)=0.
        DO I=1,NE
        A=INTELL(I,2)
        AV=0.5*(REAL(I)-1.)/NE
        TAB1(I)=A+(0.707107*(LOG(1.-(4*AV*AV))))-(0.3127313*AV*AV)
          IF (I.GE.3) THEN
          TAB2(I-1)=(TAB1(I)-TAB1(I-2))*NE
          END IF
        END DO
        TAB1(NE+1)=(3*(TAB1(NE)-TAB1(NE-1)))+TAB1(NE-2)
        TAB2(NE)=(TAB1(NE+1)-TAB1(NE-1))*NE
        TAB2(NE+1)=(2*TAB2(NE))-TAB2(NE-1)
        SEGDES STAB
 
C============================
C ON BOUCLE SUR LES NOEUDS  =
C============================
        SEGACT MCHPO1
        DO 2 IPT=1,NNNOE
        L=IPT1.NUM(1,IPT)
        IREF= (IDIM+1)*L
        XM=XCOOR (IREF-IDIM)
        YM=XCOOR (IREF-IDIM+1)
        ZM=XCOOR (IREF-IDIM+2)
        RM=((XM**2)+(YM**2))**0.5
C       write(6,*) 'RM=',RM
        if ( xm.eq.0 ) then
         if ( ym.gt.0 ) then
           theta =  xpi/2.
         else
           if ( ym.lt.0 ) then
           theta = -1.*xpi/2.
           else
           theta = 1.
           end if
         end if
        else
        THETA= ATAN(YM/XM)
        end if
 
C LECTURE DU CHAMP POINT
        BX=0.
        BY=0.
        BZ1=0.
        FLU=0.
          DO 3 I=1,MCHPO1.IPCHP(/1)
          MSOUPO=MCHPO1.IPCHP(I)
          SEGACT MSOUPO
          IPT3=IGEOC
          SEGACT IPT3
          MPOVAL=IPOVAL
          SEGACT MPOVAL
c         write (6,*) 'nbr compo =' ,nocomp(/2)
 
                DO 4 K=1,IPT3.NUM(/2)
                L1=IPT3.NUM(1,K)
                IREF1=(IDIM+1)*L1
                X=XCOOR (IREF1-IDIM)
                Y=XCOOR (IREF1-IDIM+1)
                Z=XCOOR (IREF1-IDIM+2)
                if (nocomp(/2).eq.0) then
                      H = 1.e-5
                      COUR = 1.
                end if
                if (nocomp(/2).eq.1) then
                   if (nocomp(1).eq.'E   ') then
                      H = vpocha(k,1)
                      COUR = 1.
                   else
                      if (nocomp(1).eq.'I   ') then
                      H = 1.e-5
                      COUR = vpocha(k,1)
                      else
                      call erreur (961)
                      end if
                   end if
                end if
                if ( nocomp(/2).eq.2) then
        if ((nocomp(1).eq.'E   ').and.(nocomp(2).eq.'I   ')) then
                      H=VPOCHA(K,1)
                      COUR=VPOCHA(K,2)
                   else
        if ((nocomp(1).eq.'I   ').and.(nocomp(2).eq.'E   ')) then
                      H=VPOCHA(K,2)
                      COUR=VPOCHA(K,1)
                      else
                      call erreur (961)
                      end if
                   end if
                end if
                R=((X**2)+(Y**2))**0.5
                FLU=FLU + (COUR*(FLUXBC(RM,ZM,R,Z,H,NE,TAB1)/1.E+6))
                CALL CHAMPM(RM,ZM,R,Z,H,NE,TAB1,TAB2,BR,BZ)
                IF (XM.LT.0) THEN
                BX=BX- ( COUR*(BR/1.E+6)*COS (THETA))
                BY=BY- ( COUR*(BR/1.E+6)*SIN (THETA))
                BZ1=BZ1 + (COUR*(BZ/1.E+6))
                ELSE
                BX=BX+ ( COUR*(BR/1.E+6)*COS (THETA))
                BY=BY+ ( COUR*(BR/1.E+6)*SIN (THETA))
                BZ1= BZ1 + ( COUR*(BZ/1.E+6))
                END IF
  4             CONTINUE
         SEGDES MPOVAL
         SEGDES IPT3
         SEGDES MSOUPO
  3       CONTINUE
C           WRITE(6,*) 'FLU=',FLU
C           WRITE(6,*) 'BX= ',BX,';BY= ',BY,';BZ= ',BZ1
            IGEO(IPT)=L
            BB(1,IPT)=BX
            IBIN(1,IPT)=1
            BB(2,IPT)=BY
            IBIN(2,IPT)=1
            BB(3,IPT)=BZ1
            IBIN(3,IPT)=1
            BB(4,IPT)=FLU
            IBIN(4,IPT)=1
  2     CONTINUE
        SEGDES IPT1
 
        CALL CRECHP(MTRAV,MCHPO3)
        ENT=MCHPO3
        SEGSUP MTRAV
C       WRITE(6,*) 'COUCOU1',ENT
        RETURN
        END