A discontinuous Petrov-Galerkin finite-element framework for the simulation of microwave-heated flows

1. Introduction

• Microwave heating can accelerate reactions and improve selectivity, but field sensitivity to materials and geometry and difficulty measuring reactor fields complicate industrial scale-up.The chemical industry contributes approximately 5% of global CO2 emissions, and reducing production is not considered a viable solution.Source: §1 · p. 2
• High-wavenumber Maxwell simulations require carefully chosen mesh size and polynomial degree, while the resulting indefinite systems challenge iterative-solver preconditioning.The growing degrees of freedom make direct solvers impractical at large scale; resonance modes and material discontinuities further hinder convergence.Source: §1 · p. 3

1.1. Discontinuous Petrov-Galerkin method

• DPG produces a Hermitian positive definite system suited to memory-efficient Conjugate Gradient solves and provides absolute stability and an error estimator for adaptive refinement.Its applicability to Maxwell equations is supported by structural analogy with Helmholtz and numerical evidence, while wavenumber-explicit estimates have been established only for acoustics.Source: §1.1 · p. 3

1.2. Novelty

• The work contributes an implementation guide and integrates a time-harmonic DPG electromagnetic solver into lethe for coupled flow, heat-transfer, and electromagnetism simulations.The authors describe this as the first integration of DPG within an open-source multiphysics framework outside the DPG community; discretization orders can be selected independently by physics.Source: §1.2 · pp. 3–4

2.1. Abstract problem formulation

• The practical DPG formulation uses broken trial and enriched test spaces, plus trace unknowns that weakly enforce continuity; the error representation function measures the discrete residual.The field and trace unknowns act as Lagrange multipliers enforcing orthogonality of the error representation function to the trial and trace spaces.Source: §2.1 · p. 4

2.2. Model problem and ultraweak formulation

• The Maxwell solver uses an ultraweak formulation to retain both electric and magnetic fields, accommodate magnetic losses, and impose low inter-element regularity requirements.It equips the test space with an adjoint graph norm; its Robin-boundary contribution supports convergence when Robin conditions are present.Source: §2.2 · pp. 5–6

2.3. Numerical implementation and discretization

• The implementation discretizes the Maxwell fields with hexahedral elements, uses test polynomials one degree higher than trial polynomials, and condenses the system to trace unknowns.The test-space enrichment parameter is Δp=1; after solving for trace unknowns, interior fields and the error representation function are recovered locally, and the element residual norm supplies the adaptive mesh error estimator.Source: §2.3 · pp. 6–8
• The condensed system is solved with unpreconditioned conjugate gradients; unconstrained magnetic-field degrees of freedom on Robin boundaries do not affect convergence for p≥1.The authors report that standard black-box preconditioners did not satisfactorily improve wall-clock time, leaving a formulation-specific preconditioner for future work.Source: §2.3 · p. 8

3. Multiphysics coupling

• The electromagnetic DPG solver is implemented within lethe and coupled with its existing heat-transfer and fluid-dynamics solvers; further implementation details for those solvers are referenced to Alphonius et al..Source: §3 · p. 9

3.1. Fluid dynamics

• The fluid solver in lethe solves the incompressible Newtonian Navier–Stokes equations, using SUPG and PSPG stabilization for momentum and continuity, respectively.Both stabilizations are implicit, with parameters computed from local flow conditions and mesh size; further solver and matrix-free implementation details are referenced to Prieto Saavedra et al..Source: §3.1 · p. 9

3.2. Heat transfer

• lethe solves enthalpy conservation for incompressible flow, with electromagnetic power dissipation providing the volumetric heat source.For time-harmonic fields, power dissipation is obtained from the time-averaged Poynting vector; SUPG stabilization mitigates instabilities in high-Péclet-number flows.Source: §3.2 · pp. 9–10

3.3. Time integration and coupling strategy

• Electromagnetic fields are solved initially and again when needed, while their time-averaged power dissipation is reused as a fixed heat source between solves.The strategy exploits the separation between nanosecond electromagnetic periods, millisecond fluid timescales, and second-scale heating; electromagnetic solves are therefore called orders of magnitude less often than fluid and heat solvers.Source: §3.3 · p. 10

4.1. Waveguide analytical verification

• Waveguide field errors converged at the expected p+1 rate in L2; p=3 at refinement level 3 matched p=1 at level 5 with about eight times fewer degrees of freedom.At about four degrees of freedom per wavelength, the under-resolved field showed amplitude decay and element-localized phase shifts.Source: §4.1 · pp. 12–13, Figure 4 · p. 14

4.2. Fichera oven

• The Fichera oven’s adaptive estimator refined near the re-entrant corner and edges, and global error decreased monotonically in agreement with Petrides and Demkowicz [33].The final solution was qualitatively similar to their study; final meshes differed, plausibly because their refinement procedure was unspecified.Source: §4.2 · pp. 14–15, Figure 6 · p. 15, Figure 7 · p. 16

4.3. Flow around a microwave-heated body

• The flow cases around a cylinder, square prism, and tilted square prism use SiC obstacles; the study also reproduces Peng et al.’s resonance case without flow.The flow simulations use Re=400; material properties are assumed temperature-independent because temperature-dependent data were unavailable.Source: §4.3 · pp. 15–17

4.3.1. No inlet velocity

• 1.1 K was the average temperature rise in the post after 60 seconds, agreeing with Peng et al.’s approximately 1 K experimental increase.The temperature profile followed the electromagnetic field pattern and retained its shape under the no-flow, adiabatic conditions.Source: §4.3.1 · p. 18, Figure 10 · p. 19

