FESADev/docs/linear-static-mitc4-shell/numerical-review.md

# Linear Static MITC4 Shell Numerical Review

## 1. Metadata

- feature_id: `linear-static-mitc4-shell`
- source_formulation: `docs/formulations/mitc4-shell-formulation.md`
- source_requirements: `docs/requirements/linear-static-mitc4-shell.md`
- source_research: `docs/research/linear-static-mitc4-shell-research.md`
- source_io_definition: `docs/io-definitions/linear-static-mitc4-shell-io.md`
- source_reference_inventory: `docs/reference-models/linear-static-mitc4-shell-reference-models.md`
- repository_policy: `AGENTS.md`, `docs/SOLVER_AGENT_DESIGN.md`,
  `docs/numerical-reviews/README.md`
- reviewed_head: `cf769aa` (`mathematical implementation baseline`)
- prior_pass_commit: `60b42f4` (`context-only; verdict not inherited`)
- status: `pass-for-implementation-planning`
- owner_agent: `numerical-review-agent`
- date: `2026-08-13`
- implementation_planning_authorized: `true`
- implementation_complete: `false`
- build_test_complete: `false`
- reference_comparison_complete: `false`
- physics_evaluation_complete: `false`
- release_ready: `false`

이번 재검토는 현재 HEAD의 요구조건, 연구, 정식화, I/O 및 reference-case 계약을
처음부터 상호 대조했다. 기존 review의 판정과 artifact 관찰 결과는 결론의 전제로
사용하지 않았고, 이전 finding은 현 문서의 수식으로 다시 검산한 뒤 disposition만
기록했다. 원 MITC4 local paper는 tying 위치와 covariant shear 보간을 확인하는 데
read-only로 사용했다.

이 단계에서는 Abaqus, Harness, C++ build/test, FESA 실행 및 reference comparison을
수행하지 않았다. Reference artifact를 생성, 수정, 복원 또는 정규화하지 않았다.

## 2. Review Verdict

- verdict: `pass-for-implementation-planning`
- critical_blockers: `none`
- confirmed_defects: `none in the approved current linear-static scope`
- open_blocking_questions: `none`
- reason: 현재 정식화는 24 global DOF와 20 physical DOF의 관계, MITC4 shear
  tying, plane-stress section law, 공통 `2 x 2 x 2` quadrature, residual/stiffness,
  고정 drilling 안정화, 물리 recovery 및 검증 불변식을 구현계획으로 옮길 수 있을
  만큼 명시한다. 요구조건, I/O 및 reference 계약과 모순되는 차원, 부호, 위치 또는
  pass/fail 의미도 발견되지 않았다.
- downstream_boundary: 이 판정은 Implementation Planning 진입만 허용한다. 구현,
  MSVC build/CTest, reference comparison, physics sanity 또는 release를 승인하지 않는다.

정식화 Section 15의 geometrically nonlinear residual/tangent는 future-only다. 완전한
`Phi: R24 -> R20`, map Hessian, objective drilling potential 및 finite-rotation load work가
미정인 사실은 미래 nonlinear 구현을 막지만 현재 linear-static 판정은 막지 않는다.

## 3. Critical Findings

### 3.1 Confirmed defects

현재 승인된 선형 정적 범위에서 구현계획 전에 Formulation 또는 Research로 돌려보낼
confirmed mathematical defect는 없다.

`K20`의 exact-arithmetic 대칭/positive-semidefinite 구조와 20-to-24 congruence는
일관된다. 다만 실제 구현의 rank, rigid action, patch field와 reference error는 문서
검토만으로 통과했다고 볼 수 없으며 Section 6의 downstream test evidence가 필요하다.

### 3.2 Previous finding disposition

