zkmic
C ZKMIC SOURCE FANDEUR 22/01/03 21:16:02 11136 IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) C************************************************************************ C Operateur KMAC C C OBJET : Cree un objet de type MATRIK C C SYNTAXE : RESU = KMAC INCO UN ; C C RVP : TABLE de soustype EQPR (cree par EQPR) C IMPR : impression du contenu de l'objet' C C REMARQUE : Cet objet n'est pas un objet STANDART CASTEM2000 C Il n'est donc pas listable C Il est tout juste bon a mettre dans la table RVP pour etre utilise C par les operateurs de résolution de la matrice de contrainte C C IKAS=1 KMAC calcul de C uniquement C IKAS=2 KMCT calcul de Ct C IKAS=3 KCCT calcul de C assemblage pour C et Ct C C*********************************************************************** -INC PPARAM -INC CCOPTIO -INC SMCHPOI POINTEUR IZCH2.MCHPOI,IZCCH2.MPOVAL POINTEUR IZDV.MCHPOI,IZDDV.MPOVAL,IZTU1.MPOVAL,TETAN.MPOVAL POINTEUR IZTG1.MCHPOI,IZTGG1.MPOVAL,IZBETA.MPOVAL -INC SMMATRIK -INC SMLENTI POINTEUR IZIPAD.MLENTI,MLENTI1.MLENTI,MLENTI2.MLENTI -INC SMLMOTS POINTEUR LINCO.MLMOTS -INC SMELEME POINTEUR MELEMZ.MELEME,MELEMB.MELEME,MELSTB.MELEME POINTEUR MELEM1.MELEME,MELES1.MELEME,MCTREI.MELEME POINTEUR IGEOM.MELEME,MELEMM.MELEME,MELEMA.MELEME POINTEUR IZLEMC.MELEME,MELEMS.MELEME,MELEMC.MELEME POINTEUR MELEMI.MELEME,MELEMP.MELEME CHARACTER*8 TYPE,TYPC,NOMZ,NOMP,NOMD CHARACTER*8 NOMPP,NOMDD CHARACTER*4 NOM INTEGER IPAD,IPAD2,IK REAL*8 KAUX,TETA1 DIMENSION IXV(3) C DATA IMPR/0/ C************************************************************************* CKMIC C write(6,*)' Operateur KMIC MTABX=',MTABX C C- Récupération de la table EQEX (pointeur MTAB1) C IF(MTAB1.EQ.0)THEN C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 MOTERR( 1: 8) = ' EQEX ' MOTERR( 9:16) = ' EQEX ' MOTERR(17:24) = ' KIZX ' RETURN ENDIF C C- Récupération de la table INCO (pointeur KINC) C IF(KINC.EQ.0)THEN C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 MOTERR( 1: 8) = ' INCO ' MOTERR( 9:16) = ' INCO ' MOTERR(17:24) = ' EQEX ' RETURN ENDIF C************************************************************************* C OPTIONS C KIMPL = 0 -> EXPL 1 -> IMPL 2 -> CN C KFORM = 0 -> EFSI 1 -> EF 2 -> VF 3 -> EFMC C KPRE=3 pression P0 KPRE=4 pression P1 KPRE=2 cas macro 1ère génération IAXI=0 IK=0 IF(IFOMOD.EQ.0)IAXI=2 C C- Récupération de la table des options KOPT (pointeur KOPTI) C IF (KOPTI.EQ.0) THEN C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 MOTERR( 1: 8) = ' KOPT ' MOTERR( 9:16) = ' KOPT ' MOTERR(17:24) = ' KIZX ' RETURN ENDIF IF (IERR.NE.0) RETURN C write(6,*)' Apres les options ' C************************************************************************* C C- Récupération de la table DOMAINE associée au domaine local C IF(MTABZ.EQ.0)THEN C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 MOTERR( 1: 8) = ' DOMZ ' MOTERR( 9:16) = ' DOMZ ' MOTERR(17:24) = ' KIZX ' RETURN ENDIF MACRO1=0 C write(6,*)' KMIC : MACRO1=',macro1 C write(6,*)' KMIC : MQUAD=',MQUAD IF (IERR.NE.0) RETURN MELEMI=MELEME IF(MACRO1.NE.0.AND.KPRE.NE.2)THEN C? CALL KMACRO(MACRO,MELEMM,MTABZ) C? MELEMI=MELEMM MELEMI=MACRO1 ENDIF IF(KPRE.EQ.2.AND.MACRO1.EQ.0)KPRE=3 IF(MQUAD.EQ.0.AND.MACRO1.EQ.0)KPRE=2 IF(KPRE.EQ.3)THEN MELEMP=MELEMC ELSEIF(KPRE.EQ.4)THEN ELSEIF(KPRE.EQ.2)THEN MELEMP=MELEMC ENDIF C************************************************************************* C VERIFICATIONS SUR LES INCONNUES C write(6,*)' Verification sur les inconnues ' TYPE='LISTMOTS' IF(LINCO.EQ.0)GO TO 90 SEGACT LINCO WRITE(6,*)'Operateur KMAC ' C Indice %m1:8 : contient plus de %i1 %m9:16 MOTERR( 1:8) = 'LISTINCO' INTERR(1) = 2 MOTERR(9:16) = ' MOTS ' RETURN ENDIF C On recupere PHI n et TETA n pour Cranck-Nicholson TYPE=' ' IF(TYPE.NE.'CHPOINT ')THEN WRITE(6,*)' Opérateur KMAC :' WRITE(6,*)' L objet CHPOINT ',NOMP, & ' n existe pas dans la table' C Indice %m1:8 : ne contient pas un objet de type %m9:16 MOTERR( 1: 8) = 'INC '//NOMP MOTERR( 9:16) = 'CHPOINT ' RETURN ELSE ENDIF C************************************************************************* C Le domaine de definition est donne par le SPG de la premiere inconnue C Les inconnues suivantes devront posseder ce meme pointeur C On verifie que les points de la zone sont tous inclus dans ce SPG C Inconnue Primale C write(6,*)' Verification inconnue primale ' IF(IKAS.EQ.1.OR.IKAS.EQ.3)THEN MELEMK=MELEMS ELSE MELEMK=MELEMC ENDIF IF(IRET.NE.0)THEN WRITE(6,*)' Opérateur KMAC ' WRITE(6,*)' La zone ',NOMZ,' n''est pas incluse dans le domaine' & , ' de définition de l''inconnue ',NOMP WRITE(6,*)' MELEMK=',melemk,' IGEOM0=',IGEOM0 C Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique MOTERR(1: 8) = 'INC '//NOMP MOTERR(9:16) = 'CHPOINT ' IPAS=0 RETURN ENDIF SEGSUP MLENTI C************************************************************************* TYPE=' ' IF(TYPE.NE.'CHPOINT ')THEN WRITE(6,*)' Opérateur KMAC :' WRITE(6,*)' L objet CHPOINT ',NOMD, & ' n existe pas dans la table' C Indice %m1:8 : ne contient pas un objet de type %m9:16 MOTERR( 1: 8) = 'INC '//NOMD MOTERR( 9:16) = 'CHPOINT ' RETURN ELSE ENDIF NC=TETAN.VPOCHA(/2) C************************************************************************* C Le domaine de definition est donne par le SPG de la premiere inconnue C Les inconnues suivantes devront posseder ce meme pointeur C On verifie que les points de la zone sont tous inclus dans ce SPG C Inconnue Duale C write(6,*)' IGEOM0=',igeom0 IF(IKAS.EQ.1.OR.IKAS.EQ.3)THEN MELEMK=MELEMC ELSE MELEMK=MELEMS ENDIF C write(6,*)' Verification inconnue duale ',MELEMK IF(IRET.NE.0)THEN WRITE(6,*)' Opérateur KMAC ' WRITE(6,*)' La zone ',NOMZ,' n''est pas incluse dans le domaine' & ,' de définition de l''inconnue ',NOMD WRITE(6,*)' MELEMK=',melemk,' IGEOM0=',IGEOM0 C Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique MOTERR(1: 8) = 'INC '//NOMD MOTERR(9:16) = 'CHPOINT ' IPAS=0 RETURN ENDIF SEGSUP MLENTI C************************************************************************* C Lecture du ou des coefficients C Type du coefficient : C IK1=0 CHPOINT IK1=1 scalaire IK1=2 vecteur C write(6,*)' Verification sur les coefficients ' IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)THEN C? MELEMM=MACRO TYPE=' ' SEGACT MELSTB NBELEM=MELSTB.NUM(/2)/4 NBNN=MELSTB.NUM(/1) NBSOUS=0 NBREF=0 SEGINI MELEMA MELEMA.ITYPEL=MELSTB.ITYPEL NKPE=4 IF(IDIM.EQ.3)NKPE=8 do 4878 k=1,nbelem mi=(k-1)*NKPE+1 do 4879 i=1,nbnn MELEMA.num(i,k)=melstb.num(i,mi) 4879 continue C write(6,*)k,(MELEMA.num(i,k),i=1,nbnn) 4878 continue TYPE=' ' TYPE=' ' ENDIF C 1er COEF IXV(1)=MELEMC IXV(2)=1 IXV(3)=0 & MTABX,KINC,1,IXV,IZTG1,IZTGG1,NPT1,NC1,IK1,IRET) IF(IRET.EQ.0)RETURN IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)THEN C 2ème COEF IXV(1)=0 IXV(2)=1 IXV(3)=0 & MTABX,KINC,2,IXV,IZTG2,IZBETA,NPT2,NC2,IK2,IRET) IF(IRET.EQ.0)RETURN ENDIF NRIGE=7 NKID =9 NKMT =7 NMATRI=1 IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)NMATRI=2 SEGINI MATRIK C CAS Stabilisation via MACRO IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)THEN NBME=1 NBSOUS=1 SEGINI IMATRI KSPGP=MCTREI KSPGD=MCTREI SEGACT MELSTB NBSOUS=MELSTB.LISOUS(/1) IF(NBSOUS.NE.0)THEN ENDIF C? SEGACT MELEMM NBCI=MELSTB.NUM(/2) NP =MELSTB.NUM(/1) MP =NP SEGINI IZAFM LIZAFM(1,1)=IZAFM LISDUA(1)=NOMD DO 32 J=1,NP K1=LECT(MELEMA.NUM(J,K)) ii=j do 321 i=1,np u=VPOCHA(K1,I)*IZBETA.VPOCHA(1,1) if(i.eq.1)u=abs(VPOCHA(K1,I))*IZBETA.VPOCHA(1,1) if(ii.le.np)then AM(K,II,J)=U else AM(K,II-NP,J)=U endif ii=ii+1 321 continue 32 CONTINUE 33 CONTINUE SEGSUP MLENTI ENDIF NBME=IDIM C write(6,*)'MELEMI=',MELEMI SEGACT MELEMI NBSOUS=MELEMI.LISOUS(/1) IF(NBSOUS.EQ.0)NBSOUS=1 SEGINI IMATRI IF(IKAS.EQ.2)THEN KSPGD=MELEMS KSPGP=MELEMC IRIGEL(2,1)=MELEMI IRIGEL(1,1)=MELEMP ELSE KSPGP=MELEMS KSPGD=MELEMC IRIGEL(1,1)=MELEMI IRIGEL(2,1)=MELEMP ENDIF SEGACT MELEMP C write(6,*)' ds kmac melemp=',IRIGEL(1,1) C write(6,*)' ds kmac melemd=',IRIGEL(2,1) IRIGEL(4,1)=IMATRI IF(IKAS.EQ.1)IRIGEL(7,1)=3 IF(IKAS.EQ.2)IRIGEL(7,1)=3 IF(IKAS.EQ.3)IRIGEL(7,1)=4 NK=0 DO 11 L=1,NBSOUS IPT1=MELEMI IF(NBSOUS.NE.1)IPT1=MELEMI.LISOUS(L) SEGACT IPT1 IF(IKAS.EQ.2)THEN MP=IPT1.NUM(/1) NP=MELEMP.NUM(/1) ELSE NP=IPT1.NUM(/1) MP=MELEMP.NUM(/1) ENDIF DO 12 I=1,NBME SEGINI IZAFM C write(6,*)' ni izafm np=',np,' mp=',mp,' nbel=',nbel,izafm,l,i LIZAFM(L,I)=IZAFM IF(IKAS.EQ.2)THEN WRITE(NOM,FMT='(I1,A3)')I,NOMD(1:3) LISDUA(I)=NOM//' ' ELSE WRITE(NOM,FMT='(I1,A3)')I,NOMP(1:3) LISDUA(I)=NOMD ENDIF 12 CONTINUE IPM1=LIZAFM(L,1) IPM2=LIZAFM(L,2) IPM3=LIZAFM(L,2) IF(IDIM.EQ.3)IPM3=LIZAFM(L,3) C write(6,*)' AVt KPRISS MACRO1=',MACRO1,KPRE C write(6,*)' APR KPRISS' C ============================= C Option Cranck-Nickolson C ************** TETA1 **** IF (KIMPL.NE.2) TETA1=1.0D0 C ************************* C On recupere le coeficient devant la matrice IF (KIMPL.EQ.2) THEN C write(6,*)' MELEMC=',melemc XV=IZTGG1.VPOCHA(1,1) SEGACT IPT1 SEGACT MELEMC SEGACT IPM1,IPM2,IPM3 NAT=2 NSOUPO=1 N=MELEMC.NUM(/2) C NC=1 C On initialise les segments necessaire a la conception C du second membre SEGINI MCHPO1,MSOUP1,MPOVA1 MCHPO1.IFOPOI=IFOUR MCHPO1.MOCHDE=TITREE MCHPO1.MTYPOI='SMBR' MCHPO1.JATTRI(1)=2 MCHPO1.IPCHP(1)=MSOUP1 DO LN=1,NC MSOUP1.NOCOMP(LN)=LISDUA(LN) END DO MSOUP1.IGEOC=MELEMC MSOUP1.IPOVAL=MPOVA1 SEGACT IZTU1 SEGACT MLENTI1,MLENTI2 KAUX=XV*(1.0D0-TETA1) c DO K=1,NBEL c IK=IK+1 c DO I=1,NP c IPAD=MLENTI2.LECT(MELEMC.NUM(1,IK)) C Par securite on met a zero le second membre a ajouter c MPOVA1.VPOCHA(IPAD,1)=0.0D0 c END DO c END DO IK=IK+1 DO I=1,NP DO J=1,MP C On recupere les bonnes valeurs pour la localisation dans la C matrice pour le produit matriciel. IPAD=MLENTI2.LECT(MELEMC.NUM(J,IK)) IPAD2=MLENTI1.LECT(IPT1.NUM(I,K)) C On effectue le produit matriciel MPOVA1.VPOCHA(IPAD,1)=MPOVA1.VPOCHA(IPAD,1)- & IPM1.AM(K,I,J)*IZTU1.VPOCHA(IPAD2,1)*KAUX IF (IDIM.GT.1) THEN MPOVA1.VPOCHA(IPAD,1)=MPOVA1.VPOCHA(IPAD,1)- & IPM2.AM(K,I,J)*IZTU1.VPOCHA(IPAD2,2)*KAUX END IF IF (IDIM.GT.2) THEN MPOVA1.VPOCHA(IPAD,1)=MPOVA1.VPOCHA(IPAD,1)- & IPM3.AM(K,I,J)*IZTU1.VPOCHA(IPAD2,3)*KAUX END IF END DO END DO END DO SEGDES IPM1,IPM2,IPM3 SEGSUP MLENTI1,MLENTI2 SEGDES MELEMC,IZTU1 SEGDES MCHPO1,MSOUP1,MPOVA1 C On ajoute le second membre a l'ancien (s'il y en avait un). TYPE=' ' C write(6,*)' SMBR ',type IF(TYPE.NE.'CHPOINT')THEN C write(6,*)' On cree un 1er SMBR ' ELSE CALL PRFUSE C ? CALL OPERAD C CALL DTRCHP(MCHPO1) C CALL DTRCHP(MCHPO2) C write(6,*)' On cree un SMBR ' ENDIF END IF C =================================== SEGACT IPT1,IPM1*MOD,IPM2*MOD,IPM3*MOD NP=IPM1.AM(/2) MP=IPM1.AM(/3) NK=NK+1 K1=1+(1-IK1)*(NK-1) XV=IZTGG1.VPOCHA(K1,1) DO I=1,NP DO J=1,MP IF (KIMPL.NE.2) THEN IPM1.AM(K,I,J)=IPM1.AM(K,I,J)*XV IPM2.AM(K,I,J)=IPM2.AM(K,I,J)*XV ELSE IPM1.AM(K,I,J)=TETA1*IPM1.AM(K,I,J)*XV IPM2.AM(K,I,J)=TETA1*IPM2.AM(K,I,J)*XV END IF ENDDO ENDDO IF(IDIM.EQ.3)THEN DO I=1,NP DO J=1,MP IF (KIMPL.NE.2) THEN IPM3.AM(K,I,J)=IPM3.AM(K,I,J)*XV ELSE IPM3.AM(K,I,J)=TETA1*IPM3.AM(K,I,J)*XV END IF ENDDO ENDDO ENDIF 23 CONTINUE SEGDES IPM1,IPM2,IPM3 SEGDES IPT1 11 CONTINUE SEGDES MELEMI C write(6,*)' Fin operateur KMIC' SEGDES IMATRI,MATRIK RETURN 90 CONTINUE WRITE(6,*)' Interruption anormale de KMAC ' C Option %m1:8 incompatible avec les données RETURN 1001 FORMAT(20(1X,I5)) END
© Cast3M 2003 - Tous droits réservés.
Mentions légales