4.3.2. With inlet velocity

• With flow, obstacles mostly reflect the incident electric field, while transmitted energy heats them; the tilted prism has the most heterogeneous temperature distribution, by at most about one degree.The temperature gradient is more pronounced along the flow and wave-propagation direction, where flow removes heat more efficiently at the obstacle’s front.Source: §4.3.2 · pp. 19–20
• After 60 seconds, the square obstacles’ average temperature is about 20 K higher than the cylinder’s, and none reaches thermal equilibrium despite gains exceeding 200 K.The larger obstacles absorb energy throughout their volume because the electromagnetic penetration depth is comparable to their characteristic diameter.Source: §4.3.2 · p. 20
• The tilted square prism heats approximately 0.375 K s−1 faster than the cylinder, while the square prism’s advantage grows from approximately 0.2 to 0.5 K s−1 over time.Figure 12 shows the temperature profiles and differences; the changing gap reflects nonlinear interplay between fluid flow and electromagnetic heating, even with constant physical properties.Source: §4.3.2 · pp. 20–21, Figure 12 · p. 21
• Although the square and tilted prisms have similar temperatures at 60 seconds, the square prism may reach a higher temperature later.Source: §4.3.2 · p. 21

Conclusions

• The DPG Maxwell solver adapts meshes using its built-in error estimator and couples high-order electromagnetics, heat transfer, and incompressible fluid dynamics with independently chosen polynomial orders.The framework integrates these physics within Lethe.Source: §5 · p. 21
• For the rectangular waveguide, the solver matches analytical convergence orders at p=1, 2, and 3; higher orders reduce pollution, which appears as amplitude loss rather than complete accuracy loss.Source: §5 · p. 21
• In the Fichera oven benchmark, results agree with Carstensen et al. [7] and Petrides and Demkowicz [33], while adaptive refinement near singularities monotonically reduces global energy-norm error.Source: §5 · p. 21
• For Peng et al.’s 24 mm-radius alumina obstacle in a WR430 waveguide, the solver reproduces the electric-field amplitude pattern and average obstacle temperature increase after 60 seconds.Source: §5 · p. 21
• At Reynolds number 400, three obstacle simulations show geometry-dependent heating and flow, rapid high temperatures, and large-system Maxwell solves with a limited memory footprint using CG.The simulations also demonstrate independently chosen polynomial orders to optimize degrees of freedom.Source: §5 · p. 22
• Future work includes sum-factorization to lower computational cost, preconditioning to improve iterative-solver convergence [33] [6], and coupling to Lethe’s discrete element method for microwave-heated multiphase flows.Source: §5 · p. 22

A. Dimensionless conventions

• The time-harmonic electromagnetic solver is implemented in dimensionless form using specified scalings.Source: A. Dimensionless conventions · p. 22
• E0 denotes the reference electric field, Z0, ε0, μ0, and c0 denote vacuum impedance, permittivity, permeability, and light speed, and L is the reference length scale.The subscript “dim” identifies which variable has dimension when ambiguous.Source: A. Dimensionless conventions · p. 22

• The DPG trial and test polynomial spaces require careful selection; this section restates the broken and interface spaces and refers readers to Carstensen et al. [7] for deeper discussion.Source: B. Polynomial spaces · p. 22
• The discretized broken spaces are defined elementwise for the supported 3D implementation, while L2(Ωh) equals L2(Ω) and needs no breaking.Source: B. Polynomial spaces · p. 22
• Trace spaces are defined using linear maps and trace operators on the mesh skeleton, including element traces and restrictions to element boundaries.The element outward normal is also specified.Source: B. Polynomial spaces · p. 23
• Following the finite element exterior calculus convention, Q^-_pΛ^k(□^d) denotes the trimmed tensor-product degree-p space for differential k-forms on a d-dimensional hexahedral mesh.Source: B. Polynomial spaces · p. 23
• In 3D, Q^-_pΛ^1(□^3) gives first-kind Nédélec elements conforming in H(curl), Q^-_pΛ^2(□^3) Raviart–Thomas elements in H(div), and Q^-_pΛ^3(□^3) discontinuous elements in L2.Source: B. Polynomial spaces · p. 23

C.1. Waveguide

• Table 3 summarizes the parameters for the rectangular-waveguide verification simulation.The supplied passage does not include the table’s parameter values.Source: §C.1 · p. 24, Table 3 · p. 24
• The listed parameters include s = 0.968987646 and a linear-solver relative residual tolerance of 10^-8.Source: §C.1 · p. 24

C.2. Fichera oven

• Table 4 presents the solver parameters for the Fichera oven verification simulation.Source: §C.2 · p. 25, Table 4 · p. 25

C.3. Microwave-heated flow

• Tables 5 and 6 summarize parameters for the cylinder simulation without inlet velocity and for the fully coupled simulation with inlet velocity.FD, HT, and EM denote fluid dynamics, heat transfer, and electromagnetism.Source: §C.3 · p. 26, Table 5 · p. 26, Table 6 · p. 27
• The FD solver uses GCMG preconditioning, 500 maximum Krylov vectors, PH coarsening, three inverse-diagonal smoother iterations, a direct coarse-grid solver, and enabled eigenvalue estimation.The HT method is GMRES with a relative residual tolerance of 10^-4.Source: §C.3 · p. 27