TUTORIAL #3 - Analisi termica non lineare con Code_Aster

Stampa

In questo terzo tutorial verrà illustrata la grande flessibilità del software Code_Aster nella gestione di condizioni al contorno particolari. L'esempio riguarda il raffreddamento per immersione in un bagno temprante di un componente in acciaio. Nei tutorial precedenti il raffreddamento era guidato da un coefficiente di scambio termico h [W/m^2 K ] costante; in realtà tale coefficiente subisce forti variazioni in funzione della temperatura superficiale del componente. Un generico mezzo temprante in fase liquida (sia acqua che olio) attraversa infatti 3 fasi di raffreddamento dette :

- film boiling : valori di coeff. h bassi

- nucleation/bubble boiling : valori di coeff.h molto alti

- convezione naturale : valori di coeff. h medio/bassi

 Nell'immagine sottostante un esempio di quanto il coefficiente di scambio per l'acqua e per l'olio possa variare durante le fasi di raffreddamento.

Fonte : ALD VACUUM TECHNOLOGIES

 A livello numerico la simulazione di un raffreddamento con coefficiente h variabile con la temperatura diventa pertanto non lineare in campo termico, Code_Aster in questi casi prevede l'utilizzo della funzione THER_NON_LINE.

La condizione al contorno di scambio termico tuttavia non prevede un coefficiente h funzione della temperatura h(T) ma solamente funzione del tempo : h(t); bisogna quindi sviluppare una condizione al contorno, definita dall'utente, che gestisca lo scambio di calore alla superficie.

Inizialmente vengono inseriti nel software i valori di coefficiente h in funzione delle diverse temperature di superficie, quindi vengono interpolati in modo lineare su tutto il dominio generando così la funzione "hT"; si definisce quindi la FORMULA "hec" (heat exchange coefficient) dello scambio termico funzione della temperatura T.

hec(T) = hT(T) * (Text - T) 

con Text la temperatura del bagno di tempra.

L'algoritmo di simulazione è stato testato su di un generico componente di media complessità a forma cilindrica ottenuto da GrabCad. Il componente alla temperatura di 850° in acciaio è stato immerso in un bagno di tempra a base olio. In questa simulazione per semplicità è stato considerato un tempo di immersione così ridotto da trascurare la velocità di ingresso nel bagno temprante. Di seguito i profili di temperatura rispettivamente dopo 2 e 10 secondi dall'immersione.

 

Profili di temperatura alla superficie ed in sezione dopo 2 secondi
Profili di temperatura alla superficie ed in sezione dopo 10 secondi

 CONCLUSIONI

Si può notare come, a causa dei diversi regimi di raffreddamento implementati, il gradiente di temperatura tra superficie e cuore varia in maniera non lineare rispetto al tempo: dopo 2 s il deltaT è di circa 200°C mentre dopo 10 s il deltaT è aumentato a 340 °C per poi ridursi via via che si raggiunge la temperatura di equilibrio del bagno. Questo avviene perchè alcune zone raggiungono una velocità di raffreddamento superiore ( passaggio tra film boiling e nucleation boiling) alcuni istanti prima di altre e questo fa si che la forbice di valori si allarghi notevolmente.

 

 

Di seguito il commfile:

DEBUT();

mesh=LIRE_MAILLAGE(FORMAT='MED',);

#-------- CONDUCTIVITY FUNCTION DEFINITION ---------

lam=DEFI_FONCTION(
                  NOM_PARA='TEMP',
                  VALE=(20.0 , 42 ,),
                  PROL_DROITE='CONSTANT',
                  PROL_GAUCHE='CONSTANT',);

#-------- RHO_CP FUNCTION DEFINITION ---------

rho_cp=DEFI_FONCTION(
                     NOM_PARA='TEMP',
                     VALE=(20.0 ,3588000.0 ,),
                     PROL_DROITE='CONSTANT',
                     PROL_GAUCHE='CONSTANT',);

#-------- HEAT EXCHANGE FUNCTION DEFINITION ---------
#hT is function of the wall temperature (TEMP), at the right and left (DROITE / GAUCHE) border the value is mantained constant

hT=DEFI_FONCTION(
                 NOM_PARA='TEMP',
                 VALE=(300.0 ,300.0 ,
                       385.0 ,425.0 ,
                       453.0 ,2574.0 ,
                       676.0 ,709.0 ,
                       900.0 ,138.0 ,),
                 PROL_DROITE='CONSTANT',
                 PROL_GAUCHE='CONSTANT',);

#-------- THERMAL FLUX FUNCTION DEFINITION ---------
#the same formula of heat exchange is defined here

hec = FORMULE(VALE='hT(TEMP) * (60-TEMP)',
              NOM_PARA='TEMP',);

#-------  INTERPOLATION OF THE HEC FORMULA --------
#first is necessary to define the interpolation domain via a DEFI_LIST_REEL, next the values are interpolated with CALC_FONC_INTERP

t1=DEFI_LIST_REEL(DEBUT=60.0,
                  INTERVALLE=_F(JUSQU_A=900.0,
                                NOMBRE=20,),);

hec2=CALC_FONC_INTERP(FONCTION=hec,
                      LIST_PARA=t1,
                      PROL_DROITE='CONSTANT',
                      PROL_GAUCHE='CONSTANT',);

#-------- MATERIAL DEFINITION --------
#It is important to remember that is a NON LINEAR simulation so THER_NL command has to be used

steel=DEFI_MATERIAU(THER_NL=_F(LAMBDA=lam,
                               RHO_CP=rho_cp,),);

mod_t=AFFE_MODELE(MAILLAGE=mesh,
                  AFFE=_F(TOUT='OUI',
                          PHENOMENE='THERMIQUE',
                          MODELISATION='3D',),);

mat_t=AFFE_MATERIAU(MODELE=mod_t,
                    AFFE=_F(TOUT='OUI',
                            MATER=steel,),
                    AFFE_VARC=_F(NOM_VARC='TEMP',
                                 VALE_REF=25.0,),);

#-------- BOUNDARY CONDITION DEFINITION ---------
#the heat exchange is defined as a non linear thermal flux (FLUX_NL) over the quenching surface, the interpolated hec function is used

c_ter=AFFE_CHAR_THER_F(MODELE=mod_t,
                       FLUX_NL=_F(GROUP_MA='quench',
                                  FLUN=hec2,),);

t2=DEFI_LIST_REEL(DEBUT=0.0,
                  INTERVALLE=(_F(JUSQU_A=5.0,
                                 PAS=0.1,),
                              _F(JUSQU_A=10.0,
                                 PAS=0.2,),
                              _F(JUSQU_A=20.0,
                                 PAS=0.5,),),);

tempo=DEFI_LIST_INST(DEFI_LIST=_F(METHODE='MANUEL',
                                  LIST_INST=t2,),);

#-------- THERMAL ANALYSIS ---------
#Similar to stat_non_line, to reach better convergence use REAC_ITER=1 and increase ITER_GLOB_MAXI

tnlin=THER_NON_LINE(MODELE=mod_t,
                    CHAM_MATER=mat_t,
                    COMPORTEMENT=_F(RELATION='THER_NL',),
                    EXCIT=_F(CHARGE=c_ter,),
                    INCREMENT=_F(LIST_INST=tempo,),
                    ETAT_INIT=_F(VALE=850.0,),
                    NEWTON=_F(REAC_ITER=1,),
                    CONVERGENCE=_F(ITER_GLOB_MAXI=20,),);

#-------- SAVING RESULTS ---------

IMPR_RESU(FORMAT='MED',
          RESU=_F(RESULTAT=tnlin,),);

FIN();