* fichier : rayo-2D-5.dgibi
*************************************************************
* *
* Calcul d'une plaque infinie (de largeur 2L) avec tempé- *
* rature imposée au milieu et soumise au rayonnement. *
* Il s'agit de trouver la température au bord en régime *
* permanent. *
* *
* Néanmoins nous cherchons la solution en faisant un *
* calcul transitoire avec PASAPAS pour tester celle-ci *
* en configuration température imposée + rayonnement vers *
* l'infini. *
* *
* Modélisation plane. *
* *
* Auteurs : Michel Bulik & Nadia Coulon *
* *
* Date : Avril 1997 *
* *
*************************************************************
*
*** Options ...
opti dime 2 elem seg2 echo 1 ;
*** Paramètres ...
L = 1. ;
ep = L / 10 ;
graph = faux ;
Tempini = 273. ;
Tempext = 0. * Tempini ;
cte_sb = 5.673e-8 ;
Gamma = 1. ;
si ( neg Gamma 0. ) ;
lambda = (cte_sb * Tempini * Tempini * Tempini * L) / Gamma ;
mess 'Lambda = ' lambda ;
valemis = 1. ;
sinon ;
lambda = 1. ;
mess 'Lambda choisie arbitrairement = ' lambda ;
valemis = 0. ;
finsi ;
*** Points ...
dens 0.1 ;
p1 = 0 0 ;
p2 = 0 ep ;
vechoriz = L 0 ;
*** Lignes ...
li1 = p1 d 1 p2 ;
*** Surface ...
opti elem qua4 ;
su1 = li1 tran vechoriz dini 0.1 dfin 0.01 ;
li2 = cote 3 su1 ;
p3 = li2 poin proc vechoriz ;
li3 = cote 2 su1 ;
si(graph) ;
titr 'Le maillage de la plaque' ;
trac su1 ;
finsi ;
*** Modèles ...
mocnd = mode su1 thermique ;
macnd = mate mocnd 'K' lambda 'RHO' 1. 'C' lambda ;
moray = mode li2 thermique rayonnement INFINI ;
maray = mate moray 'EMIS' valemis ;
*** Température imposée ...
blt = bloq 'T' li1 ;
ti = depi blt Tempini ;
lr1 = prog 0 100 ;
lr2 = prog 1 1 ;
ev1 = evol manu 't' lr1 'T' lr2 ;
chti = char 'TIMP' ev1 ti ;
*** Température extérieure ...
chtext = manu chpo li2 1 'T' Tempext nature diffus ;
chte = char 'TERA' ev1 chtext ;
*** Préparation de la table pour PASAPAS ...
tabnl = table ;
tabnl . MODELE = mocnd et moray ;
tabnl . CARACTERISTIQUES = macnd et maray ;
tabnl . BLOCAGES_THERMIQUES = blt ;
tabnl . CHARGEMENT = chti et chte ;
tabnl . TEMPS_CALCULES = prog 0. pas 0.1 1.
pas 0.2 2. ;
tabnl . TEMPERATURES = table ;
tabnl . TEMPERATURES . 0 = manu chpo su1 1 'T' Tempini ;
* tabnl . RAYONNEMENT = table ;
* tabnl . RAYONNEMENT . 1 = table ;
* tabnl . RAYONNEMENT . 1 . TYPE = 'INFINI' ;
* tabnl . RAYONNEMENT . 1 . MODELE = moray ;
tabnl . 'PROCEDURE_THERMIQUE' = 'DUPONT' ;
tabnl . 'CTE_STEFAN_BOLTZMANN' = cte_sb ;
tabnl . 'RELAXATION_THETA' = 1. ;
*** Appel à PASAPAS ...
pasapas tabnl ;
*** Petit post-traitement ...
listt = prog ;
listtp3 = prog ;
listtp3r = prog ;
nbpas = dime (tabnl . TEMPS) ;
repeter surpas nbpas ;
lindice = &surpas - 1 ;
listt = listt et (prog (tabnl . TEMPS . lindice)) ;
valtp3 = extr (tabnl . TEMPERATURES . lindice) T p3 ;
listtp3 = listtp3 et (prog valtp3) ;
valtp3r = (valtp3 - Tempext) / (Tempini - Tempext) ;
listtp3r = listtp3r et (prog valtp3r) ;
* ... Vérification s'il n'y a pas d'oscillations spatiales ...
si(faux) ;
chtit = chai 'Profil de temperature au temps '
(tabnl . TEMPS . lindice) ;
titr chtit ;
profil = evol chpo (tabnl . TEMPERATURES . lindice)
'T' li3 ;
dess profil ;
finsi ;
fin surpas ;
titr 'Evolution de la temperature au bord de la plaque' ;
evtp3 = evol manu 't' listt 'T' listtp3 ;
titr 'Evolution de la temperature relative au bord de la plaque' ;
evtp3r = evol manu 't' listt 'T' listtp3r ;
si(graph) ;
dess evtp3 ;
dess evtp3r ;
finsi ;
*** Vérification si c'est OK ...
************************************************************************
* *
* On va résoudre analytiquement le problème en écrivant l'égalité *
* des flux (T = température au bord de la plaque) : *
* *
* T_imp - T *
* lambda ------------- = epsilon * sigma * (T^4 - T_ext^4) *
* L *
* *
************************************************************************
A_0 = -1 * ((lambda * Tempini / L) +
(valemis * cte_sb * (Tempext**4))) ;
A_1 = lambda / L ;
A_2 = 0.0 ;
A_3 = 0.0 ;
A_4 = valemis * cte_sb ;
r1 = racp A_0 A_1 A_2 A_3 A_4 ;
* mess 'Zéros du polynôme :' r1 ' ' r2 ' ' r3 ' ' r4 ;
sol_anlt = r1 extr 2 ;
sol_calc = extr listtp3 (dime listtp3) ;
erreur = 100 * (sol_calc - sol_anlt) / sol_anlt ;
mess 'Erreur relative = ' erreur '%' ;
si(abs(erreur) > 0.01) ;
erre 5 ;
sinon ;
erre 0 ;
finsi ;
*** Bye ...
fin ;
© Cast3M 2003 - All rights reserved.
Disclaimer