TwoPhaseHeatTransfer

heat transfer between two different phases with different temperatures

In this example, the heat transfer between two touching cuboids is simulated using two different phases (chamber 1: Inconel 718, chamber 2: metal) with different temperatures. For Inconel 718, the JohnsonCook material model is used. Here should be a picture

Figure 1: Chamber information.

Since the focus in this example is on heat transfer, only the temperature is simulated, i.e. the calculations for velocity and pressure are switched off by the corresponding TimeIntegration in the KindOfProblem definition for both phases:

V:NONE T:EXPIMP(1.0)

Coupling

Both implicit (both phases solved with vp-, see IntegrationType) and explicit (one phase solved with v-- and the other with vp-, see IntegrationType) coupling of the contact surface between the cuboids are available and can be selected by setting the corresponding alias in USER_common_variables.

begin_alias{} "SOLVER_CHAMBER2" = "vp-" # "vp-" : implicit coupling (both phases solved with vp-) # "v--" : explicit coupling (phase 1 solved with vp- and phase 2 solved with v--) end_alias

Boundary conditions for temperature

The implicit coupling uses the BCON_CNTCT($...$,%ind_T%) functionality.
The explicit coupling is done by using a Neumann condition (%BND_NEUMANN%) with $ \frac{\partial T}{\partial n} = -\frac{\lambda_\mathrm{opp}}{\lambda}\cdot \frac{\partial T}{\partial n}\bigg\vert_\mathrm{opp}$ (corresponding to a heat flow condition) and a Robin condition (%BND_ROBIN%) with respect to the temperature in the opposite phase, i.e. $ T_\mathrm{opp}$.

