Coupling electromagnetics and heat transfer modules for wire heating
The purpose of this example is to present how the electromagnetics and heat transfer modules can be coupled together to simulate the heating of materials through electromagnetic effects. In particular, the example presented here involves simulating the heating profile of a copper wire when suppled with the fusing current.
Problem Description
The problem of interest involves suppling large amounts of DC current to a copper wire, such that the temperature of the wire reaches the melting point. To model this be behavior, the follow equations are utilized:
Where:
is the magnetic vector potential,
is ,
is the permeability of free space,
is the angular frequency of the system (60 Hz for this example),
is the electric conductivity of the wire,
is the supplied DC current density,
is the density of copper,
is the heat capacity of copper,
is the temperature,
is the thermal conductivity of the wire, and
is the Joule heating source.
The Joule heating, , is defined by the magnitude of the electric field, , such that:
Instead of solving for the electric field directly, this example used the magnetic vector potential because the magnetic vector potential formulation allows for the directly implementation of a suppled DC current density. This implementation of a suppled DC current density can be better displayed by looking at the relationship between the electric field and the magnetic vector potential:
The third term, is the gradient of the electrostatic potential. This electrostatic potential is the source of the DC current, which is defined as:
Since the magnitude of the electric field is the term of interest for Joule heating, AuxKernels are used to calculate the electric field magnitude from the magnetic vector potential and the supplied DC current density.
commentnote
We cannot define the supplied current density as a volumetric term when solving for the electric field directly, as this would result in an ill-posed electric field.
Boundary Conditions
The boundary conditions for the magnetic vector potential and temperature are:
Where:
is the normal unit vector at the boundary surface,
is the heat flux,
is the convective heat transfer coefficient (10 W/mK for this example), and
is the far-field temperature (293 K for this example).
Example Input File
# This test is a simpified coupled case between the electromagnetic and
# heat transfer modules. While the file microwave_heating.i is a test
# utilizing the method of manufactured solutions, where both real and
# complex components of the electromagnetic properties are provided
# (such that no term is zeroed out), this test involves only the
# real components of the electromagnetic properties. In particular,
# this test supplies the fusing current to a copper wire and simulations
# the spatial and temporal heating profile until the wire reaches its
# melting point. The PDE's of this test file are as follows:
#
# curl(curl(A)) + j*mu*omega*(sigma*A) = J
# mag(E) = mag(-j*omega*A) + mag(J/sigma)
# rho*C*dT/dt - div(k*grad(T)) = Q
# Q = 0.5*sigma*mag(E)^2
#
# Where:
# - A is the magnetic vector potential
# - j is the sqrt(-1)
# - mu is the permeability of free space
# - omega is the angular frequency of the system
# - sigma is the electric conductivity of the wire
# - J is the supplied DC current
# - E is the electric field
# - rho is the density of copper
# - C is the heat capacity of copper
# - T is the temperature
# - k is the thermal conductivity of the wire
# - Q is the Joule heating
#
# The BCs are as follows:
#
# curl(n) x curl(A) = 0, where n is the normal vector
# q * n = h (T - T_infty), where q is the heat flux,
# h is the convective heat transfer coefficient,
# and T_infty is the far-field temperature.
[Mesh<<<{"href": "../../../syntax/Mesh/index.html"}>>>]
# Mesh of the copper wire
[fmg]
type = FileMeshGenerator<<<{"description": "Read a mesh from a file.", "href": "../../../source/meshgenerators/FileMeshGenerator.html"}>>>
file<<<{"description": "The filename to read."}>>> = copper_wire.msh
[]
[]
[Variables<<<{"href": "../../../syntax/Variables/index.html"}>>>]
# The real and complex components of the magnetic vector
# potential in the frequency domain
[A_real]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = NEDELEC_ONE
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
[A_imag]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = NEDELEC_ONE
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
# The temperature of the air in the copper wire
[T]
initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 293.0 #in K
[]
[]
[Kernels<<<{"href": "../../../syntax/Kernels/index.html"}>>>]
### Physics to determine the magnetic vector potential propagation ###
# The propagation of the real component
[curl_curl_real]
type = CurlCurlField<<<{"description": "Weak form term corresponding to $\\nabla \\times (a \\nabla \\times \\vec{E})$.", "href": "../../../source/kernels/CurlCurlField.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = A_real
[]
# Current induced by the electrical conductivity
# of the copper wire
[conduction_real]
type = ADConductionCurrent<<<{"description": "Calculates the current source term in the Helmholtz wave equation using the conduction formulation of the current.", "href": "../../../source/kernels/ADConductionCurrent.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = A_real
field_imag<<<{"description": "The imaginary component of the electric field."}>>> = A_imag
field_real<<<{"description": "The real component of the electric field."}>>> = A_real
conductivity_real<<<{"description": "The real component of the material conductivity."}>>> = electrical_conductivity
conductivity_imag<<<{"description": "The imaginary component of the material conductivity."}>>> = 0.0
ang_freq_real<<<{"description": "The real component of the angular drive frequency."}>>> = omega_real
ang_freq_imag<<<{"description": "The imaginary component of the angular drive frequency."}>>> = 0.0
permeability_real<<<{"description": "The real component of the material permeability."}>>> = mu_real
permeability_imag<<<{"description": "The imaginary component of the material permeability."}>>> = 0.0
component<<<{"description": "Component of field (real or imaginary)."}>>> = real
[]
# Current supplied to the wire
[source_real]
type = VectorBodyForce<<<{"description": "Demonstrates the multiple ways that scalar values can be introduced into kernels, e.g. (controllable) constants, functions, and postprocessors. Implements the weak form $(\\vec{\\psi_i}, -\\vec{f})$.", "href": "../../../source/kernels/VectorBodyForce.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = A_real
function<<<{"description": "A function that defines a vectorValue method. This cannot be supplied with the component parameters."}>>> = mu_curr_real
[]
# The propagation of the complex component
[curl_curl_imag]
type = CurlCurlField<<<{"description": "Weak form term corresponding to $\\nabla \\times (a \\nabla \\times \\vec{E})$.", "href": "../../../source/kernels/CurlCurlField.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = A_imag
[]
# Current induced by the electrical conductivity
# of the copper wire
[conduction_imag]
type = ADConductionCurrent<<<{"description": "Calculates the current source term in the Helmholtz wave equation using the conduction formulation of the current.", "href": "../../../source/kernels/ADConductionCurrent.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = A_imag
field_imag<<<{"description": "The imaginary component of the electric field."}>>> = A_imag
field_real<<<{"description": "The real component of the electric field."}>>> = A_real
conductivity_real<<<{"description": "The real component of the material conductivity."}>>> = electrical_conductivity
conductivity_imag<<<{"description": "The imaginary component of the material conductivity."}>>> = 0.0
ang_freq_real<<<{"description": "The real component of the angular drive frequency."}>>> = omega_real
ang_freq_imag<<<{"description": "The imaginary component of the angular drive frequency."}>>> = 0.0
permeability_real<<<{"description": "The real component of the material permeability."}>>> = mu_real
permeability_imag<<<{"description": "The imaginary component of the material permeability."}>>> = 0.0
component<<<{"description": "Component of field (real or imaginary)."}>>> = imaginary
[]
### Physics to determine the heat transfer ###
# Heat transfer in the copper wire
[HeatTdot_in_copper]
type = ADHeatConductionTimeDerivative<<<{"description": "AD Time derivative term $\\rho c_p \\frac{\\partial T}{\\partial t}$ of the heat equation for quasi-constant specific heat $c_p$ and the density $\\rho$.", "href": "../../../source/kernels/ADHeatConductionTimeDerivative.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
specific_heat<<<{"description": "Property name of the specific heat material property"}>>> = specific_heat_copper
density_name<<<{"description": "Property name of the density material property"}>>> = density_copper
[]
[HeatDiff_in_copper]
type = ADHeatConduction<<<{"description": "Same as `Diffusion` in terms of physics/residual, but the Jacobian is computed using forward automatic differentiation", "href": "../../../source/kernels/ADHeatConduction.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
thermal_conductivity<<<{"description": "the name of the thermal conductivity material property"}>>> = thermal_conductivity_copper
[]
# Heating due the total current
[HeatSrc]
type = ADJouleHeatingSource<<<{"description": "Calculates the heat source term corresponding to electrostatic or electromagnetic Joule heating, with Jacobian contributions calculated using the automatic differentiation system.", "href": "../../../source/kernels/ADJouleHeatingSource.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
heating_term<<<{"description": "Material property providing the Joule heating."}>>> = 'electric_field_heating'
[]
[]
[AuxVariables<<<{"href": "../../../syntax/AuxVariables/index.html"}>>>]
# Decomposing the magnetic vector potential
# for the electric field calculations
[A_x_real]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
[A_y_real]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
[A_x_imag]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
[A_y_imag]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
# The electrical conductivity for the electric
# field calculations
[elec_cond]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
# The electric field profile determined from
# the magnetic vector potential
[E_real]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = NEDELEC_ONE
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
[E_imag]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = NEDELEC_ONE
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
[]
[]
[AuxKernels<<<{"href": "../../../syntax/AuxKernels/index.html"}>>>]
# Decomposing the magnetic vector potential
# for the electric field calculations
[A_x_real]
type = VectorVariableComponentAux<<<{"description": "Creates a field consisting of one component of a coupled vector variable.", "href": "../../../source/auxkernels/VectorVariableComponentAux.html"}>>>
variable<<<{"description": "The name of the variable that this object applies to"}>>> = A_x_real
vector_variable<<<{"description": "The variable from which to compute the component"}>>> = A_real
component<<<{"description": "The component to compute"}>>> = X
[]
[A_y_real]
type = VectorVariableComponentAux<<<{"description": "Creates a field consisting of one component of a coupled vector variable.", "href": "../../../source/auxkernels/VectorVariableComponentAux.html"}>>>
variable<<<{"description": "The name of the variable that this object applies to"}>>> = A_y_real
vector_variable<<<{"description": "The variable from which to compute the component"}>>> = A_real
component<<<{"description": "The component to compute"}>>> = Y
[]
[A_x_imag]
type = VectorVariableComponentAux<<<{"description": "Creates a field consisting of one component of a coupled vector variable.", "href": "../../../source/auxkernels/VectorVariableComponentAux.html"}>>>
variable<<<{"description": "The name of the variable that this object applies to"}>>> = A_x_imag
vector_variable<<<{"description": "The variable from which to compute the component"}>>> = A_imag
component<<<{"description": "The component to compute"}>>> = X
[]
[A_y_imag]
type = VectorVariableComponentAux<<<{"description": "Creates a field consisting of one component of a coupled vector variable.", "href": "../../../source/auxkernels/VectorVariableComponentAux.html"}>>>
variable<<<{"description": "The name of the variable that this object applies to"}>>> = A_y_imag
vector_variable<<<{"description": "The variable from which to compute the component"}>>> = A_imag
component<<<{"description": "The component to compute"}>>> = Y
[]
# The electrical conductivity for the electric
# field calculations
[cond]
type = ADMaterialRealAux<<<{"description": "Outputs element volume-averaged material properties", "href": "../../../source/auxkernels/MaterialRealAux.html"}>>>
property<<<{"description": "The material property name."}>>> = electrical_conductivity
variable<<<{"description": "The name of the variable that this object applies to"}>>> = elec_cond
execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = 'INITIAL LINEAR TIMESTEP_END'
[]
# The magnitude of electric field profile determined
# from the magnetic vector potential using:
# abs(E) = abs(-j*omega*A) + abs(supplied current / elec_cond)
# NOTE: The reason for calculating the magnitude of the electric
# field is the heating term is defined as:
# Q = 1/2 abs(E)^2 for frequency domain field formulations
[E_real]
type = ParsedVectorAux<<<{"description": "Sets a field vector variable value to the evaluation of a parsed expression.", "href": "../../../source/auxkernels/ParsedVectorAux.html"}>>>
coupled_variables<<<{"description": "Vector of coupled variable names"}>>> = 'A_x_imag A_y_imag elec_cond'
expression_x<<<{"description": "Parsed function expression to compute the x component"}>>> = 'abs(2*3.14*60*A_x_imag) + abs(60e6/elec_cond)'
expression_y<<<{"description": "Parsed function expression to compute the y component"}>>> = 'abs(2*3.14*60*A_y_imag)'
variable<<<{"description": "The name of the variable that this object applies to"}>>> = E_real
[]
[E_imag]
type = ParsedVectorAux<<<{"description": "Sets a field vector variable value to the evaluation of a parsed expression.", "href": "../../../source/auxkernels/ParsedVectorAux.html"}>>>
coupled_variables<<<{"description": "Vector of coupled variable names"}>>> = 'A_x_real A_y_real'
expression_x<<<{"description": "Parsed function expression to compute the x component"}>>> = 'abs(-2*3.14*60*A_x_real)'
expression_y<<<{"description": "Parsed function expression to compute the y component"}>>> = 'abs(-2*3.14*60*A_y_real)'
variable<<<{"description": "The name of the variable that this object applies to"}>>> = E_imag
[]
[]
[Functions<<<{"href": "../../../syntax/Functions/index.html"}>>>]
# The supplied current density to the wire
# where only the real x-component is considered
[curr_real_x]
type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../source/functions/MooseParsedFunction.html"}>>>
expression<<<{"description": "The user defined function."}>>> = '60e6' # Units in A/m^2, equivalent to 1178 A in a 5mm diameter wire
[]
# Permeability of free space
[mu_real_func]
type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../source/functions/MooseParsedFunction.html"}>>>
expression<<<{"description": "The user defined function."}>>> = '4*pi*1e-7' # Units in N/A^2
[]
# The angular drive frequency of the system
[omega_real_func]
type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../source/functions/MooseParsedFunction.html"}>>>
expression<<<{"description": "The user defined function."}>>> = '2*pi*60' # Units in rad/s
[]
# The angular frequency time permeability of free space
[omegaMu]
type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../source/functions/MooseParsedFunction.html"}>>>
symbol_names<<<{"description": "Symbols (excluding t,x,y,z) that are bound to the values provided by the corresponding items in the vals vector."}>>> = 'omega mu'
symbol_values<<<{"description": "Constant numeric values, postprocessor names, function names, and scalar variables corresponding to the symbols in symbol_names."}>>> = 'omega_real_func mu_real_func'
expression<<<{"description": "The user defined function."}>>> = 'omega*mu'
[]
# The supplied current density time permeability of free space
[mu_curr_real]
type = ParsedVectorFunction<<<{"description": "Returns a vector function based on string descriptions for each component.", "href": "../../../source/functions/MooseParsedVectorFunction.html"}>>>
symbol_names<<<{"description": "Symbols (excluding t,x,y,z) that are bound to the values provided by the corresponding items in the vals vector."}>>> = 'current_mag mu'
symbol_values<<<{"description": "Constant numeric values, postprocessor names, function names, and scalar variables corresponding to the symbols in symbol_names."}>>> = 'curr_real_x mu_real_func'
expression_x<<<{"description": "x-component of function."}>>> = 'mu * current_mag'
[]
[]
[BCs<<<{"href": "../../../syntax/BCs/index.html"}>>>]
### Temperature boundary conditions ###
# Convective heat flux BC with copper wire
# exposed to air
[surface]
type = ADConvectiveHeatFluxBC<<<{"description": "Convective heat transfer boundary condition with temperature and heat transfer coefficient given by material properties.", "href": "../../../source/bcs/ADConvectiveHeatFluxBC.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = walls
T_infinity<<<{"description": "Material property for far-field temperature"}>>> = 293
heat_transfer_coefficient<<<{"description": "Material property for heat transfer coefficient"}>>> = 10
[]
### Magnetic vector potential boundary conditions ###
# No defined boundary conditions represents
# zero curl conditions at the boundaries, such that:
# A x n = 0
[]
[Materials<<<{"href": "../../../syntax/Materials/index.html"}>>>]
[k]
type = ADGenericConstantMaterial<<<{"description": "Declares material properties based on names and values prescribed by input parameters.", "href": "../../../source/materials/GenericConstantMaterial.html"}>>>
prop_names<<<{"description": "The names of the properties this material will have"}>>> = 'thermal_conductivity_copper'
prop_values<<<{"description": "The values associated with the named properties"}>>> = '397.48' #in W/(m K)
[]
[cp]
type = ADGenericConstantMaterial<<<{"description": "Declares material properties based on names and values prescribed by input parameters.", "href": "../../../source/materials/GenericConstantMaterial.html"}>>>
prop_names<<<{"description": "The names of the properties this material will have"}>>> = 'specific_heat_copper'
prop_values<<<{"description": "The values associated with the named properties"}>>> = '385.0' #in J/(kg K)
[]
[rho]
type = ADGenericConstantMaterial<<<{"description": "Declares material properties based on names and values prescribed by input parameters.", "href": "../../../source/materials/GenericConstantMaterial.html"}>>>
prop_names<<<{"description": "The names of the properties this material will have"}>>> = 'density_copper'
prop_values<<<{"description": "The values associated with the named properties"}>>> = '8920.0' #in kg/(m^3)
[]
# Electrical conductivity (copper is default material)
[sigma]
type = ADElectricalConductivity<<<{"description": "Calculates resistivity and electrical conductivity as a function of temperature, using copper for parameter defaults.", "href": "../../../source/materials/ElectricalConductivity.html"}>>>
temperature<<<{"description": "Coupled variable for temperature."}>>> = T
block<<<{"description": "The list of blocks (ids or names) that this object will be applied"}>>> = copper
[]
# Material that supplies the correct Joule heating formulation
[ElectromagneticMaterial]
type = ElectromagneticHeatingMaterial<<<{"description": "Material class used to provide the electric field as a material property and computes the residual contributions for electromagnetic/electrostatic heating objects.", "href": "../../../source/materials/ElectromagneticHeatingMaterial.html"}>>>
electric_field<<<{"description": "The electric field vector or electrostatic potential scalar to produce the field."}>>> = E_real
complex_electric_field<<<{"description": "The complex component of the electric field vector for the time-harmonic formulation."}>>> = E_imag
electric_field_heating_name<<<{"description": "User-specified material property name for the Joule heating."}>>> = electric_field_heating
electrical_conductivity<<<{"description": "Material property providing electrical conductivity of the material."}>>> = electrical_conductivity
formulation<<<{"description": "The domain formulation of the Joule heating, time or frequency (default = time)."}>>> = FREQUENCY
solver<<<{"description": "Electrostatic or electromagnetic field solver (default = electrostatic)."}>>> = ELECTROMAGNETIC
block<<<{"description": "The list of blocks (ids or names) that this object will be applied"}>>> = copper
[]
# Coefficient for wave propagation
[mu_real]
type = ADGenericFunctionMaterial<<<{"description": "Material object for declaring properties that are populated by evaluation of Function object.", "href": "../../../source/materials/GenericFunctionMaterial.html"}>>>
prop_names<<<{"description": "The names of the properties this material will have"}>>> = mu_real
prop_values<<<{"description": "The corresponding names of the functions that are going to provide the values for the variables"}>>> = mu_real_func
[]
[omega_real]
type = ADGenericFunctionMaterial<<<{"description": "Material object for declaring properties that are populated by evaluation of Function object.", "href": "../../../source/materials/GenericFunctionMaterial.html"}>>>
prop_names<<<{"description": "The names of the properties this material will have"}>>> = omega_real
prop_values<<<{"description": "The corresponding names of the functions that are going to provide the values for the variables"}>>> = omega_real_func
[]
[]
[Preconditioning<<<{"href": "../../../syntax/Preconditioning/index.html"}>>>]
[SMP]
type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
[]
[]
[Executioner<<<{"href": "../../../syntax/Executioner/index.html"}>>>]
type = Transient
scheme = bdf2
solve_type = NEWTON
line_search = NONE
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
dt = 1.0
# NOTE: Change 'end_time' to 10s to accurately simulate the fusing current
# end_time = 10
end_time = 5
automatic_scaling = true
[]
[Outputs<<<{"href": "../../../syntax/Outputs/index.html"}>>>]
exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = true
perf_graph<<<{"description": "Enable printing of the performance graph to the screen (Console)"}>>> = true
[]