| previous item | rerun disposition | current independent basis |
| --- | --- | --- |
| `NR-C01` Jacobian/geometry inventory | `resolved` | Formulation 9.2-9.3은 center, eight stiffness points, four tying points 및 committed recovery points를 공통 fail-closed inventory로 두고 finite bases, nonzero surface measure와 `J>0`를 요구한다. 승인 범위는 calibrated distortion/warp cutoff를 요구하지 않는다. |
| `NR-C02` drilling normalization | `resolved` | Formulation 12.2는 `R+`를 오직 8 physical tangent-rotation diagonals의 finite positive 값으로 제한하므로 모든 후보의 단위가 `force*length`로 같다. |
| `NR-C03` mixed-DOF spectrum scaling | `resolved` | Formulation 12.5의 `(L_e I3,I2)` 및 `(L_e I3,I3)` congruence는 rank/condition evidence에서 translation/rotation 단위 혼합을 제거한다. Raw mixed-unit spectrum은 금지된다. |
| `NR-C04` 20/24 weak-form mismatch | `resolved` | Formulation 5.2-5.3과 7.1-7.2는 physical, drilling, external work를 모두 `V24`에서 `T_p^T`와 `T_d^T`로 결합한다. |
| `NR-C05` nonlinear 20-to-24 closure | `resolved for current scope` | Section 15는 physical chart tangent와 conditional global pullback을 분리하고 map-curvature 항을 보존하며, 미정인 global map/objective drill을 future-only blocker로 명시한다. |
| `NR-D01` drilling-direction moment | `retained and consistent` | Exact-zero moment는 별도 처리하고 nonzero moment에 `rho_M=abs(d dot M)/norm(M)<=1e-12`를 적용한다. Numerical drilling은 거부된 moment를 운반하지 않는다. |
| `NR-D02` normalized algebraic checks | `retained and consistent` | `1e-12` symmetry/frame/energy와 `1e-10` rigid/residual/equilibrium 기준은 scaled matrices와 unclamped denominators에 적용된다. |
| `NR-O01` coefficient sweep/plateau | `closed by product decision` | `k_d=1e-3 min(R+)`가 고정 계약이다. Sweep, plateau 및 coefficient optimality는 구현 gate가 아니다. |
| `NR-O02` drilling-energy ratio | `removed from scope` | Drilling energy는 내부 quadratic identity일 뿐 physical energy나 mandatory output이 아니며 ratio/warning threshold도 요구하지 않는다. |
| `NR-O03` smooth-director calibration | `removed from scope` | Pairwise positive incident-normal orientation, finite/nonzero averaging 및 duplicate-node fold modeling이 승인된 exact predicate다. 별도 angle calibration은 gate가 아니다. |
| `NR-O04` distortion/warp calibration | `removed from scope` | Basic topology, finite/nonzero surface measure 및 required-point `J>0`가 승인된 predicate다. Quality sweep이나 cutoff는 gate가 아니다. |
| `NR-O05` U/UR tolerance | `resolved` | 모든 관련 문서가 sole S4 reference에 고정 절대오차 `1.0e-5`, U blocking, UR warning-only를 동일하게 정의한다. Reference scale은 판정에 사용하지 않으며 S4R은 reference gate가 아닌 common-path evidence다. |

이전의 `needs-reference-model` 판정에 포함됐던 canonical naming, README,
`metadata.json`, provenance, expanded portfolio 및 아직 없는 comparison result는 현재
프로젝트 정책상 formulation verdict의 blocker가 아니다. 현 Reference Model 문서는
정확한 기존 input/displacement path와 row/tolerance 계약을 제공한다.

### 3.3 Open questions

- current_linear_scope: `none blocking`
- future_geometric_nonlinearity: finite global rotation coordinate, `Phi`와 그 1/2차
  미분, chart recentering, objective drilling, nodal-moment work 및 nonlinear output/state
  계약이 미정이다. 이는 별도 future formulation/review가 소유한다.
- optional_characterization: near-singular positive-J geometry의 conditioning과 original
  MITC4의 distorted-curved membrane locking을 더 넓게 정량화할 수 있으나 현재 승인된
  planning/completion gate는 아니다.
- downstream_results: implementation rank/patch evidence, S4 reference comparison과
  S4R common-path test 결과는 해당 후속 Agent가 판정한다. 부재 자체는
  pre-implementation review의 결함이 아니다.

## 4. Numerical Risk Assessment

