Works out the stiffness, mass, heat and load matrices of each 10-node tetrahedron used by the stress, vibration, heat and buckling studies on a tetrahedral mesh.
How the precision is kept
All element arrays are Float64Array and the shape functions and rule points are 64-bit constants. Stiffness, conductivity, thermal load and geometric stiffness use a 4-point rule; mass, heat capacity and heat generation use a 64-point rule; boundary faces use a 6-point triangle rule. The file headers state the rules are exact on straight-sided elements up to their degree; on elements with curved sides the result is a quadrature approximation.
Calculations under this policy
18 calculations in 2 source files. The code is read from the released app when this page loads. It is the calculation people run, not a copy kept for this page.
Builds the integration points and weights used for mass and body loads.
x = u, y = v(1 − u), z = w(1 − u)(1 − v); weight = wₐ w_b w_c ÷ 8 × (1 − u)² (1 − v)
Built once from four Gauss–Legendre nodes and weights written to 16 digits. Each point and weight is a short product of 64-bit numbers; nothing is rounded or stored in 32-bit.
tetfea.js
lines 37 to 46
PartCraft 0.8.1
/** Collapsed (Duffy / Stroud conical product) Gauss–Legendre rule on the unit tet: n³ points [x, y, z, w]. */
export function tetRule(gl = GL4) {
const out = [];
for (const [a, wa] of gl) for (const [b, wb] of gl) for (const [c, wc] of gl) {
const u = (a + 1) / 2, v = (b + 1) / 2, w = (c + 1) / 2;
// x = u, y = v(1 − u), z = w(1 − u)(1 − v); Jacobian (1 − u)²(1 − v)
out.push([u, v * (1 - u), w * (1 - u) * (1 - v), wa * wb * wc / 8 * (1 - u) ** 2 * (1 - v)]);
}
return out;
}
Shape functions of the 10-node tetrahedron
64-bit
The ten shape functions and their derivatives at a point of the element.
corner: N = L(2L − 1); midside: N = 4 LᵢLⱼ
The fourth barycentric coordinate is 1 − x − y − z and every function is a product of two 64-bit numbers. The derivatives use the integer table of dL/dξ, so no rounded factor enters.
tetfea.js
lines 51 to 57
PartCraft 0.8.1
/** T10 shape functions and their reference derivatives at (x, y, z) → { N: [10], dN: [[3] × 10] }. */
export function t10Shape(x, y, z) {
const L = [1 - x - y - z, x, y, z], N = new Array(10), dN = new Array(10);
for (let i = 0; i < 4; i++) { N[i] = L[i] * (2 * L[i] - 1); dN[i] = DL[i].map((d) => (4 * L[i] - 1) * d); }
EDGES.forEach(([i, j], k) => { N[4 + k] = 4 * L[i] * L[j]; dN[4 + k] = [0, 1, 2].map((c) => 4 * (L[i] * DL[j][c] + L[j] * DL[i][c])); });
return { N, dN };
}
Jacobian and physical derivatives
64-bit
The element's Jacobian determinant at a point and the shape-function derivatives in x, y and z.
J = Σ xₐ ∂Nₐ/∂ξ; det J by cofactors; ∂N/∂x = J⁻ᵀ ∂N/∂ξ
Each Jacobian entry is a sum over the ten nodes of an absolute coordinate times a derivative, and those derivatives add to zero, so for an element that is small compared with its distance from the origin the sum subtracts nearly equal numbers. The inverse is the cofactor matrix times 1 ÷ det; this function does not test the determinant.
tetfea.js
lines 61 to 73
PartCraft 0.8.1
/** Physical derivatives at a point of element el → { dNx: Float64Array(30) [node·3 + axis], det }. */
export function jac(X, el, dN, out) {
let j00 = 0, j01 = 0, j02 = 0, j10 = 0, j11 = 0, j12 = 0, j20 = 0, j21 = 0, j22 = 0;
for (let a = 0; a < 10; a++) { const n = 3 * el[a], d = dN[a], x = X[n], y = X[n + 1], z = X[n + 2]; j00 += x * d[0]; j01 += x * d[1]; j02 += x * d[2]; j10 += y * d[0]; j11 += y * d[1]; j12 += y * d[2]; j20 += z * d[0]; j21 += z * d[1]; j22 += z * d[2]; }
const det = j00 * (j11 * j22 - j12 * j21) - j01 * (j10 * j22 - j12 * j20) + j02 * (j10 * j21 - j11 * j20);
const id = 1 / det;
// inverse transpose: dN/dx = J⁻ᵀ dN/dξ
const i00 = (j11 * j22 - j12 * j21) * id, i01 = (j02 * j21 - j01 * j22) * id, i02 = (j01 * j12 - j02 * j11) * id;
const i10 = (j12 * j20 - j10 * j22) * id, i11 = (j00 * j22 - j02 * j20) * id, i12 = (j02 * j10 - j00 * j12) * id;
const i20 = (j10 * j21 - j11 * j20) * id, i21 = (j01 * j20 - j00 * j21) * id, i22 = (j00 * j11 - j01 * j10) * id;
for (let a = 0; a < 10; a++) { const d = dN[a]; out[3 * a] = d[0] * i00 + d[1] * i10 + d[2] * i20; out[3 * a + 1] = d[0] * i01 + d[1] * i11 + d[2] * i21; out[3 * a + 2] = d[0] * i02 + d[1] * i12 + d[2] * i22; }
return det;
}
Element stiffness
64-bit
The 30 × 30 stiffness matrix of one element and its volume.
Kₑ = Σₚ wₚ det Jₚ E · Bᵀ D B; volume = Σₚ wₚ det Jₚ
Four Gauss points, weights 1 ÷ 24 each, accumulated in a Float64Array of 900 entries. A determinant that is not above 0 throws an error instead of producing a matrix. The volume is the sum of the same four weighted determinants.
tetfea.js
lines 75 to 107
PartCraft 0.8.1
/** Element stiffness (30×30, dofs node·3 + axis) for modulus E, D = elasticD(ν); fills Ke, returns the volume. */
export function t10Stiffness(X, el, E, D, Ke = new Float64Array(900)) {
Ke.fill(0);
const g = new Float64Array(30), DB = new Float64Array(6 * 30);
let vol = 0;
for (let p = 0; p < 4; p++) {
const det = jac(X, el, SH4[p].dN, g);
if (!(det > 0)) throw new Error(t('study.tet.inverted'));
const w = GP4[p][3] * det * E; vol += GP4[p][3] * det;
// DB = D·B, B rows (xx, yy, zz, xy, yz, zx) as fea.js
for (let a = 0; a < 10; a++) {
const dx = g[3 * a], dy = g[3 * a + 1], dz = g[3 * a + 2];
// columns 3a, 3a+1, 3a+2 of B: [dx,0,0,dy,0,dz], [0,dy,0,dx,dz,0], [0,0,dz,0,dy,dx]
for (let r = 0; r < 6; r++) {
const d = D.subarray ? D : D;
DB[r * 30 + 3 * a] = d[r * 6] * dx + d[r * 6 + 3] * dy + d[r * 6 + 5] * dz;
DB[r * 30 + 3 * a + 1] = d[r * 6 + 1] * dy + d[r * 6 + 3] * dx + d[r * 6 + 4] * dz;
DB[r * 30 + 3 * a + 2] = d[r * 6 + 2] * dz + d[r * 6 + 4] * dy + d[r * 6 + 5] * dx;
}
}
for (let a = 0; a < 10; a++) {
const dx = g[3 * a], dy = g[3 * a + 1], dz = g[3 * a + 2];
for (let j = 0; j < 30; j++) {
// Bᵀ rows for dofs 3a.. times DB column j
const c0 = DB[j], c1 = DB[30 + j], c2 = DB[60 + j], c3 = DB[90 + j], c4 = DB[120 + j], c5 = DB[150 + j];
Ke[(3 * a) * 30 + j] += w * (dx * c0 + dy * c3 + dz * c5);
Ke[(3 * a + 1) * 30 + j] += w * (dy * c1 + dx * c3 + dz * c4);
Ke[(3 * a + 2) * 30 + j] += w * (dz * c2 + dy * c4 + dx * c5);
}
}
}
return vol;
}
Element mass and shape-function integrals
64-bit
The consistent mass matrix of one element and ∫N dV for each node, used for gravity.
Mₑ = Σₚ ρ wₚ det Jₚ NₐN_b; ∫Nₐ dV = Σₚ wₚ det Jₚ Nₐ
64 points are added in a fixed order into Float64Array. No lumping and no rounding.
tetfea.js
lines 109 to 118
PartCraft 0.8.1
/** Consistent scalar mass (10×10, × I₃ per node pair) for density rho, and the shape-function integrals ∫N dV. */
export function t10Mass(X, el, rho, Me = new Float64Array(100), Ni = new Float64Array(10)) {
Me.fill(0); Ni.fill(0);
const g = new Float64Array(30);
for (let p = 0; p < RULE64.length; p++) {
const det = jac(X, el, SH64[p].dN, g), w = RULE64[p][3] * det, N = SH64[p].N;
for (let a = 0; a < 10; a++) { Ni[a] += w * N[a]; const wa = rho * w * N[a]; for (let b = 0; b < 10; b++) Me[a * 10 + b] += wa * N[b]; }
}
return Me;
}
Stress and strain at a Gauss point
64-bit
The six strains and six stresses at one Gauss point from the element's displacements.
ε = B u; σ = E · D (ε − α ΔT [1 1 1 0 0 0])
Strains are sums of derivative times displacement over the ten nodes. The thermal strain, interpolated from the nodal temperature rises, is subtracted from the three normal strains before the multiplication by D. All 64-bit.
tetfea.js
lines 120 to 132
PartCraft 0.8.1
/** Strain-displacement at Gauss point p: stresses σ = E·D·B·ue (6), strains ε = B·ue. */
export function gpStress(X, el, p, ue, E, D, sig, eps, g = new Float64Array(30), th = null) {
jac(X, el, SH4[p].dN, g);
let exx = 0, eyy = 0, ezz = 0, gxy = 0, gyz = 0, gzx = 0;
for (let a = 0; a < 10; a++) {
const dx = g[3 * a], dy = g[3 * a + 1], dz = g[3 * a + 2], ux = ue[3 * a], uy = ue[3 * a + 1], uz = ue[3 * a + 2];
exx += dx * ux; eyy += dy * uy; ezz += dz * uz; gxy += dy * ux + dx * uy; gyz += dz * uy + dy * uz; gzx += dz * ux + dx * uz;
}
eps[0] = exx; eps[1] = eyy; eps[2] = ezz; eps[3] = gxy; eps[4] = gyz; eps[5] = gzx;
// th { alpha, dT (10, the element's nodal rises) }: σ = D(ε − αΔT·[1 1 1 0 0 0]); eps stays the total strain
let et = 0; if (th) { const N = SH4[p].N; for (let a = 0; a < 10; a++) et += N[a] * th.dT[a]; et *= th.alpha; }
for (let r = 0; r < 6; r++) { let s = 0; for (let k = 0; k < 6; k++) s += D[r * 6 + k] * (k < 3 ? eps[k] - et : eps[k]); sig[r] = E * s; }
}
6-point rule on a triangle face
64-bit
The points and weights used to integrate loads, areas and heat exchange over boundary faces.
The two point coordinates and the two weights are written to 15 decimals. As written, the six weights add to 0.999999999999999, so areas and load totals integrated with this rule carry a relative error near 1e-15 from the constants alone.
tetfea.js
lines 423 to 424
PartCraft 0.8.1
// 6-point degree-4 triangle rule (Dunavant): [L0, L1, L2, w] (w sums to 1)
export const TRI6 = (() => { const a = 0.445948490915965, b = 0.091576213509771, wa = 0.223381589678011, wb = 0.109951743655322; return [[a, a, 1 - 2 * a, wa], [a, 1 - 2 * a, a, wa], [1 - 2 * a, a, a, wa], [b, b, 1 - 2 * b, wb], [b, 1 - 2 * b, b, wb], [1 - 2 * b, b, b, wb]]; })();
Shape functions of the 6-node face
64-bit
The six shape functions of a boundary face and their derivatives.
corner: N = L(2L − 1); midside: N = 4 LᵢLⱼ
Products of two 64-bit barycentric coordinates; L₀ is 1 − L₁ − L₂.
The area-weighted normal and position at each rule point of a curved boundary face.
n dA = (∂x/∂L₁ × ∂x/∂L₂) × w ÷ 2
The two tangents are sums of derivative times absolute node coordinate; their cross product times the weight and the reference area 1 ÷ 2 gives the area vector. Nothing is normalized here, so no division occurs.
tetfea.js
lines 433 to 443
PartCraft 0.8.1
/** Integrate over boundary T6 triangle t: fn(N, normal (outward, |n| = area density), point) per rule point. */
function faceIntegrate(model, t, fn) {
const X = model.X, v = model.btri6.subarray(6 * t, 6 * t + 6);
for (const [, L1, L2, w] of TRI6) {
const { N, d1, d2 } = t6(L1, L2), a = [0, 0, 0], b = [0, 0, 0], p = [0, 0, 0];
for (let k = 0; k < 6; k++) for (let c = 0; c < 3; c++) { const x = X[3 * v[k] + c]; a[c] += d1[k] * x; b[c] += d2[k] * x; p[c] += N[k] * x; }
// reference triangle area 1/2
const n = [(a[1] * b[2] - a[2] * b[1]) * w / 2, (a[2] * b[0] - a[0] * b[2]) * w / 2, (a[0] * b[1] - a[1] * b[0]) * w / 2];
fn(N, n, p, v);
}
}
The area weight at each rule point of a 6-node face, for heat exchange and heat loads.
dA = |∂x/∂L₁ × ∂x/∂L₂| × w ÷ 2
Same tangent sums as the structural face integral; the magnitude is taken with Math.hypot. Uses the 15-decimal triangle rule constants.
tet-elements.js
lines 30 to 37
PartCraft 0.8.1
/** ∫ over the T6 facet Xf (18 coordinates): fn(N, dA) per rule point. */
function facetRule(Xf, fn) {
for (const [, L1, L2, w] of TRI6) {
const { N, d1, d2 } = t6(L1, L2), a = [0, 0, 0], b = [0, 0, 0];
for (let k = 0; k < 6; k++) for (let c = 0; c < 3; c++) { a[c] += d1[k] * Xf[3 * k + c]; b[c] += d2[k] * Xf[3 * k + c]; }
fn(N, Math.hypot(a[1] * b[2] - a[2] * b[1], a[2] * b[0] - a[0] * b[2], a[0] * b[1] - a[1] * b[0]) * w / 2);
}
}
Element conductivity matrix
64-bit
The 10 × 10 heat conduction matrix of one element.
Kₑ = Σₚ wₚ det Jₚ k ∇Nₐ·∇N_b
Four Gauss points, summed into a Float64Array. The determinant is not tested here.
tet-elements.js
lines 41 to 48
PartCraft 0.8.1
conductivity(X, k, out) {
out.fill(0, 0, 100);
for (let p = 0; p < 4; p++) {
const det = jac(X, ID10, SH4[p].dN, g30), w = GP4[p][3] * det * k;
for (let a = 0; a < 10; a++) for (let b = 0; b < 10; b++) out[a * 10 + b] += w * (g30[3 * a] * g30[3 * b] + g30[3 * a + 1] * g30[3 * b + 1] + g30[3 * a + 2] * g30[3 * b + 2]);
}
return out;
},
Element heat capacity
64-bit
The heat capacity matrix of one element, or its lumped diagonal for transient runs.
M = Σₚ wₚ det Jₚ ρc NₐN_b; lumped: Mₐₐ × (ΣM ÷ Σ diagonal)
The consistent matrix uses the 64-point rule. The lumped form scales the diagonal by the ratio of the matrix total to the diagonal total, one division per element, so the lumped entries add to the element's ρc × volume up to rounding.
tet-elements.js
lines 49 to 56
PartCraft 0.8.1
capacity(X, rc, out, lumped) {
const M = new Float64Array(100);
for (let p = 0; p < SH64.length; p++) { const w = RULE64[p][3] * jac(X, ID10, SH64[p].dN, g30) * rc, N = SH64[p].N; for (let a = 0; a < 10; a++) for (let b = 0; b < 10; b++) M[a * 10 + b] += w * N[a] * N[b]; }
if (!lumped) { out.set(M); return out; }
let tot = 0, dsum = 0; for (let i = 0; i < 100; i++) tot += M[i]; for (let a = 0; a < 10; a++) dsum += M[a * 11];
for (let a = 0; a < 10; a++) out[a] = M[a * 11] * tot / dsum;
return out;
},
Element heat generation load
64-bit
The nodal heat input of a uniform internal heat source.
fₐ = Σₚ wₚ det Jₚ Q Nₐ
64-point rule, 64-bit sums in a fixed order.
tet-elements.js
lines 57 to 61
PartCraft 0.8.1
heatLoad(X, Q, out) {
out.fill(0, 0, 10);
for (let p = 0; p < SH64.length; p++) { const w = RULE64[p][3] * jac(X, ID10, SH64[p].dN, g30) * Q, N = SH64[p].N; for (let a = 0; a < 10; a++) out[a] += w * N[a]; }
return out;
},
Temperature gradient at the element center
64-bit
The temperature gradient used for the heat flux of an element.
∇T = Σ ∇Nₐ Tₐ at L = (¼, ¼, ¼, ¼)
One evaluation at the centroid; the flux of a quadratic element is therefore reported as a single value per element.
tet-elements.js
lines 62 to 67
PartCraft 0.8.1
gradient(X, Te, out) {
jac(X, ID10, C0.dN, g30);
out[0] = out[1] = out[2] = 0;
for (let a = 0; a < 10; a++) for (let i = 0; i < 3; i++) out[i] += g30[3 * a + i] * Te[a];
return out;
},
Element thermal load
64-bit
The nodal forces that a temperature rise produces in one element.
f = Σₚ wₚ det Jₚ E α ΔTₚ (D₀₀ + D₀₁ + D₀₂) ∇N
The temperature rise is interpolated to each of the four Gauss points. The sum of the first row of D replaces the full product with the thermal strain vector, which is valid for the isotropic D used here.
tet-elements.js
lines 68 to 80
PartCraft 0.8.1
thermalLoad(X, E, nu, alpha, dT, out) {
out.fill(0, 0, 30);
const D = elasticD(nu);
// D·[1 1 1 0 0 0] = (D00 + D01 + D02) on the normal rows (isotropic)
const dn = D[0] + D[1] + D[2];
for (let p = 0; p < 4; p++) {
const det = jac(X, ID10, SH4[p].dN, g30), N = SH4[p].N;
let t = 0; for (let a = 0; a < 10; a++) t += N[a] * dT[a];
const s = GP4[p][3] * det * E * alpha * t * dn;
for (let a = 0; a < 10; a++) { out[3 * a] += s * g30[3 * a]; out[3 * a + 1] += s * g30[3 * a + 1]; out[3 * a + 2] += s * g30[3 * a + 2]; }
}
return out;
},
Stresses at the four Gauss points
64-bit
The stresses used to build the geometric stiffness of a buckling study.
σ = E · D (B u − α ΔT [1 1 1 0 0 0])
The same strain sums as the static stress recovery, written out per component; 24 values per element in a Float64Array.
tet-elements.js
lines 81 to 92
PartCraft 0.8.1
stressPoints(X, E, nu, ue, out, th = null) {
const D = elasticD(nu);
for (let p = 0; p < 4; p++) {
const det = jac(X, ID10, SH4[p].dN, g30); void det;
let e0 = 0, e1 = 0, e2 = 0, e3 = 0, e4 = 0, e5 = 0;
for (let a = 0; a < 10; a++) { const dx = g30[3 * a], dy = g30[3 * a + 1], dz = g30[3 * a + 2], ux = ue[3 * a], uy = ue[3 * a + 1], uz = ue[3 * a + 2]; e0 += dx * ux; e1 += dy * uy; e2 += dz * uz; e3 += dy * ux + dx * uy; e4 += dz * uy + dy * uz; e5 += dz * ux + dx * uz; }
let et = 0; if (th) { const N = SH4[p].N; for (let a = 0; a < 10; a++) et += N[a] * th.dT[a]; et *= th.alpha; }
const eps = [e0 - et, e1 - et, e2 - et, e3, e4, e5];
for (let r = 0; r < 6; r++) { let s = 0; for (let k = 0; k < 6; k++) s += D[r * 6 + k] * eps[k]; out[6 * p + r] = E * s; }
}
return out;
},
Element geometric stiffness
64-bit
The stress stiffness of one element for linear buckling.
Gₐ_b = Σₚ wₚ det Jₚ ∇Nₐᵀ σₚ ∇N_b
Four Gauss points; the 3 × 3 stress matrix is rebuilt from the six stored components at each point. 64-bit sums.
tet-elements.js
lines 93 to 105
PartCraft 0.8.1
geometricStiffness(X, sig, out) {
out.fill(0, 0, 100);
for (let p = 0; p < 4; p++) {
const det = jac(X, ID10, SH4[p].dN, g30), w = GP4[p][3] * det, s = sig.subarray ? sig.subarray(6 * p, 6 * p + 6) : sig.slice(6 * p, 6 * p + 6);
const m = [s[0], s[3], s[5], s[3], s[1], s[4], s[5], s[4], s[2]];
for (let a = 0; a < 10; a++) {
const ga = [g30[3 * a], g30[3 * a + 1], g30[3 * a + 2]];
const t0 = m[0] * ga[0] + m[3] * ga[1] + m[6] * ga[2], t1 = m[1] * ga[0] + m[4] * ga[1] + m[7] * ga[2], t2 = m[2] * ga[0] + m[5] * ga[1] + m[8] * ga[2];
for (let b = 0; b < 10; b++) out[a * 10 + b] += w * (t0 * g30[3 * b] + t1 * g30[3 * b + 1] + t2 * g30[3 * b + 2]);
}
}
return out;
},
Face areas and volume of the tetrahedral mesh
64-bit
The area of every boundary face and the mesh volume handed to the heat and buckling routines.
area = Σ ∫Nₐ dA; volume = Σₑ Σₚ wₚ det Jₚ
Each face area is the sum of the six face load integrals, held in a Float64Array. The volume is added element by element in mesh order with the 4-point rule.
tet-elements.js
lines 114 to 127
PartCraft 0.8.1
/** The MESH INTERFACE on a tetfea.js model (makeTetModel). */
export function tetMeshAdapter(model) {
const nb = model.btri6.length / 6, area = new Float64Array(nb), Xf = new Float64Array(18), L = new Float64Array(6);
const elemX = (e, X) => { for (let a = 0; a < 10; a++) { const n = model.tets10[10 * e + a]; X[3 * a] = model.X[3 * n]; X[3 * a + 1] = model.X[3 * n + 1]; X[3 * a + 2] = model.X[3 * n + 2]; } return X; };
for (let f = 0; f < nb; f++) { for (let c = 0; c < 6; c++) { const n = model.btri6[6 * f + c]; Xf[3 * c] = model.X[3 * n]; Xf[3 * c + 1] = model.X[3 * n + 1]; Xf[3 * c + 2] = model.X[3 * n + 2]; } TET10.face.load(Xf, L); area[f] = L.reduce((s, v) => s + v, 0); }
let volume = 0; const X = new Float64Array(30);
for (let e = 0; e < model.tets10.length / 10; e++) { elemX(e, X); for (let p = 0; p < 4; p++) volume += GP4[p][3] * jac(X, ID10, SH4[p].dN, g30); }
return {
type: 'tet10', element: TET10, nn: model.n, xyz: model.X, active: new Uint8Array(model.n).fill(1), ne: model.tets10.length / 10, npe: 10, conn: model.tets10, volume,
facets: { npf: 6, nodes: model.btri6, area },
facetsOn(faceIdx) { const want = new Set(faceIdx), out = []; for (let f = 0; f < nb; f++) if (want.has(model.bface[f])) out.push(f); return Int32Array.from(out); },
elemX,
};
}
In detail
Method
10-node isoparametric tetrahedron: corner functions L(2L − 1), midside functions 4 LᵢLⱼ. 4-point Gauss rule (degree 2) for stiffness; 64-point collapsed Gauss–Legendre rule, 4 points per direction, for mass and body loads; 6-point degree-4 triangle rule, named Dunavant in the code, on the 6-node faces. The 3 × 3 Jacobian is inverted by cofactors.
Number types
64-bit numbers; Float64Array for coordinates and element matrices; Int32Array for node numbers.
Tolerances
The stiffness routine stops with an error when the Jacobian determinant at a Gauss point is not above 0. There is no other tolerance.
Limits
The 4-point rule's coordinates and the Gauss–Legendre nodes are written to 16 digits; the triangle rule's points and weights are written to 15 decimals. The Jacobian is built from absolute node coordinates, so a small element far from the origin loses digits in it.
Approach: 64-bit, stated limit
The arithmetic is 64-bit, but the method or the data sets the accuracy: a solver tolerance, a mesh, a step size, a datum shift, a table. The limit is stated where the result is shown.