Thermal FEM Solver — MEMS Bolometer Pixel

Browser-based finite element analysis — static and transient heat transfer — SiN membrane + VOx + NiCr legs

ThermalFEM calculator screenshot showing 3D mesh with temperature distribution and transient response graph

Fig. 1 ThermalFEM interface: geometry and material inputs (left panel), 3D temperature distribution after solve (upper right), and transient temperature response at the pixel center (lower right).

Overview

This tool performs static and transient thermal finite element analysis (FEA) of a MEMS bolometer pixel structure directly in the browser, with no installation required. It is intended for rapid estimation of thermal conductance G, thermal capacitance C, and the thermal time constant τ = C/G of a suspended membrane pixel.

The geometry represents a square pixel membrane (block A) supported by four legs (blocks B–E), with a thin VOx sensing film (F) on top of the membrane and NiCr interconnect films (G–J) on the legs. Heat is generated uniformly over the membrane+film volume and conducted out through the legs to fixed-temperature anchors at both ends.

Device Structure

The modeled structure consists of ten rectangular blocks defined in Cartesian coordinates, all dimensions in micrometers:

Block	Material	Role
A	PECVD SiN	Square pixel membrane (substrate layer)
B, C	PECVD SiN	Bottom legs connecting membrane to anchors
D, E	PECVD SiN	Top legs connecting membrane to anchors
F	VOx (sputtered)	Sensing film on membrane surface
G, H, I, J	NiCr (sputtered)	Interconnect/resistor film on legs

Boundary conditions: all nodes at x = −d and x = l+d (the anchor ends of both leg pairs) are fixed at 20 °C. No convection or radiation is included.

Finite Element Model

Element type and mesh

The solver uses 8-node hexahedral (Hex8) elements with full 2×2×2 Gauss quadrature for both stiffness and mass matrix assembly. Each rectangular block is meshed with 10 divisions along the in-plane directions (x and y for the square, length direction for the legs) and 2 divisions through the thickness (z), giving a total of 2,000 elements and approximately 2,800 nodes for the default geometry. Adjacent blocks share nodes at their common faces (conforming mesh), ensuring continuity of temperature across material interfaces.

Static analysis

The steady-state heat equation is solved:

∇ · (k ∇T) + Q = 0

The weak form yields the global stiffness system K u = f, where u is the nodal temperature rise above the boundary condition temperature (20 °C), and f is the heat source vector. Heat is distributed uniformly by volume over the membrane region (blocks A and F combined).

The linear system is solved with a Jacobi-preconditioned conjugate gradient (PCG) method. Dirichlet boundary conditions are applied by row/column elimination, and the system is formulated in terms of temperature rise (shifted coordinates) so that numerical precision is maintained regardless of the absolute boundary temperature value.

All computations use P = 1 W internally. The displayed temperatures are then scaled by the user-specified total input power, exploiting the linearity of the problem. Thermal conductance is reported as:

G = 1 / u_max [W/K]

where u_max is the maximum temperature rise per unit watt.

Transient analysis

The time-dependent heat equation is integrated using the implicit (backward) Euler method:

(M/Δt + K) u_n+1 = (M/Δt) u_n + f

where M is the lumped thermal mass matrix assembled from nodal heat capacity contributions (ρ c_p per material). Implicit Euler is unconditionally stable, so time steps are not constrained by a stability limit; only accuracy determines the step size.

The solver runs in two phases. In the first phase it takes 100 steps of 1 ms each to estimate the thermal time constant τ. In the second phase it refines the time step to τ/40 and integrates out to 6τ, producing a smooth time–temperature curve. A new PCG solve is performed at each time step, but the system matrix (M/Δt + K) is assembled once per phase and reused.

The thermal time constant, thermal capacitance, and steady-state temperature are extracted from the response at the probe point (l/2, l/2, t_s), i.e. the center of the top surface of the membrane:

τ = time at which T reaches T_∞ × (1 − 1/e) ≈ 63.2%
C = τ × G [J/K]

Numerical note: All floating-point arithmetic uses 64-bit (double precision) throughout. The total input power is distributed as a volumetric heat source Q [W/μm³] = P / V_source computed by Gauss integration over source elements, not by element counting, to avoid rounding errors.

Interface

Left panel — inputs

Parameters are organized into four sections. Geometry values are in micrometers. Material properties (k, ρ, c_p) are entered per material; hovering over each input field shows the property name. Default values correspond to typical MEMS bolometer fabrication: PECVD SiN, sputtered VOx sensing film, and sputtered NiCr for the leg resistors.

Section	Description
Geometry	Pixel size l, leg width w and length d, substrate thickness t_s, VOx film thickness t_v. All in μm.
Material Properties	Thermal conductivity k [W/(m·K)], density ρ [kg/m³], specific heat c_p [J/(kg·K)] for each of the three materials. k affects both static and transient results; ρ and c_p affect only the transient analysis.
Thermal Load	Total input power P [W]. The static and transient results scale linearly with this value.
Mesh Statistics	Shown after Generate Mesh: node count, element count, DOF, estimated K matrix non-zeros, and current JS heap usage.

Buttons

Generate Mesh — builds the hex mesh and displays it in 3D. Must be run before any solve. Re-run whenever geometry changes.
Static Solve — solves the steady-state problem and shows the temperature distribution on the 3D mesh, Max Temperature, and G.
Transient Solve — runs the two-phase implicit Euler integration. A progress bar at the bottom of the 3D viewport shows completion. On completion, the 3D view updates to the steady-state temperature, and the time–temperature graph is shown below.

3D viewport (upper right)

Displays the mesh with element edges (black lines) and face colors indicating material (before solving) or temperature (after solving, jet colormap). The color scale bar on the right shows the temperature range. Max Temperature and thermal conductance G are shown in the top-left overlay.

Navigation: drag to rotate, scroll to zoom, Shift+drag to pan.

Transient graph (lower right)

Shows T vs. time at the membrane center. A red dashed vertical line marks τ, a yellow dashed horizontal line marks T_∞. The header bar reports τ [ms], C [J/K], T_∞ [°C], and the total simulated time t_end [ms].

Typical Results (Default Parameters)

Quantity	Value	Note
G	~ 2–5 × 10⁻⁷ W/K	Dominated by SiN and NiCr leg conductance
C	~ 3 × 10⁻⁹ J/K	Dominated by SiN membrane volume
τ = C/G	~ 6–15 ms
Mesh	2000 elements, ~2800 nodes	Default geometry
Memory	~10 MB JS heap	Single browser tab

Validation of Results

Confirmed that key simulation results are within 10% agreement with ANSYS benchmark data.
Initial parameters for the current pixel structure were adopted from the literature cited below. Using standard material properties, the calculated Heat Capacity (C) and Thermal Conductance (G) values show less than 5% deviation from the referenced study.
R. Andrew Wood "Uncooled thermal imaging with monolithic silicon focal planes", Proc. SPIE 2020, Infrared Technology XIX, (1 November 1993); https://doi.org/10.1117/12.160553

Limitations

Linear heat conduction only; temperature-dependent material properties are not supported.
No convection, radiation, or contact resistance.
Mesh density is fixed (10×10×2 per block); geometry changes require re-meshing but do not change the element count.
Single-threaded JavaScript execution; transient analysis of the default mesh takes roughly 10–60 seconds depending on the browser and hardware.

▶ Open Calculator