| risk label | assessment | required in-scope control |
| --- | --- | --- |
| `rigid_body_modes` | Physical `K20`은 six rigid modes와 expected rank 14를 가져야 한다. 24-DOF embedding은 네 drill null coordinates를 더하고 fixed drill block 뒤 expected rank 18/nullity 6이다. | 세 translation과 세 rotation을 명시적으로 구성한다. Rotation mode는 `u_I=omega x X_I`, `theta_I=omega-(omega dot d_I)d_I`, `gamma_I=0`를 사용한다. |
| `patch_test` | Bilinear membrane/bending field와 MITC tied shear는 required patch states를 표현할 계약을 갖는다. | `E11/E22/G12`, `K11/K22/K12`, `G13/G23`를 독립 시험하고 signs/component order/resultants/stress를 함께 확인한다. |
| `symmetry` | `B^T C B`, `T_p^T K20 T_p`, `T_d^T(k_d I)T_d`는 exact arithmetic에서 symmetric이다. | Scaled Frobenius check `<=1e-12`; deterministic assembly가 대칭을 깨지 않는지 확인한다. |
| `positive_definiteness` | Free element는 six-mode semidefinite이고, 충분히 구속된 nonsingular `Kff`는 positive definite가 기대된다. Geometry 또는 supports가 부적절하면 singularity가 정당하다. | Scaled spectrum/rank, non-rigid positive physical energy, constrained solve 및 singular negative cases를 분리한다. |
| `hourglass` | `1 x 1` reduced integration을 쓰지 않으므로 Abaqus-style hourglass path는 `N/A`다. | Full `2 x 2` midsurface rank test는 유지한다. S4R source label로 reduced rule을 선택하지 않는다. |
| `shear_locking` | Edge-midpoint MITC projection이 transverse-shear locking을 다루지만 모든 mesh/thickness에서 완전 제거를 주장할 수 없다. | Required shear/bending patch와 declared S4 reference를 통과한다. Broader thin/thick convergence는 nonblocking characterization이다. |
| `membrane_locking` | Original MITC4는 membrane strain을 수정하지 않아 distorted curved meshes에서 알려진 locking 위험이 남는다. | Known limitation을 유지하고 MITC4+ 성능을 주장하지 않는다. Expanded curved/distorted portfolio는 optional이다. |
| `volumetric_locking` | 승인된 homogeneous plane-stress shell에는 `N/A`다. | `C5`를 3D nearly-incompressible law로 확장하지 않는다. |
| `distortion` | 양의 `J`를 유지하는 심한 distortion/warpage는 정확도와 rank/conditioning을 악화할 수 있다. | 모든 required location의 exact predicates와 rank/finite-result checks를 시행한다. 승인되지 않은 quality cutoff를 추가하지 않는다. |
| `singular_jacobian` | Nonfinite/nonpositive `J`, zero surface measure, invalid reciprocal basis는 mapping을 무효화한다. | Center, stiffness, tying 및 recovery inventory 전체를 omission/clamp 없이 fail closed한다. |
| `conditioning` | Thin shells, near-degenerate positive-J geometry 및 작은 fixed drill scale에서 `Kff` conditioning이 나빠질 수 있다. | Spectrum/condition evidence는 오직 length-scaled matrix를 사용한다. Threshold calibration은 gate가 아니지만 factorization failure는 결정적으로 진단한다. |
| `convergence` | Current solve는 direct linear solve라 Newton convergence는 `N/A`; spatial convergence와 locking trend는 모델 의존이다. | Free residual/global equilibrium `<=1e-10`과 declared S4 reference case를 확인한다. Broader mesh sequences는 optional이다. |
| `drilling_contamination` | Numerical drill이 physical strain/recovery에 들어가면 비물리 결과가 생긴다. | `T_p`/`T_d`를 분리하고 pure drill에서 physical strain/resultant/stress/energy가 zero임을 시험한다. Full residual은 의도대로 stabilized system 전체를 포함한다. |
| `future_nonlinear_misuse` | Section 15만으로 global nonlinear element를 만들면 nonobjective drill 또는 inconsistent tangent가 된다. | Current plan에서 완전히 제외하고 별도 승인 전 실행하지 않는다. |

## 5. Consistency Checks

### 5.1 Units, dimensions, DOF order, and constrained/free system — `pass`

- Per-node global order is exactly `[UX,UY,UZ,URX,URY,URZ]`; element order is 24
  global coordinates and 20 physical coordinates
  `[uX,uY,uZ,alpha,beta]` per node plus four separately selected `gamma` coordinates.
- `T_p` is `20 x 24`, `T_d` is `4 x 24`, `K20` is `20 x 20`, and both global
  stiffness contributions are `24 x 24`.