begin_selection{"SOLVER_CHAMBER2"} case{"vp-"} # implicit coupling via BCON_CNTCT BC_T($contact_INC_METAL$) = ( %BND_ROBIN%, [&heat_trans_coeff_Contact&], [Yopp(%ind_T%)] ) BCON_CNTCT($contact_INC_METAL$,%ind_T%) = ( %BND_ROBIN%, [&heat_trans_coeff_Contact&] ) BCON_CNTCT($contact_METAL_INC$,%ind_T%) = ( %BND_ROBIN%, [&heat_trans_coeff_Contact&] ) BC_T($contact_METAL_INC$) = ( %BND_ROBIN%, [&heat_trans_coeff_Contact&], [Yopp(%ind_T%)] ) case{"v--"} # explicit coupling via heat flow condition dTdn = -LAM_opp/LAM*dTdn_opp BC_T($contact_INC_METAL$) = ( %BND_ROBIN%, [&heat_trans_coeff_Contact&], [Yopp(%ind_T%)] ) BC_T($contact_METAL_INC$) = ( %BND_NEUMANN%, [ -(Yopp(%ind_LAM%)/(Y%ind_LAM%))*equn($NormalDerivativeTemperature$) ] ) end_selection begin_equation{$NormalDerivativeTemperature$} Yopp(%ind_n(1)%)*approxY(%ind_T%, 1, 2, 2.0, 1, Y%ind_x(1)%, Y%ind_x(2)%, Y%ind_x(3)%, 0.0) + Yopp(%ind_n(2)%)*approxY(%ind_T%, 1, 2, 2.0, 2,Y%ind_x(1)%, Y%ind_x(2)%, Y%ind_x(3)%, 0.0) + Yopp(%ind_n(3)%)*approxY(%ind_T%, 1, 2, 2.0, 3,Y%ind_x(1)%, Y%ind_x(2)%, Y%ind_x(3)%,0.0) #= dT/dn (of chamber 1) end_equation BC_T($wall_METAL$) = ( %BND_NEUMANN%, 0.0 ) # --- hydrostatic pressure p --- BC_p($wall_INC$) = ( %BND_wall%, 1.0 ) BC_p($contact_INC_METAL$) = ( %BND_wall%, 1.0 ) BC_p($contact_METAL_INC$) = ( %BND_DIRICH%, 0.0 ) BC_p($wall_METAL$) = ( %BND_wall%, 1.0 ) # --- dynamic pressure p_dyn --- BCON($wall_INC$,%ind_p_dyn%) = ( %BND_NEUMANN%, 0.0 ) BCON($contact_INC_METAL$,%ind_p_dyn%) = ( %BND_NEUMANN%, 0.0 ) BCON($contact_METAL_INC$,%ind_p_dyn%) = ( %BND_DIRICH%, 0.0 ) BCON($wall_METAL$,%ind_p_dyn%) = ( %BND_NEUMANN%, 0.0 ) # --- velocity v --- BC_v($wall_INC$) = ( %BND_wall_nosl% ) BC_v($contact_INC_METAL$) = ( %BND_slip%, 0.0, 0.0, 0.0 ) BC_v($contact_METAL_INC$) = ( %BND_DIRICH%, 0.0, 0.0, 0.0 ) BC_v($wall_METAL$) = ( %BND_wall_nosl% ) ############################# # Initial Condition Section # ############################# INITDATA($INC718$,%ind_T%) = 1000.0 INITDATA($INC718$,%ind_v(1)%) = 0.0 INITDATA($INC718$,%ind_v(2)%) = 0.0 INITDATA($INC718$,%ind_v(3)%) = 0.0 INITDATA($INC718$,%ind_p%) = 0.0 INITDATA($INC718$,%ind_p_dyn%) = 0.0 INITDATA($INC718$,%ind_c%) = 0.0 INITDATA($INC718$,%ind_Sxx%) = 0.0 INITDATA($INC718$,%ind_Sxy%) = 0.0 INITDATA($INC718$,%ind_Sxz%) = 0.0 INITDATA($INC718$,%ind_Syy%) = 0.0 INITDATA($INC718$,%ind_Syz%) = 0.0 INITDATA($INC718$,%ind_Szz%) = 0.0 INITDATA($METAL$,%ind_T%) = &T_amb& INITDATA($METAL$,%ind_v(1)%) = 0.0 INITDATA($METAL$,%ind_v(2)%) = 0.0 INITDATA($METAL$,%ind_v(3)%) = 0.0 INITDATA($METAL$,%ind_p%) = 0.0 INITDATA($METAL$,%ind_p_dyn%) = 0.0 INITDATA($METAL$,%ind_c%) = 0.0 INITDATA($METAL$,%ind_Sxx%) = 0.0 INITDATA($METAL$,%ind_Sxy%) = 0.0 INITDATA($METAL$,%ind_Sxz%) = 0.0 INITDATA($METAL$,%ind_Syy%) = 0.0 INITDATA($METAL$,%ind_Syz%) = 0.0 INITDATA($METAL$,%ind_Szz%) = 0.0 ################## # Active Section # ################## ACTIVE($init_always$) = ( %ACTIVE_init%, %ACTIVE_always% ) #################### # Equation Section # #################### begin_equation{$distanceX$} abs(Y%ind_x(1)%) end_equation begin_equation{$distanceYZ$} max(abs(Y%ind_x(2)% - 0.0005), abs(Y%ind_x(3)% - 0.0005)) end_equation begin_curve{$ControlCuboidX$} depvar_default{equn{$distanceX$}} nb_functions{1} 0.0 1.0 0.001 1.0 0.0010001 0.0 0.01 0.0 end_curve begin_curve{$ControlCuboidYZ$} depvar_default{equn{$distanceYZ$}} nb_functions{1} 0.0 1.0 0.0001 1.0 0.00010001 0.0 0.001 0.0 end_curve begin_equation{$ControlRange$} if (curve{$ControlCuboidX$} > 0) :: if (curve{$ControlCuboidYZ$} > 0) :: 1.0 else :: 0.0 endif else :: 0 endif end_equation begin_equation{$TemperatureCuboid$} if (equn($ControlRange$) > 0) :: Y%ind_T% else :: 0.0 endif end_equation begin_equation{$XcoordinateCuboid$} if (equn($ControlRange$) > 0) :: Y%ind_x(1)% else :: -10.0 endif end_equation begin_curve{$Temperature_AnalyticalSolution$} depvar_default{%ind_x(1)%} include_Ucv{./AnalyticalResults/AnalyticalResultTemperature.dat} end_curve begin_equation{$TemperatureCuboidAnalytical$} if (equn($XcoordinateCuboid$) > -10.0) :: curve{$Temperature_AnalyticalSolution$} else :: 0.0 endif end_equation ####################### # Integration Section # ####################### INTEGRATION($NbPoints$) = (%PUBLICVALUE%, [real(%FLIQUID_NbParticles%)], %INTEGRATION_Header%, "number of points") INTEGRATION($TEMP_MAX_ALL$) = (%MAXIMUM_INT%, [Y%ind_T%], $INC718$, $METAL$, %INTEGRATION_Header%, "temperature max all points") INTEGRATION($TEMP_MIN_ALL$) = (%MINIMUM_INT%, [Y%ind_T%], $INC718$, $METAL$, %INTEGRATION_Header%, "temperature min all points") INTEGRATION($TEMP_AV_ALL$) = (%AVERAGE_INT%, [Y%ind_T%], $INC718$, $METAL$, %INTEGRATION_Header%, "temperature average all points") #begINTEGRATION INTEGRATION($ERROR_UpperBound$) = (%PUBLICVALUE%, [6.0], %INTEGRATION_Header%, "upper bound error") INTEGRATION($SUM_CONTROL_POINTS$) = (%SUMMATION_INT%, [equn($ControlRange$)], $INC718$, $METAL$, %INTEGRATION_Header%, "sum of all control points") INTEGRATION($SUM_ERRORS_L1$) = (%SUMMATION_INT%, [abs(equn($TemperatureCuboid$) - equn($TemperatureCuboidAnalytical$))], $INC718$, $METAL$, %INTEGRATION_Header%, "sum of all errors L1") INTEGRATION($SUM_ERRORS_L2$) = (%SUMMATION_INT%, [(equn($TemperatureCuboid$) - equn($TemperatureCuboidAnalytical$))^2], $INC718$, $METAL$, %INTEGRATION_Header%, "sum of all errors L2") INTEGRATION($ERROR_AVERAGE_L1$) = (%PUBLICVALUE%, [integ($SUM_ERRORS_L1$) / integ($SUM_CONTROL_POINTS$)], %INTEGRATION_Header%, "mean error L1") INTEGRATION($ERROR_AVERAGE_L2$) = (%PUBLICVALUE%, [sqrt(integ($SUM_ERRORS_L2$) / integ($SUM_CONTROL_POINTS$))], %INTEGRATION_Header%, "mean error L2") #endINTEGRATION

