Download couplage_thermique.dgibi
* fichier : couplage_thermique.dgibi
*
*------------ Cas test couplage_thermique.dgibi --------------------
*
*
* Tests des opérateurs de Castem-fluide
*
*
*- Options générales
*
OPTI DIME 2 ELEM QUA8 ECHO 0 ;
GRAPH = 'N' ;
EPS0 = 1.D-6 ;
*
*-------------------------------------------------------------------
* Couplage thermique entre deux domaines solides l'échange entre les
* deux domaines étant réalisés via un coefficient d'échange : on
* compare à l'état stationnaire solutions calculées et analytique.
*-------------------------------------------------------------------
* La géométrie est un barreau constitué de deux parties de longueur
* respective L1 et L2 et de diffusivité thermique D1 et D2. On note
* h le coefficient d'échange thermique à la frontière commune et
* T1(x) (resp. T2(x)) la température dans le premier barreau (resp.
* dans le deuxième).
*
* Chargement : Flux nul sur les frontières haute et basse (1D) et
* température imposées aux deux extrémités (notées T1 et T2).
*
* Solution stationnaire : Soit Tg (resp. Td) la température au niveau
* de l'interface entre les deux barreaux vue du premier barreau (resp.
* du deuxième). A l'état stationnaire, l'égalité des gradients donne :
* Tg = T1 + (T2-T1)(a2*h)/(a1*h + a2*h + a1*a2)
* Td = T2 + (T1-T2)(a1*h)/(a1*h + a2*h + a1*a2)
* avec a1=D1/L1 et a2=D2/L2.
* On a alors
* T1(x) = T1 + (Tg-T1) x/L1 pour x in [0,L1]
* T2(x) = T2 + (Td-T2) (L1+L2-x)/L2 pour x in [L1,L1+L2]
*
* Application numérique : D1=D2=h=100. ; L1=L2=0.5 ; T1=40. ; T2=20.
* T1(x) = 40 - 10x
* T2(x) = 30 - 10x
* Les valeurs considérées permettent de se placer dans la meme config
* que celle prise par Trio_U dans son dossier de test de validation
* (rapport DTP/SMTH/LMTL/96-22 de Pierre LEDAC).
*-------------------------------------------------------------------
* Auteur : F.Dabbene 01/99
*-------------------------------------------------------------------
* 06/06 : Correction d'une erreur dans le post-traitement
* Tracé de la solution via un mchaml au lieu d'un chpo
* Formulations (DISCR) LINE, MACRO et QUAF possibles
* Ajout du cas permanent
*-------------------------------------------------------------------
*
*--------------------------------------- Début de la procédure VFRON
DEBP VFRON ;
ARGU RVX*TABLE ;
*
* Restriction des champs à la frontière
*
RV = RVX . 'EQEX' ;
TINC = RV . 'INCO' ;
LPT1 = 'DOMA' RVX.'DOMZ' 'SOMMET' ;
TINC . 'T1F' = 'KCHT' RVX.'DOMZ' 'SCAL' 'SOMMET' TINC.'T1' ;
TINC . 'T2F' = 'KCHT' RVX.'DOMZ' 'SCAL' 'SOMMET' TINC.'T2' ;
s2 ma1 = 'KOPS' 'MATRIK' ;
'RESPROC' s2 ma1 ;
'FINPROC' ;
*--------------------------------------- Fin de la procédure VFRON
*
*= MAILLAGE
*
XMIN = 0. ; L1 = 0.5 ; X1 = XMIN + L1 ; L2 = L1 ; X2 = X1 + L2 ;
YMIN = 0. ; DY = 0.5 ; Y1 = YMIN + DY ; Y2 = Y1 + DY ;
NX = 20 ; NY = 1 ;
*
*------------------------------------------ Points du premier domaine
P1 = XMIN YMIN ;
P2 = X1 YMIN ;
P3 = X1 Y2 ;
P4 = XMIN Y2 ;
*------------------------------------------ Points du deuxième domaine
PD1 = P2 ;
PD2 = X2 YMIN ;
PD3 = X2 Y2 ;
PD4 = P3 ;
*------------------------------------------ Lignes du premier domaine
P1P2 = P1 'DROIT' NX P2 ;
P2P3 = P2 'DROIT' NY P3 ;
P3P4 = P3 'DROIT' NX P4 ;
P4P1 = P4 'DROIT' NY P1 ;
*------------------------------------------ Lignes du deuxième domaine
PD1PD2 = PD1 'DROIT' NX PD2 ;
PD2PD3 = PD2 'DROIT' NY PD3 ;
PD3PD4 = PD3 'DROIT' NX PD4 ;
PD4PD1 = 'INVERSE' P2P3 ;
*------------------------------------------ Maillages
*
* Afin d'éviter toute confusion de maillages, ces objets sont à écraser
* après la création des modèles (à remplacer par le maillage qui
* sous-tend chaque modèle)
DOM1 = 'DALLER' P1P2 P2P3 P3P4 P4P1 ;
DOM2 = 'DALLER' PD1PD2 PD2PD3 PD3PD4 PD4PD1 ;
DGAU = P4P1 ;
DDRO = PD2PD3 ;
DBA1 = P1P2 ;
DBA2 = PD1PD2 ;
FRON = P2P3 ;
*
*= Création des MODELES
*
DDOM1 = 'CHANGER' DOM1 'QUAF' ;
DDOM2 = 'CHANGER' DOM2 'QUAF' ;
DDOMT = DDOM1 'ET' DDOM2 ;
DDBA1 = 'CHANGER' DBA1 'QUAF' ;
DDBA2 = 'CHANGER' DBA2 'QUAF' ;
DDGAU = 'CHANGER' DGAU 'QUAF' ;
DDDRO = 'CHANGER' DDRO 'QUAF' ;
DFRON = 'CHANGER' FRON 'QUAF' ;
*
DISCR = 'LINE' ;
$DOM1 = 'MODELISER' DDOM1 'NAVIER_STOKES' DISCR ;
$DOM2 = 'MODELISER' DDOM2 'NAVIER_STOKES' DISCR ;
$DOMT = 'MODELISER' DDOMT 'NAVIER_STOKES' DISCR ;
$DBA1 = 'MODELISER' DDBA1 'NAVIER_STOKES' DISCR ;
$DBA2 = 'MODELISER' DDBA2 'NAVIER_STOKES' DISCR ;
$DGAU = 'MODELISER' DDGAU 'NAVIER_STOKES' DISCR ;
$DDRO = 'MODELISER' DDDRO 'NAVIER_STOKES' DISCR ;
$FRON = 'MODELISER' DFRON 'NAVIER_STOKES' DISCR ;
*
DOM1 = 'DOMA' $DOM1 'MAILLAGE' ;
DOM2 = 'DOMA' $DOM2 'MAILLAGE' ;
DOMT = 'DOMA' $DOMT 'MAILLAGE' ;
DGAU = 'DOMA' $DGAU 'MAILLAGE' ;
DDRO = 'DOMA' $DDRO 'MAILLAGE' ;
DBA1 = 'DOMA' $DBA1 'MAILLAGE' ;
DBA2 = 'DOMA' $DBA2 'MAILLAGE' ;
FRON = 'DOMA' $FRON 'MAILLAGE' ;
*
*= Données numériques
*
D1 = 100. ; D2 = D1 ; H = 100. ;
T1 = 40. ; T2 = 20. ; T12 = T1 + T2 / 2. ;
*
*= Solution analytique
*
DEN1 = D1/L1*H + (D2/L2*H) + (D2*D1/L1/L2) ;
TG = T2 - T1 * D2 / L2 * H /DEN1 + T1 ;
TD = T1 - T2 * D1 / L1 * H /DEN1 + T2 ;
D1X = 'EXTR' ('EVOL' 'CHPO' DBA1 ('COOR' 1 DBA1)) 'ORDO' 1 ;
D2X = 'EXTR' ('EVOL' 'CHPO' DBA2 ('COOR' 1 DBA2)) 'ORDO' 1 ;
X1X = D1X / L1 ;
NBPT = 'NBEL' ('CHAN' 'POI1' DBA1) ;
X2X = D2X - ('PROG' NBPT*L1) / L2 ;
T1X = (TG-T1)*X1X + ('PROG' NBPT*T1) ;
T2X = (T2-TD)*X2X + ('PROG' NBPT*TD) ;
EV1 = 'EVOL' 'ROUG' 'MANU' D1X T1X ;
EV2 = 'EVOL' 'ROUG' 'MANU' D2X T2X ;
EV3 = EV1 'ET' EV2 ;
*
*
*
* PREMIERE MODELISATION DU PROBLEME
*
*-------------------------------------------------------------------
* L'état stationnaire est obtenu en tant que régime asymptotique d'un
* transitoire. La résolution est implicite MAIS l'échange surfacique
* à l'interface est traité par 'ECHI' en semi-explicite. De ce fait,
* seule la solution convergée a un sens, le transitoire est faux.
*-------------------------------------------------------------------
*
*= Description du problème (table EQEX)
*
RV = 'EQEX' $DOMT 'ITMA' 20 'ALFA' 1.
'OPTI' 'EF' 'IMPL'
'ZONE' $DOM1 'OPER' 'LAPN' D1 'INCO' 'T1'
'ZONE' $DOM2 'OPER' 'LAPN' D2 'INCO' 'T2' ;
RV = 'EQEX' RV
'OPTI' 'EF' 'IMPL'
'ZONE' $FRON 'OPER' 'VFRON'
'ZONE' $FRON 'OPER' 'ECHI' H 'T2F' 'INCO' 'T1'
'ZONE' $FRON 'OPER' 'ECHI' H 'T1F' 'INCO' 'T2' ;
RV = 'EQEX' RV
'OPTI' 'EF' 'IMPL' 'CENTREE'
'ZONE' $DOM1 'OPER' 'DFDT' 1. 'T1' 10. 'INCO' 'T1'
'ZONE' $DOM2 'OPER' 'DFDT' 1. 'T2' 10. 'INCO' 'T2' ;
RV = 'EQEX' RV
'CLIM' 'T1' 'TIMP' DGAU T1
'T2' 'TIMP' DDRO T2 ;
*
*= Initialisation de la table INCO et résolution
*
RV . 'INCO' . 'T1' = 'KCHT' $DOM1 SCAL SOMMET T12 ;
RV . 'INCO' . 'T2' = 'KCHT' $DOM2 SCAL SOMMET T12 ;
EXEC RV ;
*
*= Post-traitement de la solution calculée
*
EVC1 = 'EVOL' 'CHPO' DBA1 RV.'INCO'.'T1' ;
EVC2 = 'EVOL' 'CHPO' DBA2 RV.'INCO'.'T2' ;
CT1 = 'EXTR' EVC1 'ORDO' 1 ;
CT2 = 'EXTR' EVC2 'ORDO' 1 ;
EVC3 = ('EVOL' 'VERT' 'MANU' D1X CT1)
'ET' ('EVOL' 'VERT' 'MANU' D2X CT2) ;
*
*= Tracés
*
'SI' ('EGA' GRAPH 'O') ;
CHAM1 = 'CHAN' 'CHAM' RV.'INCO'.'T1' DOM1 ;
CHAM2 = 'CHAN' 'CHAM' RV.'INCO'.'T2' DOM2 ;
'TRAC' $DOMT (CHAM1 'ET' CHAM2) ;
*
TAB1 = 'TABLE' ;
TAB1 . 'TITRE' = 'TABLE' ;
TAB1 . 'TITRE' . 2 = 'Solution castem (1)' ;
TAB1 . 'TITRE' . 1 = '---------------' ;
TAB1 . 'TITRE' . 4 = 'Solution exacte' ;
TAB1 . 'TITRE' . 3 = '---------------' ;
TAB1 . 1 = 'MARQ CROI NOLI' ;
TAB1 . 2 = 'MARQ CROI NOLI' ;
'DESS' (EVC3 'ET' EV3) 'TITR' 'Temperature en fonction de x'
'TITX' 'x' 'TITY' 'T'
'LEGE' TAB1 ;
'FINSI' ;
*
*= Controle erreur
*
SOM1 = 'ABS' (('MAXI' ('SOMM' EV3)) - ('MAXI' ('SOMM' EVC3))) ;
TEST = SOM1 '<' EPS0 ;
'SI' TEST ;
'ERRE' 0 ;
'SINO' ;
'MESS' 'Analytic Integration - Computed = Difference ' ;
'MESS' ('MAXI' ('SOMM' EV3))'-' ('MAXI' ('SOMM' EVC3)) '=' SOM1 ;
'ERRE' 5 ;
'FINSI' ;
*
*
*
* DEUXIEME MODELISATION DU PROBLEME
*
*-------------------------------------------------------------------
* L'état stationnaire est obtenu en tant que régime asymptotique d'un
* transitoire. La résolution est implicite, y compris au niveau de
* l'interface (l'échange surfacique est traité par 'MDIA' et non pas
* par 'ECHI'). Le transitoire a donc aussi un sens physique, le flux
* échangé étant continu pendant le transitoire.
*-------------------------------------------------------------------
*
*= Description du problème (table EQEX)
*
RV = 'EQEX' $DOMT 'ITMA' 5 'ALFA' 1.
'OPTI' EF IMPL
'ZONE' $DOM1 'OPER' 'LAPN' D1 'INCO' 'T1'
'ZONE' $DOM2 'OPER' 'LAPN' D2 'INCO' 'T2' ;
RV = 'EQEX' RV
'OPTI' EFM1 IMPL
'ZONE' $FRON 'OPER' 'MDIA' H 'INCO' 'T1' 'T1'
'ZONE' $FRON 'OPER' 'MDIA' (-1.*H) 'INCO' 'T2' 'T1'
'ZONE' $FRON 'OPER' 'MDIA' H 'INCO' 'T2' 'T2'
'ZONE' $FRON 'OPER' 'MDIA' (-1.*H) 'INCO' 'T1' 'T2' ;
RV = 'EQEX' RV
'OPTI' EF IMPL CENTREE
'ZONE' $DOM1 'OPER' 'DFDT' 1. 'T1' 10. 'INCO' 'T1'
'ZONE' $DOM2 'OPER' 'DFDT' 1. 'T2' 10. 'INCO' 'T2' ;
RV = 'EQEX' RV
'CLIM' 'T1' 'TIMP' DGAU T1
'T2' 'TIMP' DDRO T2 ;
*
*= Initialisation de la table INCO et résolution
*
RV . 'INCO' . 'T1' = 'KCHT' $DOM1 SCAL SOMMET T12 ;
RV . 'INCO' . 'T2' = 'KCHT' $DOM2 SCAL SOMMET T12 ;
EXEC RV ;
*
*= Post-traitement de la solution calculée
*
EVC1 = 'EVOL' 'CHPO' DBA1 RV.'INCO'.'T1' ;
EVC2 = 'EVOL' 'CHPO' DBA2 RV.'INCO'.'T2' ;
CT1 = 'EXTR' EVC1 'ORDO' 1 ;
CT2 = 'EXTR' EVC2 'ORDO' 1 ;
EVC3 = ('EVOL' 'VERT' 'MANU' D1X CT1)
'ET' ('EVOL' 'VERT' 'MANU' D2X CT2) ;
*
*= Tracés
*
'SI' ('EGA' GRAPH 'O') ;
CHAM1 = 'CHAN' 'CHAM' RV.'INCO'.'T1' DOM1 ;
CHAM2 = 'CHAN' 'CHAM' RV.'INCO'.'T2' DOM2 ;
'TRAC' $DOMT (CHAM1 'ET' CHAM2) ;
*
TAB1 = 'TABLE' ;
TAB1 . 'TITRE' = 'TABLE' ;
TAB1 . 'TITRE' . 2 = 'Solution castem (2)' ;
TAB1 . 'TITRE' . 1 = '---------------' ;
TAB1 . 'TITRE' . 4 = 'Solution exacte' ;
TAB1 . 'TITRE' . 3 = '---------------' ;
TAB1 . 1 = 'MARQ CROI NOLI' ;
TAB1 . 2 = 'MARQ CROI NOLI' ;
'DESS' (EVC3 'ET' EV3) 'TITR' 'Temperature en fonction de x'
'TITX' 'x' 'TITY' 'T'
'LEGE' TAB1 ;
'FINSI' ;
*
*= Controle erreur
*
SOM1 = 'ABS' (('MAXI' ('SOMM' EV3)) - ('MAXI' ('SOMM' EVC3))) ;
TEST = SOM1 '<' EPS0 ;
'SI' TEST ;
'ERRE' 0 ;
'SINO' ;
'MESS' 'Analytic Integration - Computed = Difference ' ;
'MESS' ('MAXI' ('SOMM' EV3))'-' ('MAXI' ('SOMM' EVC3)) '=' SOM1 ;
'ERRE' 5 ;
'FINSI' ;
*
*
*
* TROISIEME MODELISATION DU PROBLEME
*
*-------------------------------------------------------------------
* On resout directement le problème stationnaire en full implicite
*-------------------------------------------------------------------
*
*= Description du problème (table EQEX)
*
RV = 'EQEX' $DOMT 'ITMA' 1 'ALFA' 1.
'OPTI' EF IMPL
'ZONE' $DOM1 'OPER' 'LAPN' D1 'INCO' 'T1'
'ZONE' $DOM2 'OPER' 'LAPN' D2 'INCO' 'T2' ;
RV = 'EQEX' RV
'OPTI' EFM1 IMPL
'ZONE' $FRON 'OPER' 'MDIA' H 'INCO' 'T1' 'T1'
'ZONE' $FRON 'OPER' 'MDIA' (-1.*H) 'INCO' 'T2' 'T1'
'ZONE' $FRON 'OPER' 'MDIA' H 'INCO' 'T2' 'T2'
'ZONE' $FRON 'OPER' 'MDIA' (-1.*H) 'INCO' 'T1' 'T2' ;
RV = 'EQEX' RV
'CLIM' 'T1' 'TIMP' DGAU T1
'T2' 'TIMP' DDRO T2 ;
*
*= Initialisation de la table INCO et résolution
*
RV . 'INCO' . 'T1' = 'KCHT' $DOM1 SCAL SOMMET T12 ;
RV . 'INCO' . 'T2' = 'KCHT' $DOM2 SCAL SOMMET T12 ;
EXEC RV ;
*
*= Post-traitement de la solution calculée
*
EVC1 = 'EVOL' 'CHPO' DBA1 RV.'INCO'.'T1' ;
EVC2 = 'EVOL' 'CHPO' DBA2 RV.'INCO'.'T2' ;
CT1 = 'EXTR' EVC1 'ORDO' 1 ;
CT2 = 'EXTR' EVC2 'ORDO' 1 ;
EVC3 = ('EVOL' 'VERT' 'MANU' D1X CT1)
'ET' ('EVOL' 'VERT' 'MANU' D2X CT2) ;
*
*= Tracés
*
'SI' ('EGA' GRAPH 'O') ;
CHAM1 = 'CHAN' 'CHAM' RV.'INCO'.'T1' DOM1 ;
CHAM2 = 'CHAN' 'CHAM' RV.'INCO'.'T2' DOM2 ;
'TRAC' $DOMT (CHAM1 'ET' CHAM2) ;
*
TAB1 = 'TABLE' ;
TAB1 . 'TITRE' = 'TABLE' ;
TAB1 . 'TITRE' . 2 = 'Solution castem (3)' ;
TAB1 . 'TITRE' . 1 = '---------------' ;
TAB1 . 'TITRE' . 4 = 'Solution exacte' ;
TAB1 . 'TITRE' . 3 = '---------------' ;
TAB1 . 1 = 'MARQ CROI NOLI' ;
TAB1 . 2 = 'MARQ CROI NOLI' ;
'DESS' (EVC3 'ET' EV3) 'TITR' 'Temperature en fonction de x'
'TITX' 'x' 'TITY' 'T'
'LEGE' TAB1 ;
'FINSI' ;
*
*= Controle erreur
*
SOM1 = 'ABS' (('MAXI' ('SOMM' EV3)) - ('MAXI' ('SOMM' EVC3))) ;
TEST = SOM1 '<' EPS0 ;
'SI' TEST ;
'ERRE' 0 ;
'SINO' ;
'MESS' 'Analytic Integration - Computed = Difference ' ;
'MESS' ('MAXI' ('SOMM' EV3))'-' ('MAXI' ('SOMM' EVC3)) '=' SOM1 ;
'ERRE' 5 ;
'FINSI' ;
'FIN' ;
© Cast3M 2003 - All rights reserved.
Disclaimer