- Translation-translation, translation-rotation, and rotation-rotation stiffness
  blocks have units `force/length`, `force`, and `force*length`; `R+` therefore
  excludes every translational diagonal.
- The constrained/free equation is `Kff*df=Ff-Kfc*dc`. Stiffness partition and
  factorization precede load assembly, and an all-constrained valid `0 x 0 Kff` is
  not reclassified as singular.
- `r=K*d-F` fixes the internal-minus-external sign. Constrained entries are the
  required reaction rows and free entries remain residual evidence.

### 5.2 Local/global transforms, congruence, and energy — `pass`

- `(a_I,b_I,d_I)` and `(e1,e2,e3)` are deterministic right-handed orthonormal frames.
  The least-aligned-axis nodal rule avoids a fixed-axis parallel singularity.
- `[alpha,beta,gamma]^T=R_I^T theta_I^g` gives
  `delta d=beta*a-alpha*b=theta x d` with the required sign.
- `q20=T_p qg` and `gamma=T_d qg` preserve virtual work. Congruence gives
  `Kphys24=T_p^T K20 T_p` and `Kdrill24=T_d^T(k_d I4)T_d`; the corresponding local and
  global quadratic energies are identical.
- A physical rigid rotation uses only the tangent projection of `omega`, so `gamma=0`
  and drilling does not destroy the six physical rigid modes.

### 5.3 Kinematic operators and MITC tying — `pass`

- Bilinear `N_I` satisfies partition, Kronecker and derivative-sum identities.
- Membrane and bending content comes from the direct covariant small-strain operator.
  Only `epsilon_xi-zeta` and `epsilon_eta-zeta` are replaced.
- `epsilon_xi-zeta` is tied at `(0,-1,0)` and `(0,+1,0)` and interpolated in `eta`;
  `epsilon_eta-zeta` is tied at `(-1,0,0)` and `(+1,0,0)` and interpolated in `xi`.
  Each interpolation reproduces its own edge value and is constant along the edge
  direction, matching the original MITC4 construction.
- The assumed covariant tensor is reconstructed through reciprocal bases, projected
  into the stored local Cartesian frame, and converted once to engineering shear
  `gamma_ij=2 epsilon_ij`. The same projected `B_bar` drives strain, internal force,
  stiffness and recovery.

### 5.4 Constitutive and section matrices — `pass`

- `Cps=E/(1-nu^2)[[1,nu,0],[nu,1,0],[0,0,(1-nu)/2]]` uses engineering `G12`; its
  shear coefficient is exactly `G=E/[2(1+nu)]`.
- `C5=diag(Cps,(5/6)G I2)` is symmetric positive definite for `E>0` and
  `-1<nu<0.5`. `sigma33=0` and absent thickness stretch remain assumptions.
- `A=t Cps`, `B=0`, `D=t^3 Cps/12`, and `As=(5/6)Gt I2` have consistent dimensions.
  Membrane/shear strains are dimensionless, curvature is `1/length`, `N/Q` is
  `force/length`, `M` is `force`, stress is `force/length^2`, and energy is
  `force*length`.

### 5.5 Jacobian, derivative transform, and integration — `pass`

- The three-dimensional degenerated mapping uses
  `X=sum(N X_I)+(t*zeta/2)sum(N d_I)` with a unit nodal director and separate
  thickness, preventing nodal thickness double counting.
- `J=det[G_xi,G_eta,G_zeta]` and finite covariant/reciprocal bases are checked at all
  eight stiffness points, all four midsurface tying points, center, and every
  committed bottom/middle/top recovery evaluation. Failed points are not skipped,
  averaged, clamped or repaired.
- Direct natural derivatives are converted covariantly and then to the local
  Cartesian tensor through contravariant bases; no flat-element derivative shortcut
  is substituted for curved/warped accepted geometry.
- Both S4 and S4R use the common in-plane `2 x 2` points
  `+-1/sqrt(3)` with unit weights and two identical thickness points, for eight
  volume evaluations. Tied shear is evaluated at `zeta=0` and reused at both
  thickness points while the remaining mapping and direct components use the actual
  thickness point.

### 5.6 Internal force, residual, stiffness, and future tangent — `pass`

