C CAAVCT SOURCE MAGN 09/07/30 21:15:17 6458 SUBROUTINE CAAVCT C************************************************************************ C C MODIFICATIONS : IZDD -> IZD2 (chgt de diago) C TEST SUR MTABD C SIMPLIFICATION DES PARAMETRES DE PF500 (on passe C MTABD) C C************************************************************************ C C .Rajout de commentaires et de nouveaux messages d'erreurs en utilisant C la routine ERREUR de K2000 : F.D Juillet 96 C .Correction d'erreurs dans le cas des CHPO Centre : on teste le nom C des inconnues et on effectue le calcul pour toutes les composantes C des inconnues : P.G Aout 96 C************************************************************************ IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) CHARACTER*8 NOMI,TYPE,TYPE1,NOM,NOMZ,TYPC,TYP0,TYPINC CHARACTER*8 NOMIX(10) CHARACTER*9 NOMI1 CHARACTER*4 NOC,KOMP(4) C -INC CCOPTIO -INC SMLENTI POINTEUR IZIPAD.MLENTI,MLH.MLENTI -INC SMCHAML -INC SMMATRIK -INC SMCHPOI POINTEUR CLIM.MCHPOI POINTEUR IPHI.MPOVAL,IPHR.MPOVAL,IPHH.MPOVAL POINTEUR IZD.MCHPOI ,IZD0.MCHPOI ,IZG.MCHPOI ,IZGD.MCHPOI POINTEUR IZDD.MPOVAL,IZDD0.MPOVAL,IZGG.MPOVAL,IZGDD.MPOVAL POINTEUR IZD2.MPOVAL,IZSS.MPOVAL,IZS.MCHPOI,IZH.MCHPOI POINTEUR IZPHI.MCHPOI -INC SMLMOTS POINTEUR MINCOG.MLMOTS -INC SMELEME POINTEUR MAH.MELEME POINTEUR MELEMI.MELEME POINTEUR IGEOM0.MELEME -INC SMLREEL -INC SMEVOLL PARAMETER (NTB=1) DIMENSION KTAB(NTB) CHARACTER*8 LTAB(NTB) DATA KOMP/'UX ','UY ','UZ ','SCAL'/ DATA LTAB /'EQEX '/ C C- Lecture des tables transmises en arguments C ------------------------------------------ CALL LITABS(LTAB,KTAB,NTB,1,IRET) IF (IRET.EQ.0)THEN WRITE(6,*)' Opérateur AVCT :' WRITE(6,*)' On attend un ensemble de table soustypes' RETURN ENDIF KIZC = KTAB(1) C C- Lecture facultative de la "cfl" C ------------------------------ CALL LIRREE(ALFA,0,IRET) IF (IRET.EQ.0) ALFA=1.D0 C C- Lecture facultative des paramètres pour les impressions de controle C ------------------------------------------------------------------- CALL LIRCHA(NOM,0,IRET) IMPR = 0 IF (IRET.EQ.0) THEN IMPR = 0 ELSEIF (NOM.EQ.'IMPR ') THEN CALL LIRENT(IMPR,1,IRET) IF (IRET.EQ.0) RETURN ENDIF C C- Récuperation des pointeurs des tables INCO et EQEX !! C- BUG dans LITABS puisque EQEX etait facultatif --> à corriger C TYPE=' ' CALL ACMO(KIZC,'INCO',TYPE,INCO) IF(IERR.NE.0) RETURN C C- Récupération des noms d'inconnues C --------------------------------- TYPE = 'LISTMOTS' CALL ACMO(KIZC,'LISTINCO',TYPE,MLMOT2) IF (IERR.NE.0) RETURN SEGACT MLMOT2 NBINC1 = MLMOT2.MOTS(/2) C C- Récupérations des pointeurs associés aux tables KIZD et KIZG C ------------------------------------------------------------ C (contenant la matrice "masse" diagonale et l'"incrément".) C TYPE=' ' CALL ACMO(KIZC,'KIZD ',TYPE,KIZD) TYPE=' ' CALL ACMO(KIZC,'KIZK ',TYPE,KIZK) TYPE=' ' CALL ACMO(KIZC,'KIZG ',TYPE,KIZG) C C- Initialisation des données temporelles C -------------------------------------- CALL ACMF(KIZC,'TFINAL',TFINAL) TYPE=' ' CALL ACMO(KIZC,'PASDETPS',TYPE,MTABT) IF (MTABT.NE.0) THEN CALL ACMF(MTABT,'DELTAT',DT) IF (IERR.NE.0) RETURN CALL ACMF(MTABT,'TPS',TPS) IF (IERR.NE.0) RETURN CALL ACME(MTABT,'NUPASDT',IPT) IF (IERR.NE.0) RETURN ELSE DT = 1.D0 TPS = 0.D0 IPT = 1 ENDIF DT = DT * ALFA IFINAL=0 TPS1 = TPS + DT C write(6,*)' TPS=',tps,' TFINAL=',tfinal,' DT=',dt IF(TPS.GT.TFINAL)THEN IFINAL=1 GO TO 800 ELSEIF(TPS1.GT.TFINAL)THEN DT=TFINAL-TPS IFINAL=1 ENDIF TPS = TPS + DT C C- Récupération des pointeurs associés aux tables KIZG1 et KIZS C ------------------------------------------------------------ C TYPE=' ' CALL ACMO(KIZC,'KIZG1 ',TYPE,KIZG1) TYPE=' ' CALL ACMO(KIZC,'KIZS ',TYPE,KIZS) C C============================================================== C Boucle principale : Traitement de chaque inconnue de LISTINCO C============================================================== C DO 1 L=1,NBINC1 C C- NOMI est l'identifiant de l'inconnue dans chaque table C NOMI = MLMOT2.MOTS(L) C C- Activation du MPOVAL de l'inconnue au temps précédant C ----------------------------------------------------- TYPE = 'CHPOINT ' CALL ACMO(INCO,NOMI,TYPE,IZPHI) IF (IERR.NE.0) RETURN CALL LICHT(IZPHI,IPHI,TYP0,IGEOM0) NPT = IPHI.VPOCHA(/1) NC = IPHI.VPOCHA(/2) C C- Activation du terme source éventuel associé à l'inconnue NOMI C ------------------------------------------------------------- IZSS = 0 IGEOMS = 0 IF (KIZS.NE.0) THEN TYPE = 'CHPOINT ' CALL ACMO(KIZS,NOMI,TYPE,IZS) IF (IERR.NE.0) RETURN CALL LICHT(IZS,IZSS,TYP0,IGEOMS) IF (IGEOMS.NE.IGEOM0) THEN MOTERR(1: 8) = NOMI(1:4)//'KIZS' MOTERR(9:16) = NOMI(1:4)//'INCO' INTERR(1) = 1 CALL ERREUR(698) RETURN ENDIF ENDIF C C- Identification du MPOVAL de l'inconnue au nouveau pas de temps C -------------------------------------------------------------- C (pour écrasement éventuel des valeurs du pas précédant) C IPHR = IPHI C C- Correspondance entre numérotation globale et locale C --------------------------------------------------- CALL KRIPAD(IGEOM0,IZIPAD) SEGACT IGEOM0 C C- Cas ou il y a des matrices masses non diagonales C ------------------------------------------------ IF(KIZK.NE.0)THEN TYPE=' ' MATRIK=KIZK SEGACT MATRIK NBK=IRIGEL(/2) DO 411 K=1,NBK IMATRI=IRIGEL(4,K) SEGACT IMATRI NBME=LISPRI(/2) WRITE(6,*)' INCO : ',(LISPRI(ii),ii=1,nbme) 411 CONTINUE ENDIF C C- Cas où les tables KIZD et KIZG sont données C ------------------------------------------- IF (KIZD.NE.0.AND.KIZG.NE.0) THEN C C- Recherche de la matrice "masse" diagonale associée à l'inconnue NOMI C -------------------------------------------------------------------- TYPE = ' ' CALL ACMO(KIZD,NOMI,TYPE,IZD) IF (TYPE.NE.'CHPOINT ') THEN WRITE(6,*)'pas de matrice diagonale pour ',NOMI GOTO 1 ENDIF CALL LICHT(IZD,IZDD,TYPC,IGEOM) C C- Recherche de l'"incrément" associé à l'inconnue NOMI C ---------------------------------------------------- TYPE = ' ' CALL ACMO(KIZG,NOMI,TYPE,IZG) IF (TYPE.NE.'CHPOINT ') THEN GOTO 1 ENDIF CALL LICHT(IZG,IZGG,TYPC,IGEOM) C C- Imposition des conditions aux limites sur l'"incrément" C ------------------------------------------------------ TYPE=' ' CALL ACMO(KIZC,'KIZI ',TYPE,KIZI) IF (KIZI.NE.0) THEN DO 211 I=1,NC NOC = KOMP(I) IF (NC.EQ.1) NOC=KOMP(4) TYPE = ' ' CALL ACMO(KIZI,NOMI,TYPE,MCHPOI) IF (TYPE.EQ.'CHPOINT') THEN CALL LACHT(MCHPOI,MPOVAL,NOC,MELEMI) IF (MPOVAL.NE.0) THEN SEGACT MELEMI LONG = MELEMI.NUM(/2) CALL RSETX2(IZGG.VPOCHA(1,I),MELEMI.NUM, & LONG,VPOCHA,IZIPAD.LECT) SEGDES MELEMI,MPOVAL ENDIF ENDIF 211 CONTINUE ENDIF C C- Calcul de l'inconnue au nouveau pas de temps C -------------------------------------------- C- En EF C- En VF C DO 11 I=1,NC C C -------------------- Cas de des inconnues situées au Centre C -------------------------------------- IF (TYPC.EQ.'FACE') THEN TYPE=' ' CALL ACMO(KIZC,'DOMAINE',TYPE,MTABD) IF (MTABD.EQ.0) THEN MOTERR(1: 8) = 'TABLE ' MOTERR(9:16) = 'DOMAINE ' CALL ERREUR(79) RETURN ENDIF TYPE = 'CHPOINT ' cccc CALL ACMO(KIZD,'DIAC0',TYPE,MCHPO4) CALL ACMO(KIZD,NOMI,TYPE,MCHPO4) IF (IERR.NE.0) RETURN CALL LICHT(MCHPO4,IZD2,TYPINC,IGEOM2) CALL PF500(IZD2,IPHI,IZSS,IZGG,IPHR,DT,MTABD,I) C if(I .EQ. NC) SEGDES MCHPO4,IZD2 C C -------------------- Cas de des inconnues situées au Sommet C -------------------------------------- ELSEIF (TYPC.EQ.'SOMMET') THEN IF (KIZG1.EQ.0) THEN LONG=IZDD.VPOCHA(/1) CALL P500(IZDD.VPOCHA(1,I),IPHI.VPOCHA(1,I), & IZGG.VPOCHA(1,I),IPHR.VPOCHA(1,I),LONG,DT) ELSE TYPE = ' ' CALL ACMO(KIZG1,NOMI,TYPE,IZGD) IF (TYPE.NE.'CHPOINT ') THEN LONG = IZDD.VPOCHA(/1) CALL P500(IZDD.VPOCHA(1,I),IPHI.VPOCHA(1,I), & IZGG.VPOCHA(1,I),IPHR.VPOCHA(1,I),LONG,DT) ELSE CALL LICHT(IZGD,IZGDD,TYPC,IGEOM) LONG=IZDD.VPOCHA(/1) CALL P501(IZDD.VPOCHA(1,I),IPHI.VPOCHA(1,I), & IZGG.VPOCHA(1,I),IPHR.VPOCHA(1,I), & LONG,DT,IZGDD.VPOCHA(1,I)) ENDIF ENDIF ENDIF 11 CONTINUE C C- Imposition des conditions aux limites C ------------------------------------- TYPE = ' ' CALL ACMO(KIZC,'CLIM',TYPE,MCHPOI) IF (TYPE.EQ.'CHPOINT')THEN DO 111 I=1,NC IF (NC.EQ.1) THEN NOC=NOMI(1:4) ELSE WRITE(NOC,FMT='(I1,A3)')I,NOMI(1:3) ENDIF SEGACT MCHPOI NSOUPO = IPCHP(/1) DO 10111 NSP=1,NSOUPO MSOUPO = IPCHP(NSP) SEGACT MSOUPO NCOMP = NOCOMP(/2) DO 10112 NCP=1,NCOMP IF (NOCOMP(NCP).EQ.NOC) THEN MELEMI = IGEOC MPOVAL = IPOVAL SEGACT MELEMI,MPOVAL LONG = MELEMI.NUM(/2) CALL RSETX2(IPHR.VPOCHA(1,I),MELEMI.NUM, & LONG,VPOCHA(1,NCP),IZIPAD.LECT) SEGDES MELEMI,MPOVAL ENDIF 10112 CONTINUE SEGDES MSOUPO 10111 CONTINUE SEGDES MCHPOI 111 CONTINUE ENDIF C C- Mise à zero de l'ensemble des données et ménage C ----------------------------------------------- LONG = IZGG.VPOCHA(/1)*IZGG.VPOCHA(/2) CALL INITD(IZGG.VPOCHA,LONG,0.D0) SEGDES IZG,IZGG SEGDES IZD,IZDD IF (KIZG1.NE.0.AND.IZGD.NE.0) THEN CALL INITD(IZGDD.VPOCHA,LONG,0.D0) SEGDES IZGD,IZGDD ENDIF IF (KIZS.NE.0.AND.IZSS.NE.0) THEN CALL INITD(IZSS.VPOCHA,LONG,0.D0) ENDIF C -------------------------------------------------- C- Fin du cas ou les tables KIZD et KIZG sont données C -------------------------------------------------- ENDIF C C- Traitement des historiques C -------------------------- IF (MTABT.NE.0) THEN TYPE=' ' CALL ACMO(KIZC,'HIST',TYPE,KHIS) IF (KHIS.NE.0)THEN NUCOMP = IPHR.VPOCHA(/2) NOMIX(1) = NOMI NOMIX(2) = '1'//NOMI(1:7) NOMIX(3) = '2'//NOMI(1:7) NOMIX(4) = '3'//NOMI(1:7) DO 82 NUC=1,NUCOMP+1 NUCR = 1 IF (NUC.GT.1) NUCR=NUC-1 TYPE = ' ' CALL ACMO(KHIS,NOMIX(NUC),TYPE,MEVOLL) IF (TYPE.EQ.'EVOLUTIO')THEN TYPE1=' ' NOMI1='$'//NOMIX(NUC) CALL ACMO(KHIS,NOMI1,TYPE1,MAH) IF (TYPE1.EQ.'MAILLAGE') THEN SEGACT MAH ENDIF CALL ECROBJ('CHPOINT',IZPHI) CALL ECROBJ('MAILLAGE',MAH) CALL REDU CALL LIROBJ('CHPOINT',IZH,1,IRET) CALL LICHT(IZH,IPHH,TYPC,IGEOM) TYPE1=' ' CALL ACMO(KHIS,'KFIH',TYPE1,KFIH) IF (TYPE1.EQ.'ENTIER') THEN CALL ACME(KHIS,'KFIH',KFIH) ELSE KFIH = 20 ENDIF INDH = IPT - IPT/KFIH * KFIH IF (INDH.NE.0) GOTO 80 SEGACT MEVOLL NH=IEVOLL(/1) DO 81 IH=1,NH KEVOLL=IEVOLL(IH) SEGACT KEVOLL MLREE1=IPROGX IF (IH.EQ.1)THEN SEGACT MLREE1 JG=MLREE1.PROG(/1)+1 SEGADJ MLREE1 MLREE1.PROG(JG)=TPS SEGDES MLREE1 ENDIF MLREE2=IPROGY SEGACT MLREE2 SEGADJ MLREE2 if(ih.le.iphh.vpocha(/1)) then MLREE2.PROG(JG)=IPHH.VPOCHA(IH,NUCR) else MLREE2.PROG(JG)=0 endif SEGDES MLREE2,KEVOLL 81 CONTINUE SEGDES MEVOLL 80 CONTINUE SEGDES MAH ENDIF 82 CONTINUE ENDIF ELSE WRITE(6,*)' Pour des historiques il faut une table PASDETPS' ENDIF C SEGDES IPHI,IPHR SEGDES IGEOM0 SEGSUP IZIPAD 1 CONTINUE C SEGDES MLMOT2 C C- Impressions de controle C ----------------------- IF (IMPR.NE.0) THEN KFIDT = IMPR IF (MTABT.NE.0) THEN CALL ACMM(MTABT,'OPER',NOMI) CALL ACMM(MTABT,'ZONE',NOMZ) CALL ACMF(MTABT,'DTCONV',DTT1) CALL ACMF(MTABT,'DTDIFU',DTT2) CALL ACMF(MTABT,'DIAEL',DIAEL) CALL ACME(MTABT,'NUEL',NUEL) CALL ACME(MTABT,'NUPASDT',IPT) ELSE IPT=1 ENDIF IND = IPT - IPT/KFIDT * KFIDT IF (IPT.EQ.1) THEN WRITE(6,*) & ' IPT : NUMERO DU PAS DE TEMPS , NUEL : NUMERO DE L ELEMENT , ' &,' DIAEL : DIAMETRE MOYEN DE L ELEMENT ' WRITE(6,*) &' ALFA : TOLERANCE SUR LE PAS DE TEMPS , DTMAX : PAS DE TEMPS MAX' &,' DTT1 : PAS DE TEMPS DE CONVECTION , DTT2 PAS DE TEMPS DE' &,' DIFFUSION' ENDIF IF (IND.EQ.0)THEN WRITE(6,1011)NOMZ,NOMI WRITE(6,1010)IPT,NUEL,DIAEL,ALFA,DT,DTT1,DTT2 ENDIF ENDIF C C- Mise à jour de la table PASDETPS C -------------------------------- 800 CONTINUE IF (KIZD.NE.0.AND.KIZG.NE.0) THEN IF (MTABT.NE.0) THEN CALL ECMF(MTABT,'DELTAT-1',DT) CALL ECMF(MTABT,'TPS',TPS) DT=1.E30 CALL ECMF(MTABT,'DELTAT',DT) IPT = IPT + 1 CALL ECME(MTABT,'NUPASDT',IPT) ENDIF ENDIF CALL ECRENT(IFINAL) RETURN C C- Formats associés aux impression de controle C ------------------------------------------- 1010 FORMAT(2X,'N.DT',I5,' NU.EL',I5,' DIAEL=',1PE11.4,' ALFA=', & 1PE11.4,' DTMAX=',1PE11.4, & ' DT1=',1PE11.4,' DT2=',1PE11.4) 1011 FORMAT(2X,' ZONE :',A8,' OPERATEUR :',A8) END