konjp1
C KONJP1 SOURCE CB215821 20/11/25 13:32:23 10792 $ ,ICHPSU,MELEMC,MELEFE,MELLIM,IMAT) C C************************************************************************ C C PROJET : CASTEM 2000 C C NOM : KONJP1 C C DESCRIPTION : Voir KON14 C Calcul du jacobien du résidu pour la méthode de C VLH (variable primales = variables primitives, C variable duales = variables conservatives) C C Cas deux dimensions, gaz "calorically perfect" C C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec estensions CISI) C C AUTEUR : A. BECCANTINI, DRN/DMT/SFME/LTMF C C************************************************************************ C C C APPELES (Outils C CASTEM) : KRIPAD, LICHT, ERREUR C C APPELES (Calcul) : VLHJ1, VLHJ3 C C************************************************************************ C C ENTREES C C ILINC : liste des inconnues duales (pointeur d'un LISTMOTS) C C ILINP : liste des inconnues primales (pointeur d'un LISTMOTS) C C 1) Pointeurs des CHPOINT C C IRN : CHPOINT CENTRE contenant la masse volumique ; C C IUN : CHPOINT CENTRE contenant la vitesse ; C C IPN : CHPOINT CENTRE contenant la pression ; C C IGAMN : CHPOINT CENTRE contenant le gamma ; C C INORM : CHPOINT FACE contenant les normales aux faces ; C C ICHPVO : CHPOINT VOLUME contenant le volume C C ICHPSU : CHPOINT FACE contenant la surface des faces C C C 2) Pointeurs de MELEME de la table DOMAINE C C MELEMC : MELEME 'CENTRE' du SPG des CENTRES C C MELEFE : MELEME 'FACEL' du connectivité Faces -> Elts C C MELLIM : MELEME SPG des conditions aux bords C C SORTIES C C IMAT : pointeur de la MATRIK du jacobien du residu C C************************************************************************ C C HISTORIQUE (Anomalies et modifications éventuelles) C C HISTORIQUE : C C************************************************************************ C C C N.B.: On suppose qu'on a déjà controllé RO, P > 0 C GAMMA \in (1,3) C Si non il faut le faire!!! C C************************************************************************ C C C**** Variables de COOPTIO C C INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR C & ,IOTER, IOLEC, IOIMP, IOCAR, IOACQ C & ,IOPER, IOSGB, IOGRA, IOSAU, IORES C & ,IECHO, IIMPI, IOSPI C & ,IDIM, IFICLE, IPREFI C & ,MCOORD C & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE C & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU C & ,NORINC,NORVAL,NORIND,NORVAD C & ,NUCROU, IPSAUV C IMPLICIT INTEGER(I-N) INTEGER ILINC, ILINP, IRN,IUN,IPN,IGAMN,INORM,ICHPVO,ICHPSU & , IMAT, IGEOMC, IGEOMF & , NFAC, NBSOUS, NBREF, NBELEM, NBNN, NRIGE, NMATRI, NKID & , NKMT, NBME, NBEL, MP, NP & , IFAC, NGCF, NLCF, NGCG, NGCD, NLCG, NLCD, NLFL REAL*8 ROG, PG, UXG, UYG, RETG, GAMG, VOLG & , ROD, PD, UXD, UYD, RETD, GAMD, VOLD & , SURF, CNX, CNY, CTX, CTY, FUNCEL & , DFRO(4), DFRET(4), DFRUN(4), DFRUT(4) CHARACTER*8 TYPE C C**** LES INCLUDES C -INC PPARAM -INC CCOPTIO -INC SMCHPOI -INC SMMATRIK -INC SMELEME -INC SMLMOTS -INC SMLENTI POINTEUR MPRN.MPOVAL, MPUN.MPOVAL, MPPN.MPOVAL, MPGAMN.MPOVAL, & MPNORM.MPOVAL, MPVOLU.MPOVAL, MPOVSU.MPOVAL POINTEUR MELEMC.MELEME, MELEMF.MELEME, MELEFE.MELEME, & MELEDU.MELEME,MELLIM.MELEME POINTEUR MLENTC.MLENTI, MLENTF.MLENTI, MLELIM.MLENTI POINTEUR RR.IZAFM, RUX.IZAFM, RUY.IZAFM, RP.IZAFM, & UXR.IZAFM, UXUX.IZAFM, UXUY.IZAFM, UXP.IZAFM, & UYR.IZAFM, UYUX.IZAFM, UYUY.IZAFM, UYP.IZAFM, & RETR.IZAFM, RETUX.IZAFM, RETUY.IZAFM, RETP.IZAFM POINTEUR MLMINC.MLMOTS C C**** KRIPAD pour la correspondance global/local des conditions limits C C SEGACT MELLIM C C**** KRIPAD pour la correspondance global/local des centres C C C SEGACT MLENTC SEGACT MELEMC C SEGACT MELEFE C C C**** LICHT active les MPOVALs en *MOD C C i.e. C C SEGACT MPOVSU*MOD C SEGACT MPOVNO*MOD C SEGACT MPVOLU*MOD C MELEMF = IGEOMF C C SEGACT MLENTF SEGACT MELEMF C C C SEGACT MPRN*MOD C SEGACT MPPN*MOD C SEGACT MPUN*MOD C SEGACT MPGAMN*MOD C NFAC = MELEFE.NUM(/2) C C**** Maillage des inconnues primales C NBSOUS = 0 NBREF = 0 NBELEM = NFAC NBNN = 2 C SEGINI MELEDU C MELEPR = MELEDU C C**** MELEDU = 'SEG2' C MELEDU.ITYPEL = 2 C NRIGE = 7 NMATRI = 1 NKID = 9 NKMT = 7 C SEGINI MATRIK IMAT = MATRIK MATRIK.IRIGEL(1,1) = MELEDU MATRIK.IRIGEL(2,1) = MELEDU C C**** Matrice non symetrique C MATRIK.IRIGEL(7,1) = 2 C NBME = 16 NBSOUS = 1 SEGINI IMATRI MATRIK.IRIGEL(4,1) = IMATRI C C**** Variables primales (primitives) C MLMINC = ILINP SEGACT MLMINC SEGDES MLMINC C C**** Variables duales (conservatives) C MLMINC = ILINC SEGACT MLMINC SEGDES MLMINC C NBEL = NBELEM NBSOUS = 1 NP = 2 MP = 2 SEGINI RR , RUX , RUY , RP , & UXR , UXUX , UXUY , UXP , & UYR , UYUX , UYUY , UYP , & RETR , RETUX , RETUY , RETP C C**** Duale = IMATRI.LISDUA(1) = 'RN' C Primale = IMATRI.LISPRI(1) = 'RN' C -> IMATRI.LIZAFM(1,1) = RR C C Duale = IMATRI.LISDUA(2) = 'RN' C Primale = IMATRI.LISPRI(1) = 'RUXN' C -> IMATRI.LIZAFM(1,2) = RUX C ... C IMATRI.LIZAFM(1,1) = RR IMATRI.LIZAFM(1,2) = RUX IMATRI.LIZAFM(1,3) = RUY IMATRI.LIZAFM(1,4) = RP IMATRI.LIZAFM(1,5) = UXR IMATRI.LIZAFM(1,6) = UXUX IMATRI.LIZAFM(1,7) = UXUY IMATRI.LIZAFM(1,8) = UXP IMATRI.LIZAFM(1,9) = UYR IMATRI.LIZAFM(1,10) = UYUX IMATRI.LIZAFM(1,11) = UYUY IMATRI.LIZAFM(1,12) = UYP IMATRI.LIZAFM(1,13) = RETR IMATRI.LIZAFM(1,14) = RETUX IMATRI.LIZAFM(1,15) = RETUY IMATRI.LIZAFM(1,16) = RETP C DO IFAC = 1, NFAC, 1 NGCF = MELEFE.NUM(2,IFAC) NLCF = MLENTF.LECT(NGCF) IF(NLCF .NE. IFAC)THEN WRITE(IOIMP,*) 'Il ne faut pas jouer avec la table domaine' GOTO 9999 ENDIF NLFL = MLELIM.LECT(NGCF) NGCG = MELEFE.NUM(1,IFAC) NGCD = MELEFE.NUM(3,IFAC) IF(NLFL .NE. 0)THEN C C********** The point belongs on BC -> No contribution to jacobian! C MELEDU.NUM(1,IFAC) = NGCG MELEDU.NUM(2,IFAC) = NGCD ELSEIF(NGCG .NE. NGCD)THEN C C********** Les MELEMEs C MELEDU.NUM(1,IFAC) = NGCG MELEDU.NUM(2,IFAC) = NGCD C C********** Les etats G et D C NLCG = MLENTC.LECT(NGCG) NLCD = MLENTC.LECT(NGCD) C ROG = MPRN.VPOCHA(NLCG,1) PG = MPPN.VPOCHA(NLCG,1) UXG = MPUN.VPOCHA(NLCG,1) UYG = MPUN.VPOCHA(NLCG,2) GAMG = MPGAMN.VPOCHA(NLCG,1) RETG= (PG / (GAMG - 1.0D0)) + 0.5D0 * ROG * (UXG * UXG + & UYG * UYG) VOLG = MPVOLU.VPOCHA(NLCG,1) C ROD = MPRN.VPOCHA(NLCD,1) PD = MPPN.VPOCHA(NLCD,1) UXD = MPUN.VPOCHA(NLCD,1) UYD = MPUN.VPOCHA(NLCD,2) GAMD = MPGAMN.VPOCHA(NLCD,1) RETD= (PD / (GAMD - 1.0D0)) + 0.5D0 * ROD * (UXD * UXD + & UYD * UYD) VOLD = MPVOLU.VPOCHA(NLCD,1) C C********** La normale G->D C La tangente C SURF = MPOVSU.VPOCHA(NLCF,1) CNX = MPNORM.VPOCHA(NLCF,1) CNY = MPNORM.VPOCHA(NLCF,2) CTX = -1.0D0 * CNY CTY = CNX C C********** La contribution de Gauche C & DFRO,DFRUN,DFRUT,DFRET) C C C********** AB.AM(IFAC,IPRIM,IDUAL) C A = nom de l'inconnu duale (Ro,rUX,rUY,RET) C B = nom de l'inconnu primale (Ro,UX,UY,P) C IPRIM = 1, 2 -> G, D C IDUAL = 1, 2 -> G, D C i.e. C A_IDUAL = AB.AM(IFAC,IPRIM,IDUAL) * B_IPRIM + ... C C C********** Dual RN C FUNCEL = SURF * DFRO(1) RR.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG RR.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * DFRO(2) RUX.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG RUX.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * DFRO(3) RUY.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG RUY.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * DFRO(4) RP.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG RP.AM(IFAC,1,2) = FUNCEL / VOLD C C********** Dual RUXN C FUNCEL = SURF * (DFRUN(1) * CNX + DFRUT(1) * CTX) UXR.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UXR.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * (DFRUN(2) * CNX + DFRUT(2) * CTX) UXUX.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UXUX.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * (DFRUN(3) * CNX + DFRUT(3) * CTX) UXUY.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UXUY.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * (DFRUN(4) * CNX + DFRUT(4) * CTX) UXP.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UXP.AM(IFAC,1,2) = FUNCEL / VOLD C C********** Dual RUYN C FUNCEL = SURF * (DFRUN(1) * CNY + DFRUT(1) * CTY) UYR.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UYR.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * (DFRUN(2) * CNY + DFRUT(2) * CTY) UYUX.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UYUX.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * (DFRUN(3) * CNY + DFRUT(3) * CTY) UYUY.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UYUY.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * (DFRUN(4) * CNY + DFRUT(4) * CTY) UYP.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UYP.AM(IFAC,1,2) = FUNCEL / VOLD C C********** Dual RETN C FUNCEL = SURF * DFRET(1) RETR.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG RETR.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * DFRET(2) RETUX.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG RETUX.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * DFRET(3) RETUY.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG RETUY.AM(IFAC,1,2) = FUNCEL / VOLD C FUNCEL = SURF * DFRET(4) RETP.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG RETP.AM(IFAC,1,2) = FUNCEL / VOLD C C C********** La contribution de D C CNX = -1.0D0 * CNX CNY = -1.0D0 * CNY CTX = -1.0D0 * CTX CTY = -1.0D0 * CTY & DFRO,DFRUN,DFRUT,DFRET) C C C********** Dual RN C FUNCEL = SURF * DFRO(1) RR.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD RR.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * DFRO(2) RUX.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD RUX.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * DFRO(3) RUY.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD RUY.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * DFRO(4) RP.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD RP.AM(IFAC,2,1) = FUNCEL / VOLG C C********** Dual RUXN C FUNCEL = SURF * (DFRUN(1) * CNX + DFRUT(1) * CTX) UXR.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD UXR.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * (DFRUN(2) * CNX + DFRUT(2) * CTX) UXUX.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD UXUX.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * (DFRUN(3) * CNX + DFRUT(3) * CTX) UXUY.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD UXUY.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * (DFRUN(4) * CNX + DFRUT(4) * CTX) UXP.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD UXP.AM(IFAC,2,1) = FUNCEL / VOLG C C********** Dual RUYN C FUNCEL = SURF * (DFRUN(1) * CNY + DFRUT(1) * CTY) UYR.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD UYR.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * (DFRUN(2) * CNY + DFRUT(2) * CTY) UYUX.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD UYUX.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * (DFRUN(3) * CNY + DFRUT(3) * CTY) UYUY.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD UYUY.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * (DFRUN(4) * CNY + DFRUT(4) * CTY) UYP.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD UYP.AM(IFAC,2,1) = FUNCEL / VOLG C C********** Dual RETN C FUNCEL = SURF * DFRET(1) RETR.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD RETR.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * DFRET(2) RETUX.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD RETUX.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * DFRET(3) RETUY.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD RETUY.AM(IFAC,2,1) = FUNCEL / VOLG C FUNCEL = SURF * DFRET(4) RETP.AM(IFAC,2,2) = -1.0D0 * FUNCEL / VOLD RETP.AM(IFAC,2,1) = FUNCEL / VOLG C ELSE C C********** Murs (NGCG = NGCD) C C C********** Les MELEMEs C MELEDU.NUM(1,IFAC) = NGCG MELEDU.NUM(2,IFAC) = NGCD NLCG = MLENTC.LECT(NGCG) C ROG = MPRN.VPOCHA(NLCG,1) PG = MPPN.VPOCHA(NLCG,1) UXG = MPUN.VPOCHA(NLCG,1) UYG = MPUN.VPOCHA(NLCG,2) GAMG = MPGAMN.VPOCHA(NLCG,1) VOLG = MPVOLU.VPOCHA(NLCG,1) C C********** La normale sortante C SURF = MPOVSU.VPOCHA(NLCF,1) CNX = MPNORM.VPOCHA(NLCF,1) CNY = MPNORM.VPOCHA(NLCF,2) C & DFRUN) C C********** Dual RN C RR.AM(IFAC,1,1) = 0.0D0 RR.AM(IFAC,1,2) = 0.0D0 C RUX.AM(IFAC,1,1) = 0.0D0 RUX.AM(IFAC,1,2) = 0.0D0 C RUY.AM(IFAC,1,1) = 0.0D0 RUY.AM(IFAC,1,2) = 0.0D0 C RP.AM(IFAC,1,1) = 0.0D0 RP.AM(IFAC,1,2) = 0.0D0 C C********** Dual RUXN C FUNCEL = SURF * DFRUN(1) * CNX UXR.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UXR.AM(IFAC,1,2) = 0.0D0 C FUNCEL = SURF * DFRUN(2) * CNX UXUX.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UXUX.AM(IFAC,1,2) = 0.0D0 C FUNCEL = SURF * DFRUN(3) * CNX UXUY.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UXUY.AM(IFAC,1,2) = 0.0D0 C FUNCEL = SURF * DFRUN(4) * CNX UXP.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UXP.AM(IFAC,1,2) = 0.0D0 C C********** Dual RUYN C FUNCEL = SURF * DFRUN(1) * CNY UYR.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UYR.AM(IFAC,1,2) = 0.0D0 C FUNCEL = SURF * DFRUN(2) * CNY UYUX.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UYUX.AM(IFAC,1,2) = 0.0D0 C FUNCEL = SURF * DFRUN(3) * CNY UYUY.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UYUY.AM(IFAC,1,2) = 0.0D0 C FUNCEL = SURF * DFRUN(4) * CNY UYP.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG UYP.AM(IFAC,1,2) = 0.0D0 C C********** Dual RETN C RETR.AM(IFAC,1,1) = 0.0D0 RETR.AM(IFAC,1,2) = 0.0D0 C RETUX.AM(IFAC,1,1) = 0.0D0 RETUX.AM(IFAC,1,2) = 0.0D0 C RETUY.AM(IFAC,1,1) = 0.0D0 RETUY.AM(IFAC,1,2) = 0.0D0 C RETP.AM(IFAC,1,1) = 0.0D0 RETP.AM(IFAC,1,2) = 0.0D0 C C********** Dual RN C RR.AM(IFAC,2,2) = 0.0D0 RR.AM(IFAC,2,1) = 0.0D0 C RUX.AM(IFAC,2,2) = 0.0D0 RUX.AM(IFAC,2,1) = 0.0D0 C RUY.AM(IFAC,2,2) = 0.0D0 RUY.AM(IFAC,2,1) = 0.0D0 C RP.AM(IFAC,2,2) = 0.0D0 RP.AM(IFAC,2,1) = 0.0D0 C C********** Dual RUXN C UXR.AM(IFAC,2,2) = 0.0D0 UXR.AM(IFAC,2,1) = 0.0D0 C UXUX.AM(IFAC,2,2) = 0.0D0 UXUX.AM(IFAC,2,1) = 0.0D0 C UXUY.AM(IFAC,2,2) = 0.0D0 UXUY.AM(IFAC,2,1) = 0.0D0 C UXP.AM(IFAC,2,2) = 0.0D0 UXP.AM(IFAC,2,1) = 0.0D0 C C********** Dual RUYN C UYR.AM(IFAC,2,2) = 0.0D0 UYR.AM(IFAC,2,1) = 0.0D0 C UYUX.AM(IFAC,2,2) = 0.0D0 UYUX.AM(IFAC,2,1) = 0.0D0 C UYUY.AM(IFAC,2,2) = 0.0D0 UYUY.AM(IFAC,2,1) = 0.0D0 C UYP.AM(IFAC,2,2) = 0.0D0 UYP.AM(IFAC,2,1) = 0.0D0 C C********** Dual RETN C RETR.AM(IFAC,2,2) = 0.0D0 RETR.AM(IFAC,2,1) = 0.0D0 C RETUX.AM(IFAC,2,2) = 0.0D0 RETUX.AM(IFAC,2,1) = 0.0D0 C RETUY.AM(IFAC,2,2) = 0.0D0 RETUY.AM(IFAC,2,1) = 0.0D0 C RETP.AM(IFAC,2,2) = 0.0D0 RETP.AM(IFAC,2,1) = 0.0D0 C ENDIF ENDDO C SEGDES MELEMC SEGDES MELEFE SEGDES MELEMF C SEGDES MPOVSU SEGDES MPVOLU SEGDES MPNORM C SEGDES MPRN SEGDES MPPN SEGDES MPUN SEGDES MPGAMN C SEGDES MELEDU SEGDES MATRIK SEGDES IMATRI C SEGDES RR , RUX , RUY , RP , & UXR , UXUX , UXUY , UXP , & UYR , UYUX , UYUY , UYP , & RETR , RETUX , RETUY , RETP SEGSUP MLENTC SEGSUP MLENTF SEGSUP MLELIM C 9999 CONTINUE RETURN END
© Cast3M 2003 - Tous droits réservés.
Mentions légales