- Current `K20=integral(B_bar^T C5 B_bar dV)` and `f_int20=K20 q20` are mutually
  consistent and symmetric positive semidefinite in exact arithmetic.
- The complete current weak form is in `V24` and adds the numerical drilling
  gradient before subtracting the global nodal `CLOAD` vector.
- No geometric stiffness or nonlinear state enters the current product path.
  Future Section 15 correctly separates `Kmat` and the stress-dependent `Kgeo` and,
  conditionally on a future `Phi`, includes both `A^T K20 A` and the
  residual-weighted map-Hessian term.

### 5.7 Fixed drilling contract — `pass`

- `R+` contains only finite strictly positive diagonals of the eight physical
  tangent-rotation coordinates. `k_ref=min(R+)`, `k_d=1e-3 k_ref`, and
  `Kd_local=k_d I4` are dimensionally consistent and deterministic; empty `R+`
  fails validation.
- `T_d^T(k_d I4)T_d` is symmetric and positive on the four pure drilling
  coordinates. It must remove those four nonphysical null modes without changing
  the physical rank/null modes.
- Drilling is excluded from generalized strain/resultant/stress and reported
  physical strain energy. A director-parallel applied nodal moment is rejected as
  `unsupported-drilling-load`; no numerical drill load channel exists.

### 5.8 Recovery, signs, locations, units, and external comparison — `pass`

- Nodal `[U1,U2,U3,UR1,UR2,UR3]` and full-residual
  `[RF1,RF2,RF3,RM1,RM2,RM3]` are global and source-node ordered.
- Four midsurface Gauss rows recover
  `[E11,E22,G12,K11,K22,K12,G13,G23]` and
  `[N11,N22,N12,M11,M22,M12,Q13,Q23]` in the stored local frame. The definition
  `e_m(z)=epsilon0+z*kappa` fixes curvature, moment, and bottom/top stress signs.
- `[S11,S22,S12]` is evaluated directly at `zeta=-1,0,+1`; `S33=0` is documented but
  not emitted, and `S13/S23` point stress is not synthesized. Different natural or
  section locations are never averaged.
- Reference comparison first rejects missing, extra, duplicate, nonfinite or
  identity-mismatched rows. Every U/UR row then uses fixed absolute tolerance
  `1.0e-5`; no reference-scale decision term, zero clamp or row denominator is
  introduced. U1/U2/U3 is blocking and UR1/UR2/UR3 is warning-only.
- Source S4 and S4R select the same FESA MITC4 kernel/quadrature/recovery path while
  preserving source type. This is an input mapping, not an Abaqus formulation,
  integration, stabilization or recovery equivalence claim.

## 6. Verification Readiness

### 6.1 Downstream unit and invariant tests

Implementation Planning shall convert the following to `RED -> GREEN -> VERIFY`:

1. Shape identities; nodal/integration frame orthonormality, handedness and axis
   tie-break determinism.
2. `T_p`/`T_d` dimensions, orthogonal channel selection, virtual-work equality and
   nonzero transformation-energy equality.
3. Hand-calculated direct membrane/bending columns, all four covariant tying values,
   interpolation weights and engineering-shear factors.
4. `Cps/C5/A/D/As` coefficients, symmetry, positive definiteness, dimensions and
   force/length unit-rescaling invariance.
5. Common `2 x 2 x 2` point/weight order and an independent analytical or
   higher-order flat-element stiffness/recovery cross-check.
6. Required-location geometry validation: valid planar/rotated/warped cases and
   duplicate, bow-tie/self-intersecting, zero-area, reversed, nonfinite,
   nonpositive-J and opposed-normal negative cases.
7. Scaled symmetry `<=1e-12`, physical rigid action `<=1e-10`, expected physical
   rank 14, stabilized rank 18/nullity six, and positive non-rigid physical energy.
8. Exact `R+` membership, exclusion of translations, fixed coefficient, empty-`R+`
   failure, pure drill action and zero physical recovery/energy.
9. Stable COO/reduction, source/result/diagnostic order and thread-count
   repeatability.
10. `Kff/Kfc` effective RHS, nonzero prescribed values, full-residual reaction,
    singular-support negative case, and valid all-constrained `0 x 0 Kff` case.