Figure 2: Initial (left) and final (right) temperature (implicit coupling).

Validation

The simulation results are compared with the analytical solution, which has been computed with Mathematica (see file 'Analytical_Solution_Diffusion_Equation.nb' in the subfolder 'AnalyticalResults').

INTEGRATION($ERROR_UpperBound$) = (%PUBLICVALUE%, [6.0], %INTEGRATION_Header%, "upper bound error") INTEGRATION($SUM_CONTROL_POINTS$) = (%SUMMATION_INT%, [equn($ControlRange$)], $INC718$, $METAL$, %INTEGRATION_Header%, "sum of all control points") INTEGRATION($SUM_ERRORS_L1$) = (%SUMMATION_INT%, [abs(equn($TemperatureCuboid$) - equn($TemperatureCuboidAnalytical$))], $INC718$, $METAL$, %INTEGRATION_Header%, "sum of all errors L1") INTEGRATION($SUM_ERRORS_L2$) = (%SUMMATION_INT%, [(equn($TemperatureCuboid$) - equn($TemperatureCuboidAnalytical$))^2], $INC718$, $METAL$, %INTEGRATION_Header%, "sum of all errors L2") INTEGRATION($ERROR_AVERAGE_L1$) = (%PUBLICVALUE%, [integ($SUM_ERRORS_L1$) / integ($SUM_CONTROL_POINTS$)], %INTEGRATION_Header%, "mean error L1") INTEGRATION($ERROR_AVERAGE_L2$) = (%PUBLICVALUE%, [sqrt(integ($SUM_ERRORS_L2$) / integ($SUM_CONTROL_POINTS$))], %INTEGRATION_Header%, "mean error L2")

Additionally, all results can be compared in a 1D plot over the x-axis by using the python script 'evaluate_results.py' in the subfolder 'PythonScripts_ComparisonOfResults'.