11. Exact-zero and accepted/rejected `rho_M` moment projections, including proof
    that rejected drilling moments never reach stabilization.
12. Mandatory HDF5 locations/components/units, physical-only energy, nonfinite
    recovery failure and atomic finalization.

### 6.2 Required patch and sign tests

- independent constant `E11`, `E22`, and `G12` membrane fields with `N` and
  middle-stress signs;
- pure `K11` and `K22` bending with `M` order and bottom/top stress reversal;
- pure `K12` twist with `M12` sign;
- constant `G13` and `G23` transverse shear with `Q13/Q23` order;
- six physical rigid states and four pure drilling states;
- source-type-only S4/S4R variants producing identical FESA numeric rows and
  different preserved source metadata.

### 6.3 Reference and physics handoff readiness

The sole declared read-only acceptance pair is:

- `reference/shell/shell.inp` with
  `reference/shell/shell displacements.csv` for full-integration S4.

S4R remains required input coverage through source mapping, common-kernel,
deterministic-assembly and HDF5 metadata tests, but no S4R Abaqus artifact is
consumed by reference acceptance.

The Reference Model and I/O documents define deterministic HDF5-to-CSV identity,
precheck and tolerance sufficiently for later comparison. This review did not assert
that `results.h5` exists or that any row passes. Reference Verification owns numeric
U/UR outcome; Physics Evaluation owns force/moment balance, displacement direction,
symmetry, result signs, recovered-resultant consistency and physical plausibility.

### 6.4 Missing evidence classification

- blocking_for_current_formulation: `none`
- required_after_implementation: invariant, patch, MSVC build/CTest, declared
  reference comparison and physics evidence above
- nonblocking_optional: drilling coefficient sweep/energy ratio, `NR-O03`, `NR-O04`,
  canonical naming, README/metadata/provenance, expanded benchmark portfolio and
  broader mesh convergence studies
- future_only: nonlinear directional-derivative/objectivity/Newton evidence after
  its missing formulation decisions are separately approved

## 7. Required Revisions

### Formulation Agent

- None for the approved current linear-static implementation scope.
- Do not promote Section 15 to executable status until a separate formulation closes
  the nonlinear global coordinate map, objective drilling and load-work contracts.

### Research Agent

- None before current Implementation Planning.
- Optional locking/convergence characterization must remain clearly outside the
  approved implementation gate and must not imply MITC4+ or Abaqus equivalence.

### I/O Definition Agent

- None for the current numerical verdict. Preserve exact physical/full-residual
  distinction, source identity, locations, units and U-versus-UR decision rule.

### Reference Model Agent

- None for the current numerical verdict. Preserve the two declared S4 files and every
  optional existing reference artifact read-only, exclude S4R artifacts from acceptance
  comparison, and do not add administrative or portfolio gates.

## 8. Downstream Handoff

### Implementation Planning Agent

Implementation Planning is authorized and shall:

- trace the required tests in Section 6 to the approved requirement IDs before
  production work;
- keep `24 global -> 20 physical + 4 drilling` transforms, covariant MITC tying,
  common `2 x 2 x 2` integration, fixed drilling and physical recovery as explicit
  independent test seams;
- preserve stiffness assembly/partition/factorization-before-load, stable reduction,
  full-residual reaction and failure-atomic HDF5 lifecycle;
- keep future nonlinear execution, coefficient calibration, drilling output,
  `NR-O03/NR-O04`, reference-artifact mutation and Abaqus-equivalence claims outside
  the plan.

This handoff authorizes planning only. It does not authorize Harness execution,
production implementation, reference artifact changes, or completion claims.

### Reference Verification Agent

- Compare authoritative FESA HDF5 rows directly with the matching declared Abaqus
  displacement CSV after exact row-set precheck.
- Let only U1/U2/U3 affect pass/fail; report every UR1/UR2/UR3 warning without
  changing the verdict.

### Physics Evaluation Agent

- After reference verification, independently evaluate force and global moment
  balance, free residual, reaction sign, displacement direction, symmetry, positive
  physical energy and consistency of local resultants/stresses.

### Coordinator and Release Agents

- Record the Numerical Review gate as passed for planning at HEAD `a058ef7`.
- Do not infer implementation or release completion. Build/test, reference,
  physics-sanity and release-readiness gates remain pending.