docs(thermal): integrate reviewed formulation and verification contract

This commit is contained in:
KOKO\Mimi
2026-09-29 18:21:54 +09:00
7 changed files with 2063 additions and 610 deletions
+107 -148
View File
@@ -1,181 +1,140 @@
# Shell Transient Heat Transfer 결정 기록
- Feature ID: `shell-transient-heat-transfer`
- 문서 owner: 메인 assistant
- 갱신일: 2026-09-28
- Owner: 메인 assistant
- 상태: Accepted (구현 검증 완료와 구별)
- 승인 근거: ready 인터뷰의 개별 선택, 최종 “네 승인합니다”, 후속 `dryforge:go` 요청.
기준: [프로젝트 지침](../../AGENTS.md), [프로젝트 PRD](../PRD.md),
[공통 ADR](../ADR.md#adr-028). 지역 지침은 [Feature 지침](FEATURED_AGENTS.md),
요구사항과 현재 checklist는 [Feature PRD](PRD.md)에 둔다.
[공통 ADR](../ADR.md#adr-028). 요구사항과 현재 checklist는 [PRD](PRD.md),
지역 소유권은 [FEATURED_AGENTS](FEATURED_AGENTS.md)를 따른다.
아래 결정은 이전 Proposed 상태를 대체하며 실제 계산 결과의 정확도 승인을 의미하지 않는다.
수식 원본은 [formulation](formulation.md), 독립 판정은 [numerical-review](numerical-review.md)다.
각 결정의 owner는 메인 assistant다. 승인 근거는 2026-09-27~28 사용자 인터뷰의 명시 선택이다.
개별 범위 선택의 합의를 기록하며, 전체 requirement baseline과 구현 실행은 아직 승인되지 않았다.
사용자의 계획 기반 문서 작성 요청에 따라 수치 제안의 이유·대안·영향을 아래에 구체화한다.
`사용자 합의`는 명시 선택의 기록이고 `Proposed`는 검토할 우선안이다.
수식과 검산 원본은 [formulation](formulation.md), [numerical-review](numerical-review.md)를 따른다.
## STH-ADR-001: Uncoupled transient 범위
## 결정 목록
- 상태: Accepted
- 결정: 구조해석과 분리된 단일 Shell 과도 열전달 procedure를 구현한다.
- 이유: 사용자가 “지금은 uncoupled heat transfer를 구현합니다”라고 명시했다.
- 영향: 열응력·변형 피드백·양방향 연성은 장기 목표이며 이번 acceptance에서 제외한다.
| ID | 주제 | 상태 | 요구사항·검토 연결 |
| --- | --- | --- | --- |
| STH-ADR-001 | Uncoupled 범위 | 사용자 합의 | 001 |
| STH-ADR-002 | 5개 두께 온도점 | 사용자 합의 | 002~004 |
| STH-ADR-003 | 물성·초기/경계·열입력 | 사용자 합의 | 005~012 |
| STH-ADR-004 | 시간 증분·전체 결과 저장 | 사용자 합의 | 013~015 |
| STH-ADR-005 | 외부 비교·독립 검증 | 사용자 합의 | 016~017 |
| STH-ADR-006 | 단일 step | 사용자 합의 | 018 |
| STH-ADR-007 | 수치 방법론 전체와 K/C | Proposed | 정식화 3~5절 |
| STH-ADR-008 | Geometry와 접속 이상화 | Proposed | NR-01 |
| STH-ADR-009 | Backward Euler·종료 시간 | Proposed | 013, 정식화 7절 |
| STH-ADR-010 | 초기 jump와 pre-step 경계 | Proposed | 008~009, NR-02 |
| STH-ADR-011 | Heat flux·reaction 복구 | Proposed | 014~015, NR-05 |
| STH-ADR-012 | Edge 입력 단위·외부 문법 | Proposed / syntax 미결정 | 011, NR-03 |
## STH-ADR-002: 다섯 두께 온도점
요구사항 번호는 `FEM-REQ-SHELL-THERMAL-` prefix를 생략했다. 수치 기준이 미정인 NR-04는
[PRD 미결정 계약](PRD.md#구현-전-미결정-계약)에서 관리하며 승인된 tolerance를 대신하지 않는다.
- 상태: Accepted; 초기 두께 등온 선택을 대체
- 결정: 균질 centered section의 온도점은 항상 다섯 개다. 면내 및 두께 전도를 계산한다.
Source DOF 11/13/15는 bottom/midsurface/top이다.
- 이유: 기존 DS4 reference와 사용자가 선택한 두께 방향 범위를 유지한다.
- 영향: 적층·사용자 선택 온도점 수·삼각형은 제외하며 Abaqus 내부 정식화 동등성을 주장하지 않는다.
## STH-ADR-001: 현재 구현 목표를 uncoupled transient heat transfer로 한정
## STH-ADR-003: 상수 물성과 제한된 열입력
- 상태: 사용자 합의
- 결정: 현재 목표는 구조해석과 분리된 Shell 과도 열해석이다. 양방향 열–구조 연성은 장기 목표다.
- 근거: 사용자가 초기 연성 선택을 정정하며 “지금은 uncoupled heat transfer를 구현합니다”라고 명시했다.
- 영향: 이번 기능의 acceptance에 구조 변위·열응력·연성 반복이나 변형 피드백을 포함하지 않는다.
- 상태: Accepted; edge 부분은 STH-ADR-012로 대체
- 결정: 양의 유한 상수 k/rho/c와 요소 집합별 재료·일정 두께를 사용한다.
초기 온도는 node/Nset별로 지정하되 한 절점의 두께 전체에 동일하다.
지정 온도와 집중 열입력은 source DOF 11/15에만 허용하며 step 중 일정하다.
- 이유: 사용자 선택 범위에서 독립 검증 가능한 첫 열해석 기능을 구성한다.
- 영향: 온도 의존·이방성·표면 분포 flux·대류·복사·내부 발열·내부 station 직접 입력은 제외한다.
Concentrated P는 signed power이며 h나 면적을 곱하지 않는다. 중복 하중은 source 순서로 합한다.
## STH-ADR-002: 두께 방향 균일 온도에서 5개 온도점 단면으로 범위 변경
## STH-ADR-004: 모든 완료 증분과 초기 온도 저장
- 상태: 사용자 합의; 앞선 두께 방향 균일 온도 선택을 대체
- 결정: 3차원 공간의 4절점 Shell에서 면내·두께 방향 열전달을 지원한다. 균질 단면이고
두께 방향 온도점은 모든 단면에서 5개로 고정한다.
- 이유: 사용자가 기존 DS4 reference의 두께 방향 의미를 확인한 뒤 이를 지원하는 방향을 선택했다.
- 근거: reference 유지/범위 확대 질문에서 B, 균질 단면 질문에서 A, 온도점 5개 고정 질문에서 A.
- 영향: 아랫면·중앙면·윗면 온도의 identity를 구분한다. 적층·가변 온도점 수는 제외한다.
Abaqus의 내부 정식화를 복제하거나 동일한 결과 생성 절차를 약속하는 결정은 아니다.
- 상태: Accepted
- 결정: 초기 raw 온도와 모든 완료 증분의 온도·flux·반력 열률을 공식 HDF5에 저장한다.
- 이유: 한 시점의 reference 외에도 전체 과도응답을 검증해야 한다.
- 영향: Output no-op 요청이 필수 저장량을 줄이지 않는다. 초기 flux/reaction dataset은 없다.
종료 규칙은 009, 영속 schema는 013을 따른다.
## STH-ADR-003: 상수 물성·일정 입력과 제한된 두께 방향 적용 범위
## STH-ADR-005: NT11 외부 비교와 독립 물리 검증
- 상태: 사용자 합의
- 결정: 상수 등방성 열전도율·밀도·비열, 요소 집합별 재료와 일정 두께를 지원한다.
초기 온도는 절점/집합별이되 같은 절점의 두께 전체에는 균일하다. 지정 온도와 집중
열입력은 아랫면·윗면에 각각 적용하고, edge 열입력은 두께 전체에 균일하게 분포한다.
지정 온도·열입력은 시간에 무관한 일정값이다.
- 이유: 인터뷰에서 사용자가 각 첫 구현 범위를 선택했다.
- 영향: 면적 분포 표면 열유속, 내부 온도점 직접 입력, 대류·복사·내부 발열,
온도 의존·이방성 물성, 적층과 시간 이력 입력을 제외한다.
- 상태: Accepted
- 결정: 기존 CSV의 Step-1, t=37, NT11을 필수 비교하고 나머지는 해석해·수지·수렴으로 검증한다.
- 이유: Reference의 실제 coverage는 단일 시점의 아랫면 온도 729행이다.
- 영향: 공통 TOLERANCE의 row/family 기준을 유지한다. 추가 Abaqus CSV나 Abaqus 실행은 필수가 아니다.
원본 음수·단위·시간을 보정하지 않는다. Source identity 오류는 tolerance 전에 거부한다.
## STH-ADR-004: 고정 시간 증분과 전체 증분 결과 저장
## STH-ADR-006: 단일 step
- 상태: 사용자 합의
- 결정: 전체 시간과 고정 시간 증분을 사용자가 지정한다. 온도·요소 열유속·반력 열률을
모든 시간 증분에 저장하며 초기 온도 분포도 보존한다. 공식 출력은 공통 HDF5 계약을 따른다.
- 근거: 시간 제어 A, 필수 결과 C, 저장 간격 A 선택.
- 영향: 자동 시간 증분을 제외한다. 시간 적분법·열용량 행렬·마지막 증분 처리·초기 결과의
schema를 이 결정만으로 정하지 않는다.
- 상태: Accepted
- 결정: 초기 온도에서 일정한 지정 온도·열입력으로 하나의 transient step을 수행한다.
- 이유: 사용자의 명시 범위 선택이다.
- 영향: 여러 step, restart, 시간 변화 amplitude와 자동 시간 증분을 제외한다.
## STH-ADR-005: 현재 CSV의 NT11 비교와 독립 검증을 조합
## STH-ADR-007: Q4×P2 두 구간과 consistent K/C
- 상태: 사용자 합의
- 결정: 외부 필수 비교는 기존 `shell-thermal NodalTemperature.csv`의 `Step-1`, `t=37`,
`NT11`로 한정한다. 나머지 결과는 해석해·열에너지 수지·시간 및 공간 수렴성 등으로 검증한다.
- 근거: 추가 reference 결과 확보 여부 질문에서 사용자가 A를 선택했다.
- 이유: 현재 CSV는 한 시점의 아랫면 온도만 제공한다. 원본 evidence의 실제 coverage를 명시한다.
- 영향: 추가 Abaqus CSV 확보를 필수 조건으로 삼지 않지만, 내부 온도·열유속·반력 열률과
전체 과도응답의 correctness 검증은 유지한다. 공통 tolerance와 reference 원본은 변경하지 않는다.
- 상태: Accepted; 이전 수치 방법 Proposed 상태를 대체
- 결정: 요소당 20 DOF, Q4 면내 보간과 C0 두 구간 P2 두께 보간을 사용한다.
K=integral(B^T k B), C=integral(rho c H^T H), 면내 3×3 Gauss와 두께 exact Z를 채택한다.
- 이유: 전도·열저장·반력·에너지를 같은 Galerkin 보간으로 연결하며 직접 54점 적분으로 검산할 수 있다.
- 대안·영향: Lumped C, reduced integration, hourglass, 인공 확산과 온도 clamp는 추가하지 않는다.
Consistent C는 급변 경계에서 positivity를 보장하지 않으며 일반 warped patch의 3×3은 exact가 아니다.
Reference 실패에 맞춰 방법을 변경하지 않는다. 수치 상세 기준은 구현 전 독립 검토한다.
## STH-ADR-006: 단일 transient step으로 제한
## STH-ADR-008: 실제 midsurface metric과 normalized station 접속
- 상태: 사용자 합의
- 결정: 초기 온도와 일정한 경계조건·열입력으로 하나의 transient step을 수행한다.
- 근거: 단일/다중 step 지원 범위 질문에서 사용자가 A를 선택했다.
- 영향: 여러 해석 step과 step 사이의 온도 상태 전달·경계조건 변경을 제외한다.
- 상태: Accepted; 이전 물리 범위 검토 대기를 대체
- 결정: 실제 bilinear 3D midsurface의 metric/법선을 사용하고 두께 offset 곡률 보정을 생략한다.
공유 source node의 같은 normalized station을 접속하며 fold와 두께 단차도 같은 규칙을 사용한다.
- 이유: 승인한 4절점 Shell과 다섯 온도점을 일관된 product-domain으로 표현한다.
- 영향: 접합부의 상세 3D 전도·접촉 저항·추가 부피/용량을 해상하지 않는다.
방향이 일관된 manifold mesh만 허용하며 T-edge/vertex-only 접합을 거부한다.
전체 patch validity를 검증하고 수치적으로 증명할 수 없는 geometry는 진단한다.
임의 warp/aspect cutoff나 자동 normal flip·node weld로 지원 의미를 바꾸지 않는다.
## STH-ADR-007: 상세 수치 방법론 제안
## STH-ADR-009: Backward Euler와 마지막 짧은 증분
- 상태: **Proposed / 사용자 검토 대기**
- 제안: 면내 Q4 × 두께 2구간 P2의 20 DOF, midsurface metric의 thin-shell 근사,
consistent 열용량과 열전도행렬, 면내 3×3·두께 exact 적분, Backward Euler,
24개 trace 위치의 직접 Fourier flux 및 full residual 반력 열률을 사용한다.
- 근거: 사용자가 요구사항 목록을 보강하여 열전달행렬·heat flux 방법론을 포함한 상세 계획을 요청했다.
이 요청은 제안을 작성할 권한이며 위 수치 방법의 채택 승인은 아니다.
- 이유·대안: [설계 계획](design-plan.md)의 선택 표와 [정식화](formulation.md)의 상세 수식을 따른다.
- 영향: geometry·접속 이상화, 최초 온도 불일치, edge 입력 단위/syntax와 결과 위치를
검토해야 한다. 기존 reference 정확도와 positivity를 아직 보장하지 않는다.
- 검토: [독립 수치 검토](numerical-review.md)에서 명시 행렬·BE partition·에너지 수식의
일관성을 확인했다. 이는 사용자 채택 승인이나 solver 수치 검증을 대신하지 않는다.
- 행렬 선택 이유: `K=∫BᵀkB dV`, `C=∫ρcHᵀH dV`를 동일 보간에서 유도하여 전도·저장·
반력·수지를 연결한다. 두께 exact Z와 면내 3×3 적분은 독립 직접 적분으로 검산할 수 있다.
- 대안과 trade-off: 단일 4차 두께 보간, lumped C 또는 reduced integration도 별개 방법이지만
현재 우선안에 포함하지 않는다. Consistent C는 단조성·positivity를 보장하지 않으며,
일반 곡면에서 3×3은 exact가 아니다. 실패에 맞춰 방법을 조용히 바꾸지 않는다.
- 상세 선택: geometry, 시간 갱신, 초기 상태, recovery와 edge 입력은 008~012에서 분리한다.
- 상태: Accepted; 이전 정수 배수 종료 제한 제안을 대체
- 결정: A=C/delta+K를 사용하고 동일 계수 factorization을 재사용한다.
기본 dt로 진행하고 종료까지 남은 구간만 단축한다. dt>T_end이면 한 증분이다.
- 이유: 사용자의 시간 제어 선택을 유지하며 정확한 종료 시점을 저장한다.
- 영향: 마지막 delta가 실제로 달라지면 재분해하고 RHS·반력·수지 모두 같은 delta를 사용한다.
n*dt 기반 시간격자의 표현오차와 INC 사전검사를 고정한다. 무구속 transient도 유효하다.
Generic all-constrained 처리는 algebra 검증이며 공개 내부 station BC 지원이 아니다.
## STH-ADR-008: Thin-shell geometry와 두께 station 접속
## STH-ADR-010: Raw 초기값과 step BC override
- 상태: **Proposed / 물리 범위 검토 대기**
- Owner: 메인 assistant; 전문 계약은 formulation·numerical review 담당.
- 제안: 실제 3차원 Q4 midsurface의 metric·법선을 각 평가점에서 계산한다. 두께 offset의
곡률 보정을 생략한 product-domain을 사용하며 인접 요소는 같은 normalized 두께 station을 공유한다.
- 이유: 면내·두께 전도를 유지하면서 승인된 4절점 mesh와 5개 온도점을 연결할 수 있다.
- 대안·trade-off: 실제 offset-volume metric 또는 별도 3D 접합부 해상은 더 많은 geometry
계약을 요구한다. 우선안은 두께 단차의 노출면·spreading resistance·유한 접합부 부피를 해상하지 않는다.
- 영향: 곡률·왜곡·orientation·퇴화 허용 범위와 fold/이두께 접속의 물리 이상화를 구현 전에
검토해야 한다. 적분 차수를 높여도 이 모델링 오차는 제거되지 않는다.
- 근거: [정식화 2절](formulation.md#2-geometry-방향과-접속의-의미), 독립 검토 NR-01.
계획 문서화는 요청되었지만 이 이상화의 사용자 채택 승인은 아직 없다.
- 상태: Accepted; 이전 초기 상태·default 검토 대기를 대체
- 결정: t=0은 raw 초기 온도이며 미지정 값은 0이다. 첫 양의 증분에서 g를 적용한다.
첫 RHS의 C_fc*T_c_old/delta를 보존한다. Model BC를 상속하고 같은 target의 step BC가 덮어쓴다.
- 이유: 원래 초기 상태와 첫 증분의 열저장·반력·에너지 변화를 보존한다.
- 영향: 각 단계 안의 동일 BC 중복은 허용하고 상충은 거부한다. Thermal override를 구조 절차에 전파하지 않는다.
생략 BC 값은 0이며 초기값을 g로 projection하지 않는다. 순간 impulse를 별도 해상하지 않는다.
## STH-ADR-009: Backward Euler와 고정 증분 factorization 재사용
## STH-ADR-011: 직접 flux와 full residual 반력
- 상태: **Proposed / 시간 적분·종료 규칙 검토 대기**
- Owner: 메인 assistant; 전문 계약은 formulation·I/O 담당.
- 제안: `A=C/Δt+K`를 구성하고 지정 온도를 elimination으로 적용한다. 일정 계수·geometry·
구속 집합·Δt에서 factorization을 재사용한다. 전체 시간은 Δt의 정수 배로 제한하는 안을 제시한다.
- 이유: 상수 물성·일정 입력의 단일 step을 선형 solve로 구성하고 시간 1차 수렴을 검증할 수 있다.
- 대안·trade-off: 다른 시간 적분법은 별도 시간 오차·안정성 검토가 필요하다. 짧은 마지막
증분을 허용하면 factorization과 fixed-Δt 의미도 달라진다. 안정성만으로 정확도를 보장하지 않는다.
- 영향: 초기 constrained 온도 항을 포함한 RHS와 종료시간의 수치 판정을 확정한다.
무구속 transient 모델을 정적 K의 null mode만으로 singular 판정하지 않는다.
- 근거: [정식화 7절](formulation.md#7-backward-euler-지정-온도와-최초-증분), 독립 두 증분 검산.
사용자는 고정 Δt를 선택했으며 BE와 정수 배 제한은 아직 수치 제안이다.
- 상태: Accepted; 이전 recovery 검토 대기를 대체
- 결정: q_global=-k B T를 면내 네 점×두께 여섯 trace에서 직접 복구한다.
Global q와 local frame을 저장하며 local q를 중복 필수 저장하지 않는다.
반력 열률은 C*DeltaT/delta+K*T-Q의 constrained 성분이고 양수는 모델로 공급되는 열이다.
- 이유: 온도 기울기와 열저장을 포함한 반력을 서로 다른 observable로 검증한다.
- 영향: 0-/0+와 요소 경계 값을 평균하지 않는다. Free residual·열에너지는 내부 검증 evidence다.
Fourier flux 경계 적분으로 FE 반력을 대체하지 않는다. 초기 반력을 0이나 NaN으로 저장하지 않는다.
## STH-ADR-010: 초기 온도 보존과 첫 증분 경계 jump
## STH-ADR-012: Edge 분포 열입력 연기
- 상태: **Proposed / 초기 상태·default 계약 검토 대기**
- Owner: 메인 assistant; 전문 계약은 I/O·formulation 담당.
- 제안: 입력 초기 온도를 raw `T^0`로 보존하고 첫 양의 시점부터 step 지정 온도를 적용한다.
첫 RHS에 `C_fc(T_c^0−g)/Δt`를 유지하며 `t=0` 반력은 미정의로 구별한다.
- 이유: 입력 상태를 조용히 덮어쓰지 않고 첫 증분의 열저장·반력·에너지 변화를 추적할 수 있다.
- 대안·trade-off: 즉시 projection은 초기 상태와 에너지를 바꾸므로 projection jump의 저장·보고
정책이 필요하다. 두 정책은 같은 수치 결과를 보장하지 않는다.
- 영향: Reference의 step 전 `Set-10,11,11`과 step 내부 `Set-10,11,11,1000.`을 model/history
데이터로 구분한다. 생략 초기 온도·boundary 값과 재정의 의미를 확인한 뒤 한 정책을 채택한다.
- 근거: [정식화 7절](formulation.md#7-backward-euler-지정-온도와-최초-증분), 독립 검토 NR-02.
초기 온도의 면내 지정·두께 균일성은 합의됐지만 jump 처리의 채택 승인은 아직 없다.
- 상태: Deferred; 이전 edge kernel·syntax 제안을 대체
- 결정: FEM-REQ-SHELL-THERMAL-011을 이번 구현에서 제외한다. Edge kernel·keyword·공개 테스트를 추가하지 않는다.
- 이유: Ready에서 사용자가 집중 열입력만 먼저 구현하는 범위를 승인했다.
- 영향: 이전 q_in*h*L와 edge weighting 유도는 이번 제품 계약이 아니다.
후속 지원에는 별도 단위·부호·대상·중복 의미 및 수치 검증이 필요하다.
## STH-ADR-011: 위치별 Fourier flux와 residual 반력 열률
## STH-ADR-013: 입력 subset과 영속 thermal schema
- 상태: **Proposed / recovery·영속 결과 계약 검토 대기**
- Owner: 메인 assistant; 전문 계약은 formulation·I/O·numerical review 담당.
- 제안: `q=−kBT_e`를 면내 4점 × 두께 6개 trace의 24개 위치에서 계산한다.
중앙면의 `0−/0+`를 구별한다. 반력 열률은 `CΔT/Δt+KT−Q`의 구속 성분이며 모델 안으로의 공급을 양으로 둔다.
- 이유: Gradient flux와 열저장을 포함한 반력은 서로 다른 물리량이다. 각 위치의 raw 값을
유지해야 부호·열에너지 수지·경계 근처 공간 오차를 별도로 검증할 수 있다.
- 대안·trade-off: 절점 외삽·요소 간 평균은 원자료를 바꾸고 불연속을 숨길 수 있어 필수 결과에
포함하지 않는다. Global/local 중 저장 좌표계와 frame·trace identity는 아직 I/O 선택이다.
- 영향: 직접 적분 54점과 recovery 24점을 혼동하지 않는다. `t=0` 초기 온도를 보존하되
반력을 계산된 0으로 표시하지 않는다. Dataset layout·version·validity를 구현 전에 확정한다.
- 근거: [정식화 8~9절](formulation.md#8-반력-열률-에너지와-안정성), 독립 검토 NR-05.
온도·flux·reaction 저장은 합의됐으며 위 세부 출력 계약은 Proposed다.
- 상태: Accepted
- 결정: [PRD STH-C02](PRD.md#sth-c02--공개-입력-충돌과-부가-keyword)의 입력과
[STH-C05](PRD.md#sth-c05--hdf5-외부-계약)의 schema version 1을 구현한다.
- 이유: Domain/source identity에서 결과 row까지 재현 가능한 외부 계약을 고정한다.
- 영향: Temperature N×5, global flux E×4×6×3, constrained reaction P를 저장한다.
P=0도 빈 dataset을 유지하며 t=0 flux/reaction은 dataset 자체가 없다.
지정된 no-op만 warning과 함께 허용하고 unknown keyword는 오류다.
초기/완료 frame, mapping·units·finite 값을 atomic finalize 전에 self-check한다.
## STH-ADR-012: 두께 균일 edge 열입력과 공개 문법
## STH-ADR-014: Candidate 보존과 검증 evidence
- 상태: **Proposed / 단위·syntax 미결정**
- Owner: 메인 assistant; 전문 계약은 I/O 담당.
- 제안: 계산 kernel은 edge 옆면의 면적당 유입 열유속을 적분해 총률 `q_in hL`을 보존한다.
두 endpoint의 각 두께 가중치는 `(1,4,2,4,1)/24`다.
- 이유: 같은 shape function을 사용한 일관 하중으로 입력·열저장·반력을 연결한다.
- 대안·trade-off: 외부 입력이 길이당 열률이면 `q_in=λ/h`로 해석해야 하며 두께를 중복 곱하지 않는다.
현재 공개 근거만으로 DS4의 `*DFLUX, S1…S4`를 지원한다고 선언할 수 없다.
- 영향: 승인 가능한 keyword·대상 edge·단위·중복 application 의미를 확정해야 한다.
문법의 미결정을 이유로 합의된 edge 기능을 삭제하거나 표면 분포 열유속으로 대체하지 않는다.
- 근거: [계획안 5절](design-plan.md#5-열입력과-지정-온도),
[정식화 6절](formulation.md#6-가장자리절점-열입력), 독립 검토 NR-03.
두께 균일 edge 입력은 사용자 합의이며 외부 표현의 채택 승인은 아직 없다.
- 상태: Accepted
- 결정: Domain 정의, non-owning AnalysisModel, DofManager numbering, thermal mutable state를 분리한다.
Solve/recovery/검증을 마친 candidate만 반영하고 전체 성공 후 최종 HDF5를 원자 교체한다.
- 이유: 실패가 이전 accepted state와 기존 최종 결과를 오염시키지 않도록 한다.
- 영향: Deterministic COO reduction과 기존 CLI/lease/exit 분류를 유지한다.
전체 baseline을 먼저 재현하고 완료 시 전체 검사를 다시 실행한다.
재현된 무관한 HHT temporal 실패만 별도 미해결로 둘 수 있으며 새 실패·악화를 허용하지 않는다.
이를 전체 CTest green으로 보고하지 않는다. 열해석 필수 검증 실패는 면제되지 않는다.
@@ -14,10 +14,10 @@
## 계약과 작업 경계
- 인터뷰에서 합의한 지원 범위는 PRD, 채택 상태·이유·영향은 ADR에 기록한다.
`Proposed` 수치 방법과 독립 검토 결과를 사용자 채택 승인이나 구현 완료로 바꾸지 않는다.
ready에서 승인한 채택 계약과 독립 검토·실제 구현 완료 evidence를 구별한다.
- 수식 원본은 formulation, 구현 순서·코드 연결은 design-plan, 독립 지적과 검산 evidence는
numerical-review가 소유한다. 긴 유도·검토 내용을 PRD나 이 파일에 복제하지 않는다.
- 현재 checklist는 PRD 한 곳에 둔다. 문서 작성 요청만으로 solver 구현·외부 실행·커밋을 수행하지 않는다.
- 현재 checklist는 PRD 한 곳에 둔다. 현재 구현은 사용자의 `dryforge:go` 요청으로 승인되었다. 외부 solver 실행·배포 권한은 포함하지 않는다.
- 변경이 보간·물성·geometry 근사·초기 상태·적분·시간 갱신 의미를 바꾸면 관련 ADR과
정식화를 먼저 갱신하고 필요한 독립 numerical review를 받는다. Specialist는 할당 문서만 수정한다.
@@ -27,12 +27,12 @@
현재 요구는 균질 단면의 온도점 5개와 면내·두께 방향 열전달이다.
- Source DOF 11/13/15는 각각 아랫면/중앙면/윗면이다. Source DOF, 내부 equation 번호,
온도점·적분점·flux recovery 위치를 구분한다.
- 제안의 54개 직접 적분 위치와 24개 flux trace 위치는 서로 다르다. 중앙면의 두께
- 검산용 54개 직접 적분 위치와 24개 flux trace 위치는 서로 다르다. 중앙면의 두께
미분 양쪽 값을 평균하거나 요소 경계를 넘어 flux를 합치지 않는다.
- 반력 열률은 열용량 항을 포함한 full residual로 정의하는 제안이다. Fourier flux와
- 반력 열률은 열용량 항을 포함한 full residual로 정의한다. Fourier flux와
반력 열률을 같은 observable로 대체하지 않는다.
- Thin-shell 곡률 근사, 다른 두께·fold의 접속, 초기 온도 jump, edge 입력과 결과 schema는
[PRD 미결정 계약](PRD.md#구현-전-미결정-계약)을 먼저 확인한다. 미지원 의미를 임의 default로 처리하지 않는다.
- Thin-shell 곡률 근사, 다른 두께·fold의 접속, 초기 온도 jump와 결과 schema는
[PRD 상세 계약](PRD.md#sth-c01--모델-geometry와-자유도)을 따른다. Edge 입력 011은 이번 구현에서 연기한다. 미지원 의미를 임의 default로 처리하지 않는다.
- Reference는 PRD가 지정한 기존 `.inp`와 CSV를 그대로 사용한다. 음수 reference는
roundoff로 단정하거나 clamp하지 않고 [TOLERANCE](../TOLERANCE.md)를 적용한다.
- All-constrained와 단면 전체 end 온도 지정은 검산용 algebra/assembly fixture다.
@@ -44,12 +44,12 @@
| --- | --- | --- |
| [FEATURED_AGENTS.md](FEATURED_AGENTS.md) | 메인 assistant | 지역 지침·소유권·문서 연결 |
| [PRD.md](PRD.md) | 메인 assistant | 요구사항·제외·acceptance·현재 checklist |
| [ADR.md](ADR.md) | 메인 assistant | 결정·제안 상태·근거와 미해결 선택 |
| [ADR.md](ADR.md) | 메인 assistant | 채택 결정·승인 근거·이유와 영향 |
| [design-plan.md](design-plan.md) | 메인 assistant | 방법론 요약·코드 연결·구현 순서 |
| [formulation.md](formulation.md) | formulation-agent | 수식·geometry·보간·적분·시간 갱신·recovery |
| [numerical-review.md](numerical-review.md) | numerical-review-agent | 독립 검산·수치 위험·검증 설계 |
후속 I/O·구현·physics 작업의 owner는 실제 할당 시 기록한다. 필요한 PRD 항목과 관련
I/O·구현은 design-plan의 작업별 owner에게 할당하고 physics 평가는 독립 검토자에게 할당한다. 필요한 PRD 항목과 관련
전문 문서의 절, 수정 파일, 검증 기준을 전달하고 결과·변경 경로·미해결 항목을 반환한다.
## 변경 확인
+404 -179
View File
@@ -1,219 +1,444 @@
# Shell Transient Heat Transfer PRD
# Shell Transient Heat Transfer 요구사항과 승인 계약
- Feature ID: `shell-transient-heat-transfer`
- Owner: 메인 assistant
- 상태: 상세 계획을 반영한 요구사항 초안. 인터뷰 합의와 수치 제안을 구분하며 baseline 최종 검토 전이다.
- 기록일: 2026-09-28
- 상태: 사용자 승인 계약; 구현·검증 진행 중
- 승인 근거: ready 인터뷰의 명시 선택, 최종 “네 승인합니다”, 후속 `dryforge:go` 실행 요청.
기준: [프로젝트 지침](../../AGENTS.md), [프로젝트 PRD](../PRD.md),
[공통 ADR](../ADR.md#adr-028), [ARCHITECTURE](../ARCHITECTURE.md),
[TOLERANCE](../TOLERANCE.md). 지역 지침은 [Feature 지침](FEATURED_AGENTS.md),
결정 근거는 [Feature ADR](ADR.md)을 따른다.
[CODINGSTYLE](../CODINGSTYLE.md), [TOLERANCE](../TOLERANCE.md).
지역 지침은 [FEATURED_AGENTS](FEATURED_AGENTS.md), 결정 근거는 [ADR](ADR.md),
수식은 [formulation](formulation.md), 구현 순서는 [design-plan](design-plan.md),
독립 수치 판정은 [numerical-review](numerical-review.md)를 따른다.
## 제품 목적과 문서 범위
개인 연구용 솔버가 일차 목적이며 외부 배포와 양방향 열–구조 연성은 장기 목표다.
현재 목표는 3차원 공간의 4절점 shell에 대한 단일 step uncoupled transient heat transfer다.
입력 의미, 열 방정식, 시간 적분, 결과 복구, HDF5 저장과 독립 검증을 함께 구현한다.
채택 승인은 정확도 통과 evidence가 아니다. 아래는 승인된 실행 명세를 영속 문서에 반영한 계약이다.
초기 제안 중 edge 입력 지원과 정수 배수 종료 제한은 아래 계약으로 대체한다.
3차원 공간에 배치된 4절점 Shell의 시간에 따른 온도장과 열전달을 계산한다.
양방향 열–구조 연성은 장기 목표이며, 이번 기능은 구조해석과 분리된 uncoupled transient
heat transfer다. 해석 사용자는 승인된 Abaqus 입력 subset으로 모델을 정의하고, 개발·검증
담당자는 온도·열유속·반력 열률과 열에너지 수지를 재현할 수 있어야 한다.
## 요구사항과 경계
이 PRD는 [상세 설계·구현 계획안](design-plan.md), [정식화 제안](formulation.md),
[독립 수치 검토](numerical-review.md)를 연결하는 요구사항 원본이다. 현재 작업은 이 계획을
바탕으로 `FEATURED_AGENTS.md`, `PRD.md`, `ADR.md`를 구성하는 문서 작업이며 구현 실행은 포함하지 않는다.
`FEM-REQ-SHELL-THERMAL-`를 아래 표의 ID prefix로 사용한다. 기존 ID를 재번호화하지 않는다.
초기의 두께 방향 균일 온도 선택은 reference 확인 후 사용자가 명시적으로 변경했다.
현재 범위는 균질 단면의 두께 방향 온도점 5개와 면내·두께 방향 열전달을 포함한다.
## 상태와 채택의 구분
| 구분 | 현재 상태 | 의미 |
| --- | --- | --- |
| 기능 범위 | 인터뷰 합의 | 아래 001~018의 지원·제외 범위를 사용자가 선택했다. |
| 수치 방법 | Proposed | 20 DOF 보간·K/C·적분·BE·flux recovery의 우선안을 작성했다. 문서 작성 요청을 방법별 채택 승인으로 확대하지 않는다. |
| 독립 수치 검토 | 문서 검토·정확 산술 검산 수행 | 명시 행렬과 두 증분의 일관성 evidence가 있으며 geometry·초기 상태·I/O 선택은 남아 있다. |
| 구현·제품 검증 | 미실행 | Solver code, build/CTest, FESA/reference 비교와 실제 결과의 physics 검증 완료를 주장하지 않는다. |
## 사용자와 합의한 요구사항
아래 acceptance는 확인해야 할 관찰 가능한 동작과 검증 의도다. 테스트가 실행되었거나
정식화·I/O 세부 계약이 승인되었다는 뜻은 아니다.
| ID | 요구사항 | Acceptance intent |
| --- | --- | --- |
| FEM-REQ-SHELL-THERMAL-001 | 구조 변위와 연성하지 않는 과도 열전달 해석을 제공한다. | 구조해석 없이 초기 온도에서 시간별 열해석 결과를 생성한다. |
| FEM-REQ-SHELL-THERMAL-002 | 3차원 공간의 4절점 사각형 Shell을 지원한다. | 평판·접힌 판·곡면 mesh의 열전달과 모델 회전에 대한 결과 변환을 검증한다. 허용 요소 왜곡은 후속 계약에서 정한다. |
| FEM-REQ-SHELL-THERMAL-003 | 모든 단면에 두께 방향 온도점 5개를 사용하고 면내·두께 방향 열전달을 계산한다. | 두께 방향으로 온도가 달라지는 독립 사례와 균일 온도 사례를 구분하여 확인한다. |
| FEM-REQ-SHELL-THERMAL-004 | 각 단면은 두께 전체에 동일한 재료를 사용하는 균질 단면이다. | 적층 단면은 지원 범위에 포함하지 않는다. |
| FEM-REQ-SHELL-THERMAL-005 | 열전도율·밀도·비열은 온도에 무관한 상수이며 열전도율은 등방성이다. | 상수 물성의 전도·열저장 응답을 검증하고 온도 의존·이방성 입력은 지원으로 해석하지 않는다. |
| FEM-REQ-SHELL-THERMAL-006 | 여러 재료를 정의하고 요소 집합별로 배정할 수 있다. | 재료 배정과 서로 다른 물성 영역 사이의 열전달을 확인한다. |
| FEM-REQ-SHELL-THERMAL-007 | 요소 집합별 두께를 지원하며 각 요소 내부의 두께는 일정하다. | 서로 다른 두께가 전도·열저장·가장자리 열입력에 반영되는지 검증한다. |
| FEM-REQ-SHELL-THERMAL-008 | 초기 온도는 절점 또는 절점 집합별로 지정할 수 있다. | 비균일 면내 초기 온도를 읽고 해당 절점의 5개 두께 방향 온도점에 같은 값을 적용한다. |
| FEM-REQ-SHELL-THERMAL-009 | 아랫면·윗면에 각각 독립적인 지정 온도를 적용한다. | 한쪽 면만 구속한 경우와 두 면에 서로 다른 온도를 구속한 경우를 확인한다. 내부 온도점 지정 온도는 제외한다. |
| FEM-REQ-SHELL-THERMAL-010 | 절점 집중 열입력은 아랫면·윗면 온도점에 각각 독립적으로 적용한다. | 선택한 면의 열입력이 조립·열에너지 수지에 반영되는지 확인한다. 내부 온도점 직접 입력은 제외한다. |
| FEM-REQ-SHELL-THERMAL-011 | Shell 가장자리에 두께 전체로 균일하게 분포하는 열입력을 적용한다. | 대상 edge, 두께 방향 분포와 총 열률의 일관성을 독립 적분값으로 확인한다. |
| FEM-REQ-SHELL-THERMAL-012 | 지정 온도와 열입력은 해석 중 일정하다. | 시간 이력 표·진폭 함수에 의한 변화 없이 초기 온도에서 과도응답을 계산한다. 시작 시점의 적용 의미는 후속 계약에서 명시한다. |
| FEM-REQ-SHELL-THERMAL-013 | 사용자가 전체 해석 시간과 고정 시간 증분을 지정한다. | 지정된 시간 격자와 종료 시점을 확인한다. 나누어떨어지지 않는 종료 시간의 처리는 후속 계약에서 명시한다. |
| FEM-REQ-SHELL-THERMAL-014 | 온도, 요소 열유속, 지정 온도 경계의 반력 열률을 필수 결과로 제공한다. | 온도점·열유속 위치와 성분·반력 열률의 물리적 의미를 확인한다. 두께 방향 범위 확대에 따른 열유속 성분은 후속 출력 계약에서 명시한다. |
| FEM-REQ-SHELL-THERMAL-015 | 모든 시간 증분의 결과를 HDF5로 저장하고 초기 온도 분포도 보존한다. | 입력의 lexical sibling인 `<stem>_results.h5`에서 초기 온도와 완료 증분별 결과·시간·source identity를 확인한다. |
| FEM-REQ-SHELL-THERMAL-016 | 기존 reference CSV의 `Step-1`, `t=37`, `NT11`을 외부 비교 대상으로 한다. | source node와 두께 위치·시점을 대응시킨 뒤 공통 tolerance 정책을 적용한다. |
| FEM-REQ-SHELL-THERMAL-017 | 외부 CSV에 없는 결과도 독립적으로 검증한다. | 해석해·열에너지 수지·시간 및 공간 수렴성으로 전체 온도 이력, 내부 온도, 열유속과 반력 열률을 검증한다. 추가 Abaqus CSV 확보를 필수 조건으로 삼지 않는다. |
| FEM-REQ-SHELL-THERMAL-018 | 하나의 transient step만 지원한다. | 초기 온도와 일정한 경계조건·열입력으로 지정 시간까지 계산하고, 여러 해석 step의 입력은 지원 범위 밖으로 진단한다. |
## 제외 범위
- 열–구조 순차 연계와 양방향 연성, 변형에 따른 열전달 변화
- 3절점 요소, 적층 단면, 이방성·온도 의존 열물성
- 대류·복사·내부 발열, Shell 윗면·아랫면에 직접 적용하는 분포 열유속
- 두께 방향 비균일 초기 온도, 내부 온도점의 지정 온도·집중 열입력
- 두께 방향 비균일 가장자리 열입력, 절점별 가변 두께
- 사용자 선택 두께 방향 온도점 수, 자동 시간 증분, 시간 변화 경계조건·열입력
- 여러 해석 step 및 step 사이의 온도 상태 전달·경계조건 변경
- Abaqus 내부 요소·적분·열용량 처리·결과 생성 절차의 동등성
면의 온도점에 적용하는 절점 집중 열입력은 지원하며, 면적에 걸쳐 적용하는 분포 표면
열유속은 제외한다. 두 입력을 혼동하지 않는다.
## Reference 정적 관찰
사용자가 지정한 논리적 경로는 `reference/shell-thermal`이다. 2026-09-28 확인 시
현재 `ThermalAnalysis` worktree에는 없고 동일 저장소의 `dev` checkout에서 다음 원본을
읽었다. 원본 복사·수정·solver 실행은 수행하지 않았다.
| Artifact | 확인한 내용 |
| ID | 구현할 동작 |
| --- | --- |
| `reference/shell-thermal/shell-thermal.inp` | `DS4` 676개, 절점 729개, instance `Part-1-1`, 단일 재료, `*Shell Section` 데이터 `0.5, 5` |
| 같은 `.inp`의 재료 | Conductivity `54`, Density `7850`, Specific Heat `460`. 입력 숫자만 기록하며 단위계를 추정하지 않는다. |
| 같은 `.inp`의 step | `*Heat Transfer, end=PERIOD`와 `1., 100., , ,`; `Set-10, 11, 11, 1000.` 지정 온도. 명시적인 초기 온도·집중 열입력·가장자리 열입력은 없다. |
| 같은 `.inp`의 step 전 경계 | `*Boundary`의 `Set-10, 11, 11`이 존재한다. 생략된 값과 step 내부 재지정의 의미를 I/O에서 확인하며 model-level 경계를 무시하지 않는다. |
| `reference/shell-thermal/shell-thermal NodalTemperature.csv` | `Step`, `Frame`, `Node Label`, `NT11`; `Step-1`, `Increment 37: Step Time = 37.00`의 729행. 다른 시간·열유속·반력 열률은 없다. |
| 001 | 구조해석과 분리된 과도 열전달 procedure |
| 002 | 3차원 공간의 DS4 4절점 shell; 평판·유효한 warped patch·곡면 mesh·접힌 판 |
| 003 | 절점마다 두께 온도점 5개, 면내와 두께 방향 전도 |
| 004 | 단면 전체가 같은 재료인 centered homogeneous shell |
| 005 | 상수 등방성 양의 유한 열전도율 k, 밀도 rho, 비열 c |
| 006 | 요소 집합별 서로 다른 재료 |
| 007 | 요소 집합별 양의 유한 두께 h; 요소 내부에서는 일정 |
| 008 | 절점/절점 집합의 초기 온도; 동일 절점의 5점에는 같은 값 |
| 009 | 아랫면·윗면의 독립적인 지정 온도 |
| 010 | 아랫면·윗면의 독립적인 절점 집중 열입력 |
| 011 | **이번 구현에서 연기:** 모서리 분포 열입력. Kernel·keyword·공개 테스트 지원을 추가하지 않음 |
| 012 | 지정 온도와 열입력은 step 중 일정 |
| 013 | 사용자 지정 총 시간과 기본 시간 간격; 마지막 증분만 필요 시 단축 |
| 014 | 온도·요소 전역 열유속·지정 온도점의 반력 열률 |
| 015 | 초기 온도와 모든 완료 증분을 공식 HDF5에 저장 |
| 016 | 원본 reference Step-1, t=37, NT11의 필수 외부 비교 |
| 017 | 전체 온도·flux·반력에 대한 해석해·수지·시간 및 공간 수렴 검증 |
| 018 | 하나의 transient step만 지원 |
CSV에는 instance 열이 없으므로 `.inp`의 단일 instance로 row identity를 복원할 수 있는지
I/O 단계에서 명시적으로 검증한다. 파일명·행 순서만으로 다중 instance를 추정하지 않는다.
이 사례 하나로 다중 재료·두께, 가장자리·집중 열입력 또는 전체 시간 이력의 coverage를 주장하지 않는다.
독립 검토의 CSV 정적 통계에서 729개 label이 고유하고 값이 유한했으며 NT11 범위는
`[-15.5742216, 1000]`이었다. 초기 경계 급변·공간 이산과 관련된 undershoot 가능성을
조사해야 하며 reference 음수를 단순 roundoff로 취급하지 않는다.
추가 상세 계약은 이 문서의 `STH-C01`–`STH-C08`이다. 위 요구사항의 구체화이며 별도 기능을
추가하지 않는다. 제외: 열–구조 연성, 삼각형, 적층, 이방성·온도 의존 물성, 잠열,
대류·복사·내부 발열, 윗면/아랫면 분포 열유속, 모서리 분포 입력, 두께 비균일 초기 온도,
내부 온도점 직접 BC/CFLUX, 절점별 두께, 사용자 선택 온도점 수, 시간 이력·amplitude,
적응 시간 간격, 여러 step, restart, Abaqus 내부 정식화 동등성.
Abaqus의 공개 열전달 Shell 규약에서 자유도 11은 아랫면, 5개 온도점 단면의 13은 중앙면,
15는 윗면에 해당한다. 기하학적 절점 위치와 온도 자유도의 두께 위치는 구별한다.
[Abaqus 문서: Heat Transfer Shell Elements](https://docs.software.vt.edu/abaqusv2025/English/SIMACAEELMRefMap/simaelm-c-shellelem.htm).
원본은 CAE 2024 생성 표시가 있으며 공개 의미 확인에는 접근 가능한 2025 문서를 사용했다.
## STH-C01 — 모델, geometry와 자유도
## 제안하는 계산 방법과 상세 계획
### Source 의미
아래는 사용자가 선택한 기능 범위를 구현할 수 있도록 작성한 수치 방법 제안이다.
인터뷰로 확정한 요구사항과 구분하며 최종 채택 승인을 받은 것으로 기록하지 않는다.
- Shell geometry 절점은 midsurface에 있다. Source DOF `[11,12,13,14,15]`는 각각
`zeta=[-1,-1/2,0,1/2,1]`이며 11은 bottom, 13은 midsurface, 15는 top이다.
- 요소 절점 순서는 자연좌표 `(-1,-1),(1,-1),(1,1),(-1,1)`이다. 요소 로컬 DOF는
node-major `p=5*(i-1)+a`, i=1..4, a=1..5의 20개다.
- Source instance/label, Domain index, full DOF, equation index는 서로 다른 identity다.
집합을 전개해도 source row와 target identity를 진단까지 보존한다.
- 모든 열 방정식은 활성 DS4에 연결된 온도 DOF에 대해 구성한다. 비활성/미연결 절점을
열 방정식에 넣어 0 용량 행을 만들지 않는다. 그런 절점을 명시적으로 대상으로 하는 열
initial/BC/load는 유효한 열 target이 아니므로 모델 오류다. 일반 source model의 미사용
절점은 보존할 수 있으며 결과의 활성 절점 mapping으로 구분한다.
- 한 열해석 모델에 structural B33/S4/S4R 또는 coupled 요소를 혼합하지 않는다.
각 DS4에는 정확히 하나의 shell section과 thermal material이 연결되어야 한다.
Section 중복 지정은 값이 같아도 오류이며, set 내부의 동일 member는 한 번만 전개한다.
- 면내 Q4와 두께 두 구간의 C0 2차 보간으로 요소당 온도 자유도 20개를 구성한다.
- 실제 3차원 midsurface의 metric과 법선을 이용해 `∇T=B T_e`를 구성한다.
Thin-shell product-domain 근사와 접힌 판·다른 두께의 normalized-station 접속 의미를 명시한다.
- `K_e=∫BᵀkB dV`, `C_e=∫ρcHᵀH dV`를 사용한다. 두께의 명시 Z 행렬과
면내 3×3 Gauss 적분을 조합하고 직접 54점 적분으로 검산한다.
- Backward Euler의 `(C/Δt+K)T^(n+1)=Q+C T^n/Δt`를 free/constrained 공간으로
분할한다. 첫 증분의 이전 constrained 온도 항까지 보존하고 고정 계수 factorization을 재사용한다.
- `q=−kBT_e`로 면내·두께 방향 열유속을 직접 계산한다. 면내 4점과 두께 6개 trace의
총 24개 위치를 제안하며 중앙면의 양쪽 flux를 구별한다.
- 반력 열률은 `C(T^(n+1)−T^n)/Δt+K T^(n+1)−Q`의 구속 성분에서 계산한다.
### Geometry와 접속
[계획안](design-plan.md)은 각 선택의 이유, 단계별 수정 파일과 검증을,
[정식화](formulation.md)는 shape function·미분·행렬·적분·에너지 수식의 원본을 소유한다.
기존 `Element`/`DofManager`의 구조 자유도 전제를 확인했으며 별도 thermal kernel/state와
DOF layout 확장이 필요하다. 구체적인 수치 제안·I/O 의미·geometry 허용 범위를 검토한 뒤 구현한다.
`r=sum(N_i*x_i)`, `a1=r_,xi`, `a2=r_,eta`, `g=A^T A`, `A=[a1,a2]`,
`Js=|a1 cross a2|`, `n=(a1 cross a2)/Js`를 평가점별로 사용한다.
Warped 요소를 중심 평면에 투영하지 않는다.
## 입력·결과 계약의 범위
`grad_s(N_i)=A*g^-1*[N_i,xi,N_i,eta]^T`,
`grad_*T=grad_s(T)+n*(2/h)*T_,zeta`, `dV*=Js*(h/2)*dxi*deta*dzeta`다.
실제 offset surface의 곡률에 따른 metric·면적 변화는 포함하지 않는다.
이는 두께 방향 등온 가정이 아니다. 얇은 shell product-domain 모델의 한계를 문서화하며,
큰 h/곡률 또는 초기 두께 경계층에 대한 정량 정확도를 검증 없이 보장하지 않는다.
입력의 topology·material·section·load·boundary·step은 Domain에 보존하고, 실행 view와
mutable thermal state를 분리한다. `DS4`, `*SHELL SECTION`, 열물성과 `*HEAT TRANSFER`를
읽는 것만으로 기능 완료를 선언하지 않는다. Source label/instance는 방정식과 결과까지 보존한다.
공유 source node의 같은 normalized zeta 온도는 공유한다. 서로 다른 두께에서도 같으며,
실제 물리 z 좌표 일치, 접촉 열저항, 접합부의 유한 부피·추가 용량을 의미하지 않는다.
접힌 판·두께 단차의 국소 3D 전도장은 이번 모델의 해상 대상이 아니다.
| 물리량·입력 | 의미와 위치 | 단위 차원·확정 상태 |
| --- | --- | --- |
| 온도 | 절점마다 두께 위치 5개, initial 및 모든 완료 증분 | Temperature. 값의 단위계를 추정하지 않는다. |
| 열전도율·밀도·비열 | 요소 집합이 참조하는 상수 재료 | Power/(length·temperature), mass/volume, energy/(mass·temperature) |
| 절점 집중 열입력 | 아랫면 또는 윗면 온도점의 입력 | Power. 면적이나 두께를 중복 곱하지 않는 제안 |
| 가장자리 열입력 | 두께 전체에 균일한 edge 입력 | 외부 단위·syntax 미확정. 계산안의 면적당 flux는 power/area, 적분 후 총 열률은 power |
| Heat flux | 요소별·면내 위치·두께 trace·component·시점 | Power/area. 제안은 24개 trace 위치, global/local 저장 선택은 미확정 |
| 반력 열률 | 지정 온도 자유도·시점 | Power. 제안은 모델 안으로 공급되는 열을 양으로 정의 |
| 열에너지·free residual | 검증 evidence | Energy / power. 별도 필수 공개 dataset으로 추가 승인된 것은 아님 |
Topology는 방향이 일관된 manifold shell 표면이다. 하나의 edge에 최대 두 요소가 접속하고
공유 edge를 따라 절점 순서는 반대여야 한다. 불일치는 오류이며 자동으로 뒤집지 않는다.
T자 분기 edge, vertex-only 접합과 nonmanifold vertex fan은 거부한다.
분리된 여러 shell 영역은 허용하며 좌표만 같은 서로 다른 source node를 자동 병합하지 않는다.
서로 다른 instance도 명시적인 공유 source identity가 없으면 연결되지 않는다.
공식 출력은 입력 경로의 lexical sibling인 `<stem>_results.h5`다. Dataset layout·schema
version·초기 frame validity·좌표계는 I/O 계약에서 정한다. `t=0` 초기 온도 저장을
시간 차분이 없는 반력 열률의 계산 완료로 표시하지 않는다.
유한 좌표, 서로 다른 4절점, 양의 Js와 비특이 metric, 자기 교차 없는 patch를 요구한다.
적분점 몇 곳의 검사만으로 전체 patch 유효성을 증명했다고 하지 않는다.
Roundoff를 고려한 전체 patch의 유효성 판정과 geometry 오차 분리는 STH-C07의 선행
수치 계약 작업에서 고정한다. 양의 작은 Jacobian을 clamp하거나 구조 shell의 drilling/
director/두께 비율 제한을 검토 없이 전용하지 않는다.
## 검증 기준과 후속 계약 경계
## STH-C02 — 공개 입력, 충돌과 부가 keyword
외부 `NT11` 비교는 [TOLERANCE](../TOLERANCE.md)의 공통 정책을 적용한다. Row identity,
필수 component·시간·두께 위치의 대응과 missing·extra·duplicate·nonfinite 사전검사를
tolerance 전에 수행한다. 음수 reference를 zero-clamp하거나 결과에 맞춰 보정하지 않는다.
공통 ASCII·대소문자·source ordering·single physical file·part/identity-instance·set 규칙은
[공통 입력 계약](../contracts.md)의 기존 계약을 따른다. 열 입력의 추가 subset은 다음과 같다.
독립 검증은 최소한 균일 온도 보존, 열입력·열저장·반력 열률의 에너지 수지, 면내 및
두께 방향의 해석 가능한 전도 사례, 열유속·반력의 부호, 시간 및 공간 수렴성을 다룬다.
구체적인 사례·mesh·시간 간격·판정 기준은 정식화와 독립 numerical review에서 확정한다.
외부 비교 tolerance를 내부 에너지 수지나 수렴성 기준으로 대신 사용하지 않는다.
| 입력 | 허용 형식과 의미 |
| --- | --- |
| `*ELEMENT, TYPE=DS4` | 4절점; 기존 선택적 ELSET 및 topology 구문 사용 |
| `*SHELL SECTION, ELSET=..., MATERIAL=...` | 단일 row `h` 또는 `h,5`; 온도점 수 생략은 5. 명시 값은 5만 허용; offset/composite/nodal thickness 제외 |
| `*CONDUCTIVITY` | parameter 없는 상수 k 한 row; 선택적 마지막 빈 field 하나 허용 |
| `*DENSITY` | 기존 상수 rho 구문; 양의 유한 값 |
| `*SPECIFIC HEAT` | parameter 없는 상수 c 한 row; 선택적 마지막 빈 field 하나 허용 |
| `*INITIAL CONDITIONS, TYPE=TEMPERATURE` | model-level, `node-or-nset,T`; 각 target의 5점에 같은 유한 T |
| `*BOUNDARY` | model 또는 step, `target,dof,dof[,T]`; dof는 11 또는 15, T 생략은 0 |
| `*CFLUX` | step-level, parameter 없이 `target,dof,P`; dof 11 또는 15, 유한한 signed P |
| `*STEP` | 기존 NAME, INC; NLGEOM은 생략 또는 NO, PERTURBATION 제외 |
| `*HEAT TRANSFER` | END 생략 또는 PERIOD; 한 row `dt,T_end` 뒤 최대 3개의 빈 field 허용 |
| `*END STEP` | 단일 step 종료 |
정식화·계획안에서 제시한 방안은 아직 사용자 미승인이다. 다음 계약은 제안을 검토하여
확정할 필요가 있으며 기존 인터뷰 범위를 조용히 변경하지 않는다.
추가 material temperature/field columns, ANISOTROPIC/ORTHOTROPIC, STEADY STATE,
DELTMX, 시간증분 최소/최대의 명시 숫자, OP/AMPLITUDE 등 미지원 option은 거부한다.
`*BOUNDARY`의 11..15 범위 표기는 내부 station을 포함하므로 거부한다.
물성·두께·시간은 양의 유한 값이며 온도는 signed finite scalar다. 섭씨/절대온도 단위를
숫자로 추정하지 않는다. Reference 파일의 상수 값을 변환하지 않는다.
- 정식화: 곡면·왜곡 geometry의 허용 범위, 서로 다른 두께 단면의 접속 의미, 온도 보간,
면내·두께 적분, 열용량 행렬, 시간 적분법, 초기/지정 온도 적용 순서와 반력 열률 정의
- I/O: 승인 keyword subset, edge 열입력의 단위·부호·대상 표기, 생략된 초기 온도 등 default,
중복/상충 집합 진단, source identity, 열유속 위치·좌표계·성분, 초기 상태와 증분의 HDF5 schema
- Numerical Review: reference frame 대응과 comparison family, analytical/invariant 사례와
각 판정 기준, 고정 5개 온도점 조건에서의 공간 수렴성·두께 해상도 한계
- 구현 전 runtime 계약: 별도 thermal procedure/state, 실패 candidate와 기존 최종 HDF5 보존,
deterministic assembly·결과 순서, 기존 구조해석의 회귀 검증 범위
미지정 초기 온도는 0이다. 같은 초기 선언에서 집합이 겹치면 동일 값은 허용하고 다른
값은 오류다. BC는 model 단계와 step 단계 각각에서 target을 완전히 전개하여 충돌을
확인한 뒤 step의 값이 동일 node/DOF에 대한 model 값을 덮어쓴다. 각 단계 안의 동일
값 중복은 허용하며 다른 값은 오류다. 이 override는 thermal에만 적용한다.
## 기능 완료 기준
집중 입력은 source row 순서로 더한다. 양수 P는 모델 안으로 공급되는 열률이며 음수는
열 제거다. P에 h·면적·shape weight를 곱하지 않는다. 지정 온도와 겹치는 입력도 full Q에
유지하여 반력에 반영한다. 온도 또는 열입력이 지정되지 않은 외부 경계는 자연 단열이다.
| 검증 묶음 | 요구사항 연결 | 필요한 evidence |
| --- | --- | --- |
| 입력과 모델 의미 | 001~013, 018 | 양성 입력의 semantic mapping과 미지원·상충·불량 값 진단; pre-step boundary와 단일 transient step 구분 |
| Element·조립 | 002~007, 010~011 | Shape·K/C·edge 일관 하중 검산, source DOF scatter, 결정적인 조립·결과 순서 |
| 시간 적분과 물리 | 008~013, 017~018 | 초기 경계 불일치, 해석해·에너지 수지·시간/공간 수렴, 실패 candidate 보존 |
| 결과와 저장 | 014~015 | 모든 온도점·완료 증분, flux trace, reaction의 identity·단위·좌표계 및 HDF5 round trip·실패 보존 |
| 외부 비교 | 016 | 원본 NT11 t37 row 대응 사전검사와 공통 tolerance의 blocking 비교 |
| 통합 회귀 | 전체 | 구성된 MSVC x64 Debug build/CTest, 신규 warning 없음, 영향받는 구조해석 회귀와 독립 검토 결과 |
Reference에 있는 다음 부가 구문만 열 procedure의 no-op으로 허용한다. Data row가
필요 없는 구문에 data를 추가하거나 미지원 option을 붙이면 오류다. 허용한 각 no-op
block에 source 위치를 갖는 warning을 하나씩 남기며 필수 결과를 줄이지 않는다.
단일 kernel 통과, reference 한 시점의 일치 또는 문서의 정확 산술 검산만으로 위 묶음을
대체하지 않는다. 검증 상세는 [계획안 9~11절](design-plan.md#9-구현-순서와-각-단계의-완료-조건)과
[독립 검토 3절](numerical-review.md#3-구현-후-반증-가능한-검증-계획)을 따른다.
- `*PREPRINT, ECHO=NO, MODEL=NO, HISTORY=NO, CONTACT=NO`.
- `*RESTART, WRITE, FREQUENCY=0`: restart 파일을 생성하지 않는다.
- `*OUTPUT, FIELD, VARIABLE=PRESELECT` 및 `*OUTPUT, HISTORY, FREQUENCY=0`.
Child NODE/ELEMENT OUTPUT 및 그 밖의 output 설정은 제외한다.
- 하중이 참조하지 않는 assembly `*SURFACE, TYPE=ELEMENT, NAME=...`:
기존 target 및 SPOS/SNEG facet 의미를 검증하고 보존하되 열하중으로 사용하지 않는다.
- Model-level `*PHYSICAL CONSTANTS, ABSOLUTE ZERO=<finite>`:
이번 상수 물성·비복사 모델에서 온도 변환이나 에너지 영점 이동에 쓰지 않는다.
## 구현 전 미결정 계약
Unknown/unsupported keyword를 일반적으로 무시하는 fallback은 없다.
| 항목 | 우선 제안·남은 선택 | 결정·상세 계약 owner | 근거 |
## STH-C03 — 행렬과 시간 적분
두께는 `[-1,0]`, `[0,1]`의 C0 quadratic 구간 둘로 구성한다. 각 구간에서
`s=2*zeta+1` 또는 `s=2*zeta-1`로 두고
`l=[s(s-1)/2, 1-s*s, s(s+1)/2]`를 각각 station `(1,2,3)`과 `(3,4,5)`에 놓는다.
`dL/dzeta=2*dl/ds`이며 중립면의 온도는 하나지만 미분의 두 trace는 별도다.
`H_ia=N_i*L_a`, `B_ia=L_a*grad_s(N_i)+n*N_i*(2/h)*L_a,zeta`다.
```text
M_ij = integral_surface N_i N_j dA
S_ij = integral_surface grad_s(N_i) dot grad_s(N_j) dA
Z0_ab = integral_-1^1 L_a L_b dzeta
Z1_ab = integral_-1^1 L_a,zeta L_b,zeta dzeta
C_e = rho*c*(h/2) * (M tensor Z0)
K_e = k * [(h/2)*(S tensor Z0) + (2/h)*(M tensor Z1)]
Z0_segment = [[4,2,-1],[2,16,2],[-1,2,4]] / 30
Z1_segment = [[7,-8,1],[-8,16,-8],[1,-8,7]] / 3
```
각 segment의 3×3 행렬을 공유 station 3에 조립한다. 면내 production 적분은 3×3
Gauss, 두께는 위 exact Z다. 동일 보간의 직접 54점 적분과 교차 검산한다.
Affine 평판의 exact 면내 2×2 검산과 production 규칙을 혼동하지 않는다.
Structural 질량을 C로 재사용하지 않으며 lumping, reduced integration, hourglass,
인공 확산·온도 clamp를 추가하지 않는다. Consistent C는 급변 BC에서 undershoot를
만들 수 있으므로 온도 범위 초과만으로 성공/실패를 결정하지 않는다.
기본 dt>0, 총 시간 T_end>0에서 종료까지의 증분을 구성한다. 마지막 구간이 남으면 그
길이 delta로 단축한다. dt>T_end이면 delta=T_end인 한 증분이다. 적응 제어는 없다.
일반 시점은 누적 덧셈 대신 n*dt로 만들고 마지막 저장 시점은 입력 T_end다.
정수 배수의 부동소수점 표현 오차 때문에 거의 0인 가짜 증분을 만들지 않는다.
시간 비교의 roundoff bound는 입력·곱셈·나눗셈 연산을 근거로 STH-C07에서 고정한다.
정상 간격과 표현 오차 수준에서 같은 마지막 간격에는 같은 계수를 사용하며, 실제로
단축된 마지막 간격은 모든 방정식·반력·에너지에 동일한 delta를 사용한다.
시점이 유한하고 엄격히 증가해야 하며 표현할 수 없는 격자는 진단으로 거부한다.
INC는 기존 STEP 기본 100 또는 명시 양의 uint32다. 필요한 완료 증분 수가 INC를
넘으면 실행 전에 오류로 처리하며 중간 결과를 성공 결과로 저장하지 않는다.
각 증분의 Backward Euler 식은 다음과 같다.
```text
A = C/delta + K
A_ff*T_f_new = Q_f + C_ff*T_f_old/delta + C_fc*T_c_old/delta - A_fc*g_new
T_c_new = g_new
```
입력 초기 T_raw를 t=0에 보존하며 첫 old state로 그대로 쓴다. t=0에 g로 projection하지
않는다. 따라서 첫 증분의 `C_fc*(T_c_old-g_new)/delta`는 생략할 수 없다.
이 정책은 초기 불연속을 첫 유한 증분에 반영하며 순간 impulse를 해상한다는 뜻은 아니다.
동일 계수의 A_ff factorization은 재사용하고 마지막 delta가 바뀌면 다시 구성·분해한다.
무구속 단열 모델도 C가 양정치여서 transient A는 양정치다. 연결 component별 K의 상수
null mode를 과도계의 singular error로 오인하지 않는다. Generic 0×0 free partition은
solve 없이 처리하되 공개 입력은 내부 station을 구속하지 못하므로 all-constrained는
algebra fixture로만 검증한다.
## STH-C04 — 열유속, 반력과 수지
각 완료 온도장으로 `q_global=-k*B*T_e`를 직접 평가한다. 단위 차원은 power/area다.
`e1=a1/|a1|`, `e2=n cross e1`, `R=[e1,e2,n]`이고 후처리 local 성분은 `R^T*q_global`다.
Local q3는 bottom에서 top으로 흐를 때 양수다. Local 성분을 중복 필수 저장하지 않는다.
면내 점 순서는 `(-g,-g),(g,-g),(g,g),(-g,g)`, g=1/sqrt(3)이다.
각 점에서 두께 순서는 `(-1,+),(-1/2,interior),(0,-),(0,+),(1/2,interior),(1,-)`로
고정한다. 총 24개/요소다. 부호는 해당 두께 위치로 접근하는 쪽이며 바깥 법선 부호가 아니다.
중립면 양쪽, 재료 경계와 인접 요소의 값을 평균하지 않는다. Nodal smoothing/extrapolation을
공식 flux로 추가하지 않는다. Point source 부근의 mesh 독립적인 유한 flux를 보장하지 않는다.
```text
r = C*(T_new-T_old)/delta + K*T_new - Q
reaction = r의 constrained 성분 [power, 양수: 모델에 공급]
free residual = r의 free 성분 [내부 검증 evidence]
DeltaE = 1^T*C*(T_new-T_old) [energy]
DeltaE = delta*(sum(Q)+sum(reaction)+sum(free residual))
```
물리 수지 defect와 free residual 크기를 별도로 검사한다. 위 등식이 대수적으로 맞는다는
사실만으로 free equation 정확도를 주장하지 않는다. 에너지 영점을 자동 추정하지 않으며
온도 기준 이동에 대한 에너지 변화의 불변성을 검사한다. Fourier flux의 경계 적분은
일반적으로 FE 반력과 같지 않으며 반력 대신 사용하지 않는다. BE의 C-norm 소산도
실제 열에너지 보존과 별도 검증한다. t=0 반력은 미정의다.
## STH-C05 — HDF5 외부 계약
CLI는 기존 `fesa.exe <model.inp>`이고 output은 lexical sibling `<stem>_results.h5`다.
새 thermal 형식은 기존 structural 형식을 변경하지 않는다. 아래 모든 필드는 mandatory다.
F64/U64/U8/I8은 고정 폭 little-endian, 문자열은 UTF-8다. Compound는 field/type/extent로
읽으며 메모리 padding을 계약으로 삼지 않는다. Float 결과는 모두 finite다.
### Metadata와 model
`/metadata` scalar attributes:
| 이름 | 값/형식 |
| --- | --- |
| schema_version | U64, 1 |
| feature_id | `shell-transient-heat-transfer` |
| analysis_procedure | `uncoupled-transient-heat-transfer` |
| source_step_name | 입력의 원문 step 이름, 생략 시 Step-1 |
| solver_version | 기존 build version |
| source_input_identity | 실제 읽은 입력 byte의 기존 identity 계약 |
| unit_system_label | `user-consistent-unspecified` |
| internal_formulation | `FESA-Q4-P2x2-thermal-shell` |
| integration_rule | `surface-3x3-gauss; thickness-exact-P2x2` |
| time_integration | `backward-euler` |
| nominal_time_increment | F64, 입력 dt |
| time_period | F64, 입력 T_end |
Model의 기본 row 순서는 기존 Domain의 stable source declaration/instance expansion 순서다.
| Dataset | Shape/field와 의미 |
| --- | --- |
| `/model/nodes` | 기존 node compound: internal_node_id U64, instance_name/string, source_label/string, coordinates F64[3] |
| `/model/elements` | thermal DS4 row: internal_element_id U64, instance_name/string, source_label/string, source_element_type=`DS4`, internal_formulation/string, node_internal_ids U64[4], shell_section_internal_id U64, material_internal_id U64 |
| `/model/thermal_materials` | internal_material_id U64, name/string, conductivity F64, density F64, specific_heat F64 |
| `/model/shell_sections` | internal_section_id U64, source_file/string, source_line U64, source_elset/string, material_internal_id U64, thickness F64 |
| `/model/thermal_node_ids` | U64[N], 활성 온도 절점의 internal_node_id, 원래 node 순서를 유지 |
| `/model/thermal_element_ids` | U64[E], flux 요소 행에 대응하는 internal_element_id, 원래 element 순서 유지 |
| `/model/temperature_dofs` | U64[5]=[11,12,13,14,15] |
| `/model/temperature_zeta` | F64[5]=[-1,-0.5,0,0.5,1] |
| `/model/thermal_constraint_mask` | U8[N,5], free=0/constrained=1 |
| `/model/thermal_reaction_dofs` | compound[P]: internal_node_id U64, source_dof U64; mask의 1인 항목을 thermal node 행/DOF 순서로 나열 |
| `/model/thermal_flux_inplane_coordinates` | F64[4,2], 위 xi/eta 순서 |
| `/model/thermal_flux_zeta` | F64[6]=[-1,-0.5,0,0,0.5,1] |
| `/model/thermal_flux_trace_side` | I8[6]=[1,0,-1,1,0,-1]; 0=interior, +/-1=위의 접근 방향 |
| `/model/thermal_flux_local_frames` | F64[E,4,3,3], 마지막 두 축은 global component와 local axis; R의 열이 e1,e2,n |
ID는 해당 model table에 존재해야 하며 중복·누락 mapping은 self-check 실패다.
Reaction mapping은 source DOF 11/15만 포함한다. P=0이면 길이 0인 dataset을 유지한다.
Material/section의 source provenance는 기존 모델 writer convention을 보존한다.
Structural elastic E/nu나 displacement mask를 thermal material/DOF 대신 기록하지 않는다.
### Frames
Storage slot은 `/steps/Step-1`로 고정하며 실제 step 이름은 metadata로 구별한다.
`/steps/Step-1/accepted_frames`는 compound[N_step+1]이고 field는
`frame_index U64`, `increment_index U64`, `time F64`, `time_increment F64`, `is_initial U8`다.
Frame/index는 0..N_step이고 서로 같다. 초기 행은 time=0, time_increment=0, is_initial=1,
나머지는 완료 시간과 해당 방정식에 사용한 양의 delta, is_initial=0이다.
각 `/steps/Step-1/frames/<decimal-index>` group에도 위 다섯 scalar attribute를 저장한다.
| Frame 상대 경로 | Shape | 초기 frame | 완료 frame |
| --- | --- | --- | --- |
| 보간·행렬·적분 | Q4×두 P2 구간, consistent C, 면내 3×3·두께 exact Z의 채택 | 메인/사용자, formulation·numerical review | STH-ADR-007 |
| Geometry·접속 | Thin-shell 곡률 근사와 normalized station 접속의 적용 범위·퇴화 기준 | 메인/사용자, formulation·numerical review | STH-ADR-008, NR-01 |
| 시간 제어 | Backward Euler 및 종료시간을 Δt 정수 배로 제한할지 | 메인/사용자, formulation·I/O | STH-ADR-009 |
| 초기 상태 | Raw 초기값의 첫 증분 사용, pre-step 경계·생략값·step 재정의 의미 | 메인/사용자, I/O·formulation | STH-ADR-010, NR-02 |
| 결과 계약 | 24개 trace, global/local 저장, 반력 부호, t=0 validity | 메인/사용자, I/O·numerical review | STH-ADR-011 |
| Edge 입력 | 면적당 또는 길이당 단위·공개 keyword·대상 및 중복 의미 | 메인/사용자, I/O | STH-ADR-012, NR-03 |
| 내부 수치 판정 | Geometry margin, 적분 오차, residual·수렴성 기준 | numerical review, 해당 계약 owner | NR-04; 외부 tolerance와 구분 |
| `nodal/temperature` | F64[N,5] | raw 초기 온도 | 완료 온도 |
| `nodal/reaction_heat_rate` | F64[P] | **dataset 없음** | constrained full residual, reaction mapping 순서 |
| `element/heat_flux` | F64[E,4,6,3] | **dataset 없음** | global qx,qy,qz |
Owner 역할은 후속 할당 책임이며 아직 수행하지 않은 I/O 작업이나 새 승인을 뜻하지 않는다.
NR 항목의 상세 근거는 [독립 수치 검토](numerical-review.md#1-영향이-큰-판단과-필요한-결정)에 있다.
Temperature attributes는 `component_order=NT11,NT12,NT13,NT14,NT15`,
`units_label=temperature`, `coordinate_system=scalar`, `location=nodal-thickness-stations`다.
Reaction은 `component_order=RFL`, `units_label=power`, `coordinate_system=scalar`,
`location=prescribed-temperature-dofs`, `positive_direction=into-model`이다.
Flux는 `component_order=QX,QY,QZ`, `units_label=power/area`,
`coordinate_system=global-cartesian`, `location=element-inplane-thickness-traces`다.
Model geometry의 좌표·두께는 length, normalized coordinates/frames/masks/IDs는
dimensionless다. `/model/nodes`는 `units_label=length`,
`/model/shell_sections`는 `thickness_units_label=length`를 가진다.
`/model/thermal_materials`의 세 attribute는
`conductivity_units_label=power/(length*temperature)`,
`density_units_label=mass/volume`, `specific_heat_units_label=energy/(mass*temperature)`다.
`/steps/Step-1/accepted_frames`에는 `time_units_label=time`을 저장하며 time과
time_increment 모두에 적용한다. 서로 다른 단위의 compound field를 하나의
모호한 units_label로 합치지 않는다.
## 사용자 검토 항목
초기 frame의 누락 flux/reaction은 정해진 schema이며 NaN이나 계산된 0으로 대체하지 않는다.
Free residual·열에너지는 검증 보고서의 evidence로 유지하고 필수 공개 dataset을 추가하지 않는다.
Warnings는 기존 `/diagnostics` compound convention의 source/severity/code/message를 유지한다.
형식 discriminator, schema, mandatory field·shape·identity·finite value·frame 완전성을
atomic finalization 이전에 self-check한다. Unknown optional additive field를 무시할 수
있으나 잘못된 mandatory schema를 허용하지 않는다.
개별 범위 선택 인터뷰를 마쳤다. 이 초안의 범위·제외·acceptance intent를 사용자에게
요약하여 제시하고 문서 전체의 requirement baseline을 확정한다. 구현 실행은 별도 요청 범위다.
## STH-C06 — 런타임 소유권과 실패
후속 정식화·I/O에서 physics나 공개 의미의 추가 선택이 필요하면 근거와 함께 사용자에게
돌려보낸다. 위 전문 항목을 이미 결정되었거나 검증된 것으로 취급하지 않는다.
Domain은 thermal material/section/initial/BC/load/step 정의를 소유하며 parsing 이후
불변이다. AnalysisModel은 non-owning active view, DofManager는 열 5개와 구조 6개의
layout, numbering, scatter, constraint 및 sparse pattern을 전담한다. Equation ID를
Node/Element에 분산시키지 않는다. Thermal kernel/procedure/state/recovery를 구분하고
구조 displacement API에 temperature를 위장하거나 structural mass로 C를 취급하지 않는다.
미래 coupling state를 미리 할당하지 않는다. Vendor type은 기존 adapter 경계 안에 둔다.
Assembly는 element-local contribution, stable COO ordering/reduction을 사용한다.
병렬 worker가 global CSR를 직접 갱신하지 않는다. 같은 입력은 같은 sparse structure,
축약 순서, diagnostic order와 result row order를 만든다.
증분마다 solve, full temperature 복구, residual/flux/finite/수지 검증을 끝낸 candidate만
accepted state에 반영한다. 실패 시 마지막 accepted state를 오염시키지 않는다.
최종 output은 전체 해석 성공과 writer self-check·close 후 원자적으로 교체한다.
실패 시 기존 최종 HDF5를 보존하며 로그는 실패 진단으로 갱신될 수 있다.
Output lease와 CLI exit 2/3/4/5/6의 기존 분류를 유지한다. 초기 해석 성공이나 이전
파일 존재만으로 이번 실행 성공을 표시하지 않는다. Failed trial을 accepted_frames에 넣지 않는다.
## STH-C07 — 독립 검증과 수치 계약의 선행 조건
채택된 수식의 일관성에 관한 기존 독립 검토와 실제 solver 정확도 검증을 구별한다.
수치 전문 작업은 kernel 구현 전에 [정식화](formulation.md)와 [독립 수치 검토](numerical-review.md)에 아래 세부 계약을 문서화하고 독립 numerical review를
통과해야 한다. 이는 기존 PRD가 numerical owner에게 맡긴 상세 설계 작업이다.
1. 전체 bilinear patch 유효성·orientation 판정, 비다양체 판정, scale-aware degeneracy
bound와 검증할 geometry군. 전역 자기 교차 검출이나 contact 기능을 추가하지 않는다.
2. Roundoff metric: 행렬 norm·연산 수·machine epsilon·conditioning을 근거로 exact
fixture 및 free residual의 bound를 도출한다. 예: backward error의 분모는
`||A_ff||*||T_f||+||b_f||`; 분모 0은 별도 처리한다. `max(1,scale)` floor 금지.
3. 에너지 defect는 저장·외력·반력 항의 절댓값 합을 기준으로 정의하고 free residual
평가와 분리한다. 이산 수지와 BE 소산 각각의 roundoff 전파를 문서화한다.
4. 3×3/5×5/7×7 quadrature 비교, 시간/공간 convergence의 norm·mesh·시간·series
절단 오차·판정 구간을 독립 oracle와 예상 asymptotic order에서 도출한다.
5. 고정 5점의 slab 오차 plateau와 초기 열침투 깊이 sqrt(alpha*t)를 보고하는 기준,
비교 시점과 각 물리량의 허용 이산 오차를 테스트 실행 전에 고정한다.
6. 시간격자 roundoff 동등성·increment count·overflow와 representability 판정을 고정한다.
이 수치 기준을 구현 결과에 맞춰 완화하거나 외부 reference 5%/1%를 내부 검사에
전용하지 않는다. 검토에서 3×3 방법 변경, 지원 geometry 축소, 추가 물리 또는 accuracy
acceptance 변경이 필요하면 구체적 근거를 제시하고 사용자 결정으로 돌아간다.
수치 상세 기준의 독립 검토가 끝나기 전 production kernel 착수를 통과로 표시하지 않는다.
필수 검증은 다음과 같다. 별도 허용을 뜻하는 fixture와 공개 CLI 사례를 구별한다.
| 검증 | 독립 oracle와 실패 검출 |
| --- | --- |
| Shape/K/C | sum(H)=1, sum(B)=0, 정확 M/S/Z, C SPD, K rank19/상수 null, symmetry, 양의 에너지, h 및 1/h scaling, 분리식/직접54점 일치 |
| 방향/flux | constant/면내선형/두께선형 patch, global rotation covariance, normal 반전과 station 재매핑, R 직교성 |
| 중립면 trace | kernel state (0,0,0,1,2)에서 q3(0-)=0, q3(0+)=-4k/h; 공개 초기조건 사례로 오인하지 않음 |
| 초기 jump | C=[[2,1],[1,2]], K=[[1,-1],[-1,1]], delta=1/2, old=(0,0), 두 번째 g=1: 첫 free=-1/5, 다음=1/25, 반력24/5와36/25 |
| 시간격자 | 정수배, 10/3 마지막1, dt>T, decimal roundoff, INC 초과, 비표현 가능 시간, 마지막 factorization 변경 및 실제delta 반력/수지 |
| 열저장/하중 | 단열 균일 온도 유지, no-BC transient, P 및 -P, overlapping source load 합산, constrained DOF에 P 유지, disconnected component별 보존 |
| Partition | 0×0 all-constrained generic algebra; 공개 내부 station BC 금지 유지 |
| 면내 해석해 | 단열 직사각형 T=T0+A*cos(pi*x/L)*exp(-alpha*(pi/L)^2*t), alpha=k/(rho*c); 온도와 미분 flux, 면내 L2 온도2차/gradient1차 |
| 시간 해석해 | 같은 고정 spatial discretization의 독립 eigenpair Kv=lambda Cv 및 exp(-lambda*t), BE 증폭률; dt refinement의 1차 오차 |
| 두께 slab | 균일 초기온도, bottom step 온도/top 단열 Fourier series; 5점 온도·flux·열에너지·반력, series 절단 bound 및 두께 오차 plateau |
| 재료/두께/접속 | full-end prescribed strip의 sum L/(k*h*w) 저항은 assembly fixture; 실제 bottom/top-only CLI는 해당 조건의 독립 oracle 사용. Fold 등온 및 조립 열률 수지 |
| Geometry | 왜곡/warped patch 적분 비교, scaling/translation, zero/negative/uncertain metric, 자기 교차, 방향 충돌, T edge와 vertex-only junction 음성 사례 |
| IO/실패 | 모든 양성/음성 subset, source identity, HDF5 round trip, 초기 dataset 부재, P=0 반력, 마지막delta, no-op warnings, candidate/lease/write/self-check/close 실패 시 기존결과 보존 |
구현 단계에는 열해석 관련 검사를 먼저 수행하고 전체 MSVC x64 Debug build/CTest와
필수 reference 비교를 수행한다. 새 compiler warning을 추가하지 않는다.
구현 전 전체 검사를 실행하여 기존 HHT 시간수렴 실패 2건의 실제 재현·이름·metric을
기록한다. 기록만으로 현재 실패를 가정하지 않는다. 구현 후 전체 검사에서 제외하지 않고,
새 실패·기존 실패 악화를 허용하지 않는다. 재현된 기존의 무관한 HHT 실패만 별도
미해결로 유지할 수 있으며 전체 CTest green이라고 보고하지 않는다.
예상하지 못한 baseline failure는 자동 면제하지 않고 원인과 영향 범위를 확인한다.
열해석의 필수 검증 실패는 이 예외에 포함되지 않는다.
## STH-C08 — Reference 원본과 완료 evidence
대상 경로는 다음 두 개다.
- `reference/shell-thermal/shell-thermal.inp`
- `reference/shell-thermal/shell-thermal NodalTemperature.csv`
구현 preflight에서 동일 저장소의 기존 dev worktree에 있는
동일 상대 경로에서 byte 그대로 복사한다. `git worktree list --porcelain`으로 dev
worktree를 찾고 파일별 SHA256을 아래와 대조한다. 이 방식은 머신 절대경로를 계약에
고정하지 않는다. 원본이 달라졌거나 없으면 임의 생성·보정·대체하지 않는다.
| 파일 | 확인한 SHA256 |
| --- | --- |
| shell-thermal.inp | 955C8BEEC35CE53DA5F5DD48785CF8C12FCE726EBFDE7EAE40CF69DBF00D6EA4 |
| shell-thermal NodalTemperature.csv | 765B0D0F64C24CDE6D033D6EF80E0267275A496C97E28D932466EB877810D3E5 |
목적 파일이 이미 있으면 hash 일치를 확인해 재사용하고 다르면 덮어쓰지 않는다.
복사 전후 hash와 사용한 경로를 evidence에 남긴다. Abaqus를 실행하거나 원본 파일명,
내용, metadata를 수정하는 권한은 포함하지 않는다. FESA HDF5/log는 복사한 입력 옆에
생성한다. 보조 deterministic CSV view는 공식 결과가 아니다.
외부 family는 `(shell-thermal, source step Step-1, time 37, bottom temperature NT11,
scalar, user-consistent temperature units, blocking)`이다. 원본은 단일 Part-1-1,
676 DS4, 729 node이며 CSV의 `Step,Frame,Node Label,NT11`에서 instance 열 부재는
입력의 단일 instance 및 label 집합 검증으로만 복원한다. HDF5의 source_step_name,
정확한 t=37 accepted frame, source DOF11을 일대일 대응한다. Nearest-time matching,
station 평균, 행 순서 추정, 잘못된 instance 추정은 금지한다.
Missing/extra/duplicate/nonfinite/identity 오류는 tolerance 전에 실패다.
[공통 tolerance](../TOLERANCE.md)의 S=max|reference|, B=0.01S, near-zero absolute<=B,
그 밖의 relative<=0.05, family RMS/S<=0.01, S=0 exact-zero 규칙을 그대로 적용한다.
관찰된 min=-15.5742216, max=1000은 정적 관찰이며 solver 통과 evidence가 아니다.
음수·미소값을 clamp하거나 온도 기준·단위를 조정해 비교를 통과시키지 않는다.
필수 비교 실패 시 원인을 분리 조사하고 방법/acceptance 충돌은 사용자 결정으로 돌린다.
완료 보고서는 구현 diff, 실제 검증 명령/환경·결과, 독립 numerical/physics 검토,
외부 비교 row/family 판정, 기존 HHT 실패, 남은 물리 한계를 구분한다.
추가 Abaqus CSV 확보는 필수 조건이 아니지만 다른 물리량 검증을 생략하는 근거도 아니다.
## 현재 checklist
- [x] 물리·요소·재료·두께·초기 온도·경계조건·열입력 범위 인터뷰
- [x] 고정 시간 증분·모든 증분 결과 저장·필수 결과 선택
- [x] 기존 reference 원본 정적 확인과 외부 비교 범위 합의
- [x] 핵심 문서 3개에 합의와 미결정 항목 기록
- [x] 단일 transient step으로 범위 확정
- [x] 열전도·열용량행렬, 열유속·반력, 시간 적분과 코드 연결을 포함한 상세 계획 제안
- [x] [독립 수치 검토](numerical-review.md)의 문서 지적 반영; 미결정 물리·I/O 계약은 별도로 유지
- [x] 상세 계획을 바탕으로 지역 FEATURED_AGENTS·PRD·ADR 구성
- [ ] 사용자 검토 후 requirement baseline 확정
- [ ] 사용자 검토 후 수치 방법·I/O 계약·geometry 허용 범위와 실행 계획 확정
- [x] 요구사항 인터뷰와 단일 uncoupled step 범위 확정
- [x] Q4×P2 두 구간, consistent C, 3×3/exact Z, BE 채택 승인
- [x] Geometry 접속 이상화, raw 초기값, 마지막 짧은 증분, 입력·HDF5 계약 승인
- [x] 모서리 분포 입력 011 연기 및 reference 복사 승인
- [x] 실행 계획 승인 및 구현 요청
- [x] 정식화 정합과 상세 수치 기준의 [독립 검토](numerical-review.md) 완료
- [x] Production 변경 전 전체 baseline과 reference hash evidence 기록
- [ ] Thermal semantic model·DOF·geometry·K/C·flux 구현 및 검증
- [ ] Thermal 입력 subset·진단·결정적 조립·시간 적분 구현 및 검증
- [ ] HDF5·CLI·candidate와 기존 결과 보존 검증
- [ ] 해석해·수지·수렴성과 필수 NT11 외부 비교 통과
- [ ] 전체 MSVC Debug build/CTest 및 기존 실패 대비 확인
- [ ] 독립 physics 검토와 완료 evidence 기록
Solver build/test, Abaqus 실행, FESA/reference 수치 비교와 실제 결과의 physics 검토는
수행하지 않았다. [독립 수치 검토](numerical-review.md)는 두께 Z 행렬, 평판 20 DOF 행렬과
두 증분 BE 에너지 수지를 정확 유리수로 검산하고 잔여 물리·I/O 결정을 기록했다.
미실행 검사를 통과로 표시하지 않는다. 전문 문서의 수치 유도·검산, 실제 solver 검사,
외부 reference 비교와 physics 검토는 서로 다른 evidence다.
@@ -1,29 +1,29 @@
# Shell Transient Heat Transfer 상세 설계·구현 계획안
- 상태: 사용자 검토용 제안. 요구사항 인터뷰를 바탕으로 수치 방법론을 구체화했으며 구현은 시작하지 않았다.
- 상태: 사용자 승인 실행 계획. 구현·검증 진행 중이며 수치 통과 evidence와 구별한다.
- Owner: 메인 assistant
- 기준: [Feature PRD](PRD.md), [프로젝트 지침](../../AGENTS.md), [ARCHITECTURE](../ARCHITECTURE.md),
[CODINGSTYLE](../CODINGSTYLE.md), [TOLERANCE](../TOLERANCE.md)
- 수식의 상세 원본: [정식화 제안](formulation.md). 현재 checklist는 [PRD](PRD.md#현재-checklist)에만 둔다.
- 수식의 상세 원본: [정식화](formulation.md). 현재 checklist는 [PRD](PRD.md#현재-checklist)에만 둔다.
- 독립 검토: [numerical-review.md](numerical-review.md). 수식의 일관성 검토이며 구현·정확도 통과 판정은 아니다.
## 1. 구현할 해석과 제안 방법
## 1. 구현할 해석과 채택 방법
입력의 4절점 Shell mesh에서 온도만 풀며, 한 절점에 두께 방향 온도 5개를 둔다.
요소당 20개 온도 자유도로 면내·두께 방향 전도와 열저장을 함께 계산한다.
기존 MITC4의 변위·회전 강성행렬을 열행렬로 전용하지 않는다.
| 결정할 항목 | 우선 제안 | 선택 이유와 대안 |
| 항목 | 채택 방법 | 선택 이유와 한계 |
| --- | --- | --- |
| 면내 보간 | Bilinear Q4 | 승인된 4절점 topology와 연결되며 상수 온도·평판 선형 온도장을 재현 |
| 두께 보간 | 2개 구간의 C0 2차 Lagrange 보간 | 5개 온도점을 사용하고 중앙 온도점을 공유. 전체 두께 단일 4차 보간보다 구간별 gradient와 적분을 분리하기 쉬움 |
| geometry | 실제 3차원 midsurface의 metric과 법선을 이용하는 thin-shell product-domain | 접힌 판·곡면을 같은 scalar 열전도 문제로 구성. 두께 offset에 따른 곡률 보정은 생략하는 모델임을 명시 |
| 열용량 | Consistent heat-capacity matrix | Galerkin 약형에서 직접 유도하고 열저장·반력·에너지 계산에 같은 행렬 사용. Lumped 대안은 초기 급변의 진동 및 reference 불일치 조사 대상으로만 남김 |
| 열용량 | Consistent heat-capacity matrix | Galerkin 약형에서 직접 유도하고 열저장·반력·에너지 계산에 같은 행렬 사용. Lumping·인공 확산·온도 clamp는 승인 범위 밖 |
| 적분 | 면내 3×3 Gauss, 두께의 두 P2 구간은 정확 적분한 Z 행렬 사용 | 직접 54점 적분과 대조. Affine 평판은 2×2 면내 적분도 정확하며 곡면·왜곡은 고차 적분 대조로 오차 확인 |
| 시간 적분 | Backward Euler, 고정 Δt | 일정 계수 선형 문제에서 factorization을 재사용할 수 있고 첫 구현의 시간 수렴 검증이 명확 |
| 열유속 | 보간 온도의 기울기에서 직접 복구 | 요소별·위치별 Fourier 열유속과 전체 residual 반력을 구별하여 검증 |
이 표는 FESA가 채택할 제안이며 Abaqus 내부 알고리즘 동등성의 약속이 아니다.
이 표는 사용자가 채택한 FESA 방법이며 Abaqus 내부 알고리즘 동등성의 약속이 아니다.
특히 consistent C가 현재 reference의 과도응답을 공통 tolerance 안에서 재현하는지는
아직 검증하지 않았다. 시간 적분의 안정성이 양의 온도 보존이나 정확도를 보장하지도 않는다.
@@ -73,9 +73,9 @@ product-domain 모델의 교차항은 0이다. 평판에서 각각 대략 `kh`,
`G`의 퇴화, 뒤집히거나 자기 교차하는 요소, 불일치하는 surface side는 입력 geometry
진단으로 처리한다. `|det J|`로 invalid 요소를 정상화하거나 임의의 최소 두께를 대입하지 않는다.
중앙면 접속은 같은 normalized thickness station을 공유하는 이상화로 제안한다.
중앙면 접속은 같은 normalized thickness station을 공유하는 이상화다.
두께가 다른 요소의 단차·급한 fold의 실제 3D 접합부 열장은 이 모델의 해상 범위 밖이다.
이 접속 의미와 얇은 Shell 곡률 근사는 사용자 검토가 필요한 물리 선택이다.
이 접속 의미와 얇은 Shell 곡률 근사는 사용자가 채택한 물리 모델이며 상세 3D 전도를 보장하지 않는다.
## 4. 열용량행렬 C와 적분 절차
@@ -105,22 +105,15 @@ Affine 평판은 면내 2×2 × 두께 6점, 총 24점으로도 정확하다. 5
## 5. 열입력과 지정 온도
열입력의 양의 부호는 모델 안으로 유입되는 방향으로 제안한다.
공개 입력은 source DOF 11/15에 대한 일정한 절점 집중 power P만 지원한다.
양수는 모델로 공급되는 열이며 h·면적·shape weight를 곱하지 않는다. Source row 순서로
더하고 constrained DOF의 입력도 full Q에 유지한다. Edge 분포 입력 011은 연기되었으며
이번에 kernel이나 keyword를 추가하지 않는다. 지정되지 않은 외부 경계는 자연 단열이다.
- **절점 집중 열입력** `P`는 단위 시간당 에너지다. 선택한 node의 아랫면 또는 윗면 자유도에
직접 더하며 면적·두께를 다시 곱하지 않는다.
- **가장자리 열입력**은 외곽 edge의 옆면에 작용하는 면적당 열유속 `q_in`으로 제안한다.
`Q_e=∫edge∫[-h/2,h/2] Hᵀ q_in dz ds`로 조립한다. 직선 edge 길이가 `ℓ`이면
전체 열률은 `q_in hℓ`다. 두께 방향으로 균일해도 5개 온도점에 단순 1/5씩 배분하지 않는다.
- **지정 온도**는 elimination으로 적용한다. 구속 온도점도 C·K·외력 full-space에서 유지하여
반력 열률과 전체 에너지 수지를 복구한다.
- 열입력이나 지정 온도가 없는 외부 경계에는 자연 경계조건인 단열 조건을 적용한다.
실제 source 문법은 지원 의미를 검증한 뒤 고정한다. 특히 DS4의 가장자리 열유속을
`*DFLUX, S1...S4`로 받을 수 있다고 가정하지 않는다. 확인한 공개 shell library의
element-based flux 목록은 BF·SNEG·SPOS이며, 가장자리 입력을 위한 공개 문법 대응은
아직 확정하지 못했다. 이 항목은 계산 kernel과 별개의 I/O 결정이다.
[공개 shell load 정의](https://docs.software.vt.edu/abaqusv2025/English/SIMACAEELMRefMap/simaelm-r-shellgeneral.htm).
초기 온도는 node/Nset별 signed finite 값으로 5점에 동일 적용하고 미지정 값은 0이다.
BC 생략값은 0이며 model 단계와 step 단계 각각에서 중복 충돌을 검사한 뒤 step이
동일 target의 model 값을 덮어쓴다. 동일 값 중복은 허용하고 다른 값은 오류다.
정확한 subset·미지원 옵션·no-op warning은 [PRD STH-C02](PRD.md#sth-c02--공개-입력-충돌과-부가-keyword)를 따른다.
## 6. Transient equation과 해법
@@ -144,15 +137,17 @@ free/constrained partition에서
\]
를 푼다. 초기 구속점 온도와 step 지정 온도가 다르면 첫 증분의
`C_fc(T_c^0−g)/Δt` 항을 빠뜨리면 안 된다. 제안은 입력 초기 상태를 `t=0`에 보존하고
`C_fc(T_c^0−g)/Δt` 항을 빠뜨리면 안 된다. 입력 초기 상태를 `t=0`에 보존하고
첫 양의 시점부터 step 지정 온도를 부과하는 것이다. 초기 값을 덮어쓴 다른 해석과
에너지·reference 결과가 달라질 수 있으므로 이 시작 의미를 명시적으로 검토한다.
에너지·reference 결과가 달라지므로 raw 초기값을 조용히 projection하지 않는다.
고정 물성·geometry·구속 집합·Δt에서는 A_ff를 한 번 factorize하고 매 증분 RHS만 갱신한다.
전체 시간이 Δt로 나누어떨어지는 입력을 우선 지원하는 것으로 제안한다. 마지막 증분을
조용히 줄이지 않는다. 시간은 누적 덧셈보다 increment index와 Δt에서 결정적으로 생성한다.
전체 시간이 Δt로 나누어떨어지지 않으면 마지막 증분만 종료까지 단축한다. Δt가 전체 시간보다
크면 한 증분이다. 시간은 n*Δt로 생성하고 마지막 저장 시점은 입력 종료 시간이다. 표현 오차
수준으로 같은 마지막 간격은 같은 계수를 사용하고 실제로 짧은 간격만 재분해한다.
시간격자·INC·표현 가능성의 상세 판정은 정식화와 독립 numerical review에서 고정한다.
절연체의 K는 상수 온도 null mode를 가지지만 양의 C 때문에 A는 양정정일 수 있다.
절연체의 K는 상수 온도 null mode를 가지지만 양의 C 때문에 A는 양정치다.
정적 solver처럼 “온도 구속이 없으므로 singular”로 판정하지 않는다. 일반 partition 연산은
전체 구속 시 0×0 free solve를 건너뛰도록 검산한다. 다만 이번 공개 입력은 bottom/top만
구속하므로 내부 3개 온도점까지 모두 구속한 사례는 CLI 지원 모델이 아닌 algebra fixture다.
@@ -171,8 +166,7 @@ local 직교 frame `R=[e_1,e_2,n]`이면 `q_local=Rᵀq_global`이고,
- 기본 authoritative flux는 면내 2×2 Gauss 위치와 두께 trace
`[-1+, -1/2, 0−, 0+, 1/2, 1−]`의 곱인 24개 위치에서 복구한다. 두께 trace는 구간
안쪽에서 평가하는 위치/접근 방향이다. Global 성분과 local frame 또는 local 성분의
저장 선택을 I/O에서 고정하여 좌표 의미를 보존한다.
안쪽에서 평가하는 위치/접근 방향이다. Global 성분과 local frame을 저장하고 local 성분은 Rᵀq로 후처리한다.
- 위 trace 위치에서 B를 직접 평가한다. C0 보간의 중앙면에서 두께 도함수는 한쪽 극한
두 개일 수 있으므로 `0−`와 `0+`를 평균으로 합치지 않는다. 행렬의 적분 위치와 복구 위치는 독립이다.
- 재료 경계와 요소 경계를 넘어 nodal averaging한 값을 공식 heat flux로 사용하지 않는다.
@@ -187,15 +181,15 @@ local 직교 frame `R=[e_1,e_2,n]`이면 `q_local=Rᵀq_global`이고,
free 성분은 평형 잔차, constrained 성분은 지정 온도 유지에 필요한 반력 열률이다.
양의 반력은 모델에 추가로 공급한 열이다. transient 반력에는 열저장 항이 포함되어야 한다.
초기 온도는 저장하지만 유한 차분 시간구간이 없는 `t=0` 반력에 계산된 값처럼 0을 쓰지 않는다.
초기 온도는 저장하지만 `t=0` flux와 reaction dataset은 생성하지 않는다. 반력을 계산된 0이나 NaN으로 대체하지 않는다.
## 8. 기존 코드와 연결할 책임·파일
현재 `Element::Recover`는 displacement를 소비하고 `DofManager`는 `kDofsPerNode=6`을
전제로 한다. 따라서 parser에 DS4만 추가해서 완료할 수 없다. 다음은 실제 저장소 경로와
새 파일의 제안 배치다. 세부 C++ signature는 승인된 I/O·정식화 계약을 기준으로 구현 계획에서 고정한다.
신규 파일 배치다. 세부 C++ signature는 승인된 I/O·정식화 계약과 기존 객체 책임을 따른다.
| 작업 경계 | 기존 수정 대상 | 신규 파일 제안과 책임 |
| 작업 경계 | 기존 수정 대상 | 신규 파일과 책임 |
| --- | --- | --- |
| Semantic model | `include/fesa/model/domain.h`, `src/fesa/model/domain.cpp`, `include/fesa/model/model_types.h`, `include/fesa/materials/material.h` | `include/fesa/materials/isotropic_thermal_material.h`, `src/fesa/materials/isotropic_thermal_material.cpp`: 양의 k·ρ·c와 source identity |
| Thermal element | `include/fesa/elements/element_definition.h`, `src/fesa/elements/element_factory.cpp`의 definition/runtime 생성 경계 | `include/fesa/elements/thermal_shell.h`, `src/fesa/elements/thermal_shell.cpp`: 20×20 K/C와 직접 flux 복구. 열 kernel을 구조 Element의 displacement 계약에 억지로 넣지 않음 |
@@ -212,49 +206,55 @@ Domain이 정의를 소유하고 AnalysisModel이 non-owning view를 제공하
## 9. 구현 순서와 각 단계의 완료 조건
아래 순서는 향후 구현 작업을 위한 계획이다. 현재 실행하거나 완료한 목록이 아니다.
아래 dependency는 승인된 실행 계획이다. 같은 행의 독립 작업만 병행하며 완료 evidence 없이
후속 수치 kernel로 진행하지 않는다. 공통 CMake 등록은 메인이 각 묶음 통합 시 소유한다.
| 순서 | 작업·연결 요구사항 | 추가/변경할 검증과 단계 완료 조건 |
| --- | --- | --- |
| 1 | 물리·I/O 선택 확정 — 001~018 | Thin-shell 접속·초기 jump·edge syntax·flux identity의 결정을 문서화. 독립 수치 검토 지적 해결 |
| 2 | Thermal semantic model과 input — 001~013, 018 | `tests/unit/io/abaqus/thermal_domain_test.cpp`: DS4·물성·5점 단면·initial/boundary/CFLUX·single step mapping, 잘못된 DOF·불량 물성·상충 집합·미지원 옵션 진단 |
| 3 | Element kernel — 002~007, 014 | `tests/unit/elements/thermal_shell_test.cpp`: H partition of unity, gradient sum zero, K/C 대칭성·rank, exact flat matrices, plane/thickness linear field와 flux, rotation covariance |
| 4 | DOF·assembly·load — 006~011 | `tests/unit/fem/dof_manager_test.cpp`, `tests/unit/assembly/thermal_assembler_test.cpp`: 20 DOF scatter, 인접 node station 공유, edge 총 열률, source 순서, 병렬/직렬 동등성과 structural 회귀 |
| 5 | Transient procedure — 008~013, 018 | `tests/unit/analysis/transient_heat_procedure_test.cpp`: 두 증분 직접 계산, C_fc 초기 항, no-Dirichlet, generic all-constrained algebra fixture, factorization 재사용, 후보 실패 보존 |
| 6 | Recovery·HDF5·CLI — 014~015 | `tests/unit/io/hdf5/thermal_hdf5_results_writer_test.cpp`, `tests/integration/app/thermal_application_test.cpp`: 모든 온도점·증분·flux 위치·reaction의 round trip, 초기 반력 미정의, 파일 finalize 실패 시 이전 결과 보존 |
| 7 | Physics·reference — 016~017 | `tests/integration/analysis/thermal_physics_test.cpp`, `tests/reference/thermal_reference_comparison_test.cpp`, `tests/reference/thermal_reference_cases_test.cpp`: 아래 사례와 기존 CSV 비교·identity 음성 검사 |
| 8 | 통합 확인 | `src/fesa/CMakeLists.txt`, `tests/CMakeLists.txt` 등록 후 MSVC x64 Debug 전체 build/CTest와 필수 reference 비교. 신규 warning 없음 |
| 작업 | 선행 작업 | 소유 대상·목적 | 검증 |
| --- | --- | --- | --- |
| T1 | 없음 | 메인: 지역 지침/PRD/ADR/계획; 정식화 담당: formulation; 독립 numerical 담당: numerical-review | 전체 patch·roundoff·시간격자·수렴·slab 오차의 구체적 기준 확정 |
| T2 | 없음 | 구현 담당: material/model/Domain/AnalysisModel의 thermal 의미·소유권 | 양성/불량 값·source identity·structural 회귀; parser/kernel/DOF 제외 |
| T3 | T1,T2 | 구현 담당: thermal geometry·20 DOF K/C·24 trace kernel·Domain geometry 연결·공유 시간격자 계산 | Exact M/S/Z, 직접54점, rank/SPD/null, scaling·rotation·trace·geometry 음성·시간격자 경계 |
| T4 | T2 | 구현 담당: DofManager thermal layout/scatter/pattern | 구조6·열5 분리, 활성 node와 stable 번호, partition·공유 node |
| T5 | T3 | 구현 담당: Abaqus thermal subset·semantic mapping·진단 | 승인 input/default/override/no-op와 모든 미지원·충돌 사례 |
| T6 | T3,T4 | 구현 담당: deterministic assembly·thermal state/result record | K/C·source-ordered Q·초기 상태·순서와 공유 recovery 인터페이스 |
| T7 | T6 | 구현 담당: BE procedure·복구·full residual·candidate | Exact 두 증분·no BC·짧은 마지막·factor reuse·수지·실패 state 보존 |
| T8 | T6 | 구현 담당: thermal HDF5 writer·self-check | Synthetic state round trip·schema/identity/초기 dataset 부재·P=0·atomic 실패 |
| T9 | T5,T7,T8 | 구현 담당: prepared_analysis와 CLI 연결 | 양성 CLI·lexical output·exit·lease·기존 결과 보존 |
| T10 | T9 | 구현 담당: analytical physics·reference comparator/integration | 독립 oracle·온도/flux/reaction·수렴과 NT11 source mapping/tolerance |
| T11 | T10 | 검증 담당·독립 physics 검토자; 메인 완료 문서 | 전체 MSVC Debug build/CTest·필수 reference·baseline 대비·coverage 한계 |
각 행동 변경은 먼저 기대 결과를 고정하는 검사를 작성하고 실패를 확인한 뒤 구현한다.
기존 coverage가 충분한 지원 코드 정리에는 구현 구조만 따라 하는 테스트를 추가하지 않는다.
독립 수치 검토와 실제 결과의 물리 검토를 구현 담당자의 자기 승인으로 대체하지 않는다.
행동 변경은 기대 결과를 고정한 검사에서 실패를 확인한 뒤 구현한다. 구현 구조를 복제한
expected 값으로 정확도를 주장하지 않는다. 독립 numerical/physics 검토를 구현자의 자기 승인으로 대체하지 않는다.
각 담당은 할당 파일만 수정하고 source/reference artifact와 공통 프로젝트 계약을 임의 변경하지 않는다.
입력의 INC 사전검사와 procedure는 같은 시간격자 계산을 사용해야 한다. T1에서 검토한
시간격자 함수를 T3의 공통 thermal validation 경계에서 구현하고 T5와 T7이 재사용한다.
이 파일 소유권 배분은 승인된 dependency를 유지하며 parser와 procedure에 같은 수식을
각각 구현하거나 아직 존재하지 않는 후속 procedure에 parser를 의존시키지 않기 위한 것이다.
## 10. 방법론을 반증할 검증 사례
| 사례 | 기대되는 관찰과 오류 검출 목적 |
| --- | --- |
| 완전 단열·균일 초기 온도 | 모든 온도점이 시간에 따라 불변, flux=0, 열에너지 불변. K의 상수 null mode 검증 |
| 평판의 면내 선형 온도 | Kernel에서 `q=−k∇T`를 정확히 재현. Global 회전 시 scalar T 불변·vector q 회전 |
| 두께 방향 선형 온도 | 하·상면 온도차의 정상 한계에서 `q_n=−k(T_top−T_bottom)/h`, 두 면 반력 합 0 |
| 두께 전체 단열의 면내 Fourier mode | 초기 온도가 두께 방향 균일한 허용 입력으로 analytical decay와 비교. BE 시간 1차 수렴과 면내 mesh 수렴 분리 |
| 한쪽 면 온도 step, 반대 면 단열 slab | 허용된 초기 조건에서 1D transient analytical series와 비교. 5점 두께 근사의 한계와 초기 jump 처리 검증 |
| 두 재료·두 두께 strip | 단면 전체 end 온도를 구속한 assembly fixture에서 총 열저항 `Σ L/(k h w)` 확인. 내부 station 구속은 검산 전용이며 CLI 기능이 아님. 공개 bottom/top 경계 사례는 경계층을 포함한 별도 oracle과 비교 |
| Edge flux와 집중 열입력 | 입력한 총 열률과 `Δ(1ᵀCT)/Δt` 및 reaction 합의 일치. 단위와 면적/두께 중복 곱 검출 |
| 작은 열망의 직접 해 | `A_ff`, RHS, 첫 두 증분, full residual을 독립 hand calculation과 비교. 전체 구속 특수 조건은 generic partition fixture로 한정 |
| 실제 reference | 원본 `.inp`와 t37의 NT11 729행 비교. 누락/추가/중복/nonfinite·instance ambiguity를 tolerance 전에 거부 |
필수 사례의 원본은 [PRD STH-C07](PRD.md#sth-c07--독립-검증과-수치-계약의-선행-조건)이다.
정확 fixture·norm·시간/mesh 목록·oracle 절단 및 order 판정은 [독립 수치 검토](numerical-review.md)가 소유한다.
전체 열에너지는 `E=1ᵀCT`로 두고 기준 온도의 상수 차이가 에너지 변화에 영향을 주지 않음을
확인한다. 검증할 discrete balance는
`E_(n+1)−E_n = Δt[1ᵀQ + Σ constrained reaction]`이며 free residual의 영향도 보고한다.
점별 Fourier flux와 FE 반력은 다른 observable이므로 각각의 검증을 수행한다.
- Shape partition, affine exact K/C, constant null mode, C SPD, 양의 에너지와 회전 변환.
- 일정 온도·면내 선형·두께 선형 flux 및 중앙면 0-/0+를 구별하는 직접 미분.
- 두 자유도 C=[[2,1],[1,2]], K=[[1,-1],[-1,1]], delta=.5, old=0, 두 번째 g=1에서
첫 free=-.2, 두 번째 .04, reaction=4.8와1.44. 첫 C_fc 이전 상태 항을 검출한다.
- No-BC 균일 보존, P/-P·constrained P·source 중복과 component별 에너지 수지.
- Semidiscrete 독립 eigenmode로 시간1차, Fourier 직사각형의 L2 온도2차·gradient1차.
- Bottom step/top 단열 slab series의 온도·flux·에너지·반력; 고정5점 두께 오차 plateau와 초기 침투 깊이 보고.
- 다중 재료/두께 strip 저항은 full-end assembly fixture로 확인한다. 내부 station BC를 CLI로 확대하지 않는다.
- 원본 Step-1 t37 NT11 729행의 identity precheck와 공통 tolerance blocking 비교.
시간 수렴은 공간 오차를 충분히 낮추고 Δt를 줄이며, 면내 수렴은 충분히 작은 Δt에서 mesh를
세분화한다. 제품의 두께 온도점 5개가 고정되므로 면내 mesh를 줄이는 것만으로 두께 오차까지
0으로 간다고 주장하지 않는다. 세부 norm·기대 order·내부 판정 기준은 독립 검토가 정한다.
에너지 defect의 scale은 저장·외력·반력 항의 절댓값 합이며 free residual을 별도로 검사한다.
E=1ᵀCT의 기준 온도 이동에 대한 DeltaE 불변성과 BE C-norm 소산은 서로 다른 검사다.
외부 5%/1% 기준으로 내부 수지를 판정하거나 max(1,scale)로 작은 모델의 오류를 숨기지 않는다.
## 11. 검증 실행과 완료 판단
구현 시 MSVC x64 developer 환경과 GoogleTest/MKL/TBB/HDF5 경로를 준비하고 다음을 실행한다.
Production 변경 전에 전체 baseline을 실행하고 실제 실패 이름·metric·환경을 기록한다.
구현 중 관련 검사를, 기능 완료 시 아래 전체 검사를 실행한다.
```powershell
cmake --preset windows-msvc-ninja-debug
@@ -263,39 +263,25 @@ ctest --preset windows-msvc-ninja-debug --show-only=json-v1
ctest --preset windows-msvc-ninja-debug --output-on-failure
```
선택된 `.harness/config.json`이 동일 실행을 구성하면 `python scripts/execute.py`를 사용할 수 있다.
현재 계획 작성에서는 이 개인 설정 파일을 확인하지 못했으므로 runner의 preset 선택을 가정하지 않는다.
Python Harness를 수정하지 않는 이 기능에 Harness Python suite를 형식적으로 추가 실행하지 않는다.
MSVC x64 developer 환경과 승인 local GoogleTest/MKL/TBB/HDF5를 사용한다.
`.harness/config.json`이 동일 실행을 구성할 때만 `python scripts/execute.py`를 대안으로 사용한다.
개인 설정이 없다고 runner가 preset을 선택한다고 가정하지 않는다. 새 compiler warning을 추가하지 않는다.
Reference 원본은 현재 별도 dev checkout에서 읽었으며 현 worktree에는 없다. 실행 단계에서
read-only reference 위치와 결과 output 경로를 명시적으로 연결한다. 파일을 찾지 못하면
comparison을 통과나 coverage 완료로 기록하지 않는다. reference 재생성·이름 변경·보정은 포함하지 않는다.
Reference는 [PRD STH-C08](PRD.md#sth-c08--reference-원본과-완료-evidence)의 두 원본을
승인된 dev checkout에서 byte 그대로 복사·hash 확인한다. Abaqus 실행·재생성·이름 변경·보정은 제외한다.
실제 명령/exit, hash, baseline 및 최종 결과는 implementation-verification 문서에 기록한다.
완료는 solver 실행 성공뿐 아니라 모든 필수 결과·증분 저장, 에너지 수지·해석해·수렴,
기존 구조해석 회귀, 공통 tolerance의 NT11 비교와 필요한 독립 검토 evidence가 있어야 한다.
재현된 무관한 기존 HHT 시간수렴 실패는 전체 검사에서 제외하지 않고 마지막에도 비교한다.
새 실패나 기존 실패 악화를 허용하지 않으며 기존 실패가 남으면 전체 CTest green이라고 보고하지 않는다.
열해석의 필수 physics/reference 실패에는 이 예외를 적용하지 않는다.
현재 문서 검토에서는 독립 reviewer가 Z 계수, 평판 20 DOF의 C 양정치/K rank 19,
첫 두 증분의 C_fc 항·반력·에너지 수지를 정확 유리수 연산으로 확인했다.
구현 코드·FESA 결과의 정확도나 곡면·접합부 지원 한계가 검증된 것은 아니다.
## 12. 채택 모델의 한계
## 12. 구현 전에 검토할 제안과 한계
Product-domain은 실제 offset 곡률과 상세 접합부를 해상하지 않는다. 고정5점은 초기 두께
경계층을 충분히 해상하지 못할 수 있고 면내 refinement가 두께 오차를 제거하지 않는다.
Consistent C/BE의 안정성은 positivity나 정확도를 보장하지 않는다. Point source 부근의
mesh 독립적인 유한 flux는 보장하지 않는다. 이 한계를 임의 clamp나 tolerance 완화로 숨기지 않는다.
1. Q4×두 구간 P2, consistent C, 면내 3×3/두께 exact 적분, Backward Euler의 선택
2. Thin-shell product-domain과 normalized thickness station 접속; 실제 접합부/곡률 두께 효과 한계
3. 입력 초기 상태를 보존하는 첫 증분 경계 jump와 t=0 반력 미정의
4. 가장자리 열입력의 면적당 단위와 승인 가능한 외부 입력 문법
5. 모든 온도점 온도·3성분 flux·반력의 위치/좌표계와 HDF5 저장 계약
6. 고정 Δt로 나누어떨어지는 종료 시간만 우선 지원하는 제한
이 항목은 사용자 인터뷰에서 이미 승인한 범위와 새 수치·I/O 제안을 구별하기 위한 것이다.
현재 문서는 상세 방법론을 검토할 수 있는 계획안이며, 모든 구현 계약이 확정되었다는 선언이 아니다.
Reference에는 step 전 `*Boundary / Set-10, 11, 11`과 step 내부
`Set-10, 11, 11, 1000.`이 함께 있다. Step 전 경계를 지우거나 무시하지 않고 model/history
정의로 구분하여 읽는다. 생략 초기 온도·생략 boundary 값의 default, step에서 같은 대상에
새 값을 적용하는 의미를 원문과 대조해야 최초 C_fc 항과 reference 비교를 확정할 수 있다.
CSV 정적 통계의 NT11 최솟값은 `−15.5742216`, 최댓값은 `1000`이다. 공통 scale `S=1000`,
near-zero band `B=10`에서 위 음수는 near-zero 밖이므로 상대오차 검사 대상이다.
이를 0으로 보정하거나 단순 roundoff로 설명하지 않는다. 초기 조건과 열용량·보간 선택에
따른 undershoot의 원인을 조사하고 analytical accuracy와 reference 일치를 각각 평가한다.
Reference의 NT11 최소값 -15.5742216은 S=1000, B=10의 near-zero 범위 밖이다.
공통 상대오차 검사를 유지하고 초기 jump·공간 이산·열용량의 영향을 독립적으로 평가한다.
현재 채택 계약을 바꾸어야 하는 충돌이 확인되면 근거를 제시하고 사용자 결정으로 돌린다.
+388 -77
View File
@@ -1,17 +1,20 @@
# Shell 과도 열전달 정식화 제안
# Shell 과도 열전달 정식화
- 상태: **Proposed / 사용자 미승인 / 독립 numerical review 별도 보고서 참조**. 구현 실행 지시가 아니다.
- 상태: **2026-09-29 승인된 물리·방법을 반영한 구현용 정식화**.
독립 검토의 판정과 검증 기준은 [numerical-review](numerical-review.md)를 따른다.
- 담당 범위: `FEM-REQ-SHELL-THERMAL-001`–`018`의 수식·물리 의미와 검증 함의.
- 기준: [Feature PRD](PRD.md), [Feature ADR](ADR.md), [프로젝트 지침](../../AGENTS.md),
[프로젝트 PRD](../PRD.md), [공통 ADR](../ADR.md#adr-028), [ARCHITECTURE](../ARCHITECTURE.md).
- 아래 수식은 FESA의 제안이다. Abaqus 내부 요소·적분·열용량·결과 생성과 동등하다는 주장이 아니다.
- 아래 수식은 승인된 FESA 모델이다. Abaqus 내부 알고리즘과 동등하다는 주장이 아니다.
이전 Proposed 보고서의 검산은 제한된 과거 evidence이며 현재 solver 검증을 대신하지 않는다.
## 1. 적용 가정과 차원
고정된 3차원 midsurface의 4절점 Q4 요소에 두께 방향 온도점 5개를 둔다.
요소마다 두께 `h > 0`, 등방성 열전도율 `k > 0`, 밀도 `ρ > 0`, 비열 `c > 0`가
상수이며 단면 전체는 같은 재료다. 서로 다른 요소에는 다른 값을 배정할 수 있다.
변형·대류·복사·내부 발열·적층·면적 분포 표면 열유속은 [PRD 제외 범위](PRD.md#제외-범위)를 따른다.
요소마다 유한한 두께 `h > 0`, 등방성 열전도율 `k > 0`, 밀도 `ρ > 0`, 비열 `c > 0`가
상수이며 centered 단면 전체는 같은 재료다. 요소 집합별 재료·두께를 배정할 수 있다.
단일 uncoupled transient step만 지원한다. 변형·연성·대류·복사·내부 발열·적층·분포
표면 열유속·분포 edge 입력·시간 변화 입력·다중 step은 제외한다.
일관된 입력 단위계를 사용한다. SI 예시는 `h [m]`, `T [K]`, `t [s]`,
`k [W/(m·K)]`, `ρ [kg/m³]`, `c [J/(kg·K)]`, 열유속 `q [W/m²]`이다.
@@ -47,31 +50,117 @@ dV_* = J_s\frac h2\,d\xi\,d\eta\,d\zeta.
정확한 offset-volume의 접선 `aα+z n_,α`, 두께별 metric·면적 변화는 생략한다.
따라서 `h |κ1|`, `h |κ2|`가 작고 열장이 변하는 면내 길이에 비해 두께가 작은 영역을
대상으로 한다. 적분 차수 증가로 이 모델링 오차가 없어지지는 않는다. 곡률의 정량 허용치,
왜곡 허용치와 validation margin은 아직 미정이다. 구조 shell의 두께 비율 기준을 전용하지 않는다.
곡률을 반영하는 별도 shell 열모델이 있다는 점은 [COMSOL의 곡면 열전달 예제](https://www.comsol.com/model/heat-transfer-in-curved-layered-shells-104241)로 확인할 수 있다.
대상으로 하는 물리 근사다. 이는 승인되지 않은 curvature/thickness cutoff를 입력에 추가하는
조건이 아니다. 평판·warped patch·곡면 mesh·fold를 지원하며 적분 차수를 높여도 생략한
offset-volume 모델링 오차는 없어지지 않는다.
Geometry의 필수 조건은 유한 좌표, 서로 다른 4절점, `J_s>0`, 비특이 metric, 자기 교차 없는
유효 patch와 일관된 orientable connectivity이다. 적분점 검사만으로 전체 patch 유효성을
증명하지 않는다. 검사 위치·scale-aware degeneracy 기준·비다양체 허용 범위는 I/O와
numerical review가 확정해야 하며 작은 양수 Jacobian을 clamp하지 않는다.
증명하지 않는다. 아래 whole-patch certificate를 사용하며 작은 양수 Jacobian을 clamp하지 않는다.
**접속 제안:** 인접 요소는 같은 source node의 같은 `ζ` 온도 DOF를 공유한다.
**승인된 접속:** 인접 요소는 같은 source node의 같은 `ζ` 온도 DOF를 공유한다.
매끈한 면에서 이는 온도 연속성과 약형 열률 평형을 준다. 접힌 판에서는 요소별 normal을
유지하고 공유 edge의 `ζ`별 온도만 같게 묶는 이상적 무저항·무용량 접합이다.
서로 다른 두께에서는 같은 normalized `ζ`를 연결하므로 실제 `z` 좌표가 같다는 뜻이 아니다.
실제 step 단차의 노출면, 유한 접합부 부피, 국소 spreading resistance는 해상하지 않는다.
이 모델은 사용자 결정이 필요한 물리 이상화다. 동의 없이 실제 3D 접속을 지원한다고 표시하지 않는다.
이는 승인된 이상화이며 상세 3D 접합부 온도장·열유속의 재현을 의미하지 않는다.
공유 edge의 요소 순서는 서로 반대여야 하고 bottom/top의 연속된 의미를 유지해야 한다.
온도는 scalar이므로 normal 사이의 회전행렬을 곱하지 않는다. 법선 반전은 해당 두께 station
`a ↔ 6-a`와 면 하중/경계조건을 함께 바꿔야 같은 물리 문제가 된다.
### 2.1 Whole-patch geometry certificate
먼저 네 source node identity와 좌표가 각각 서로 다르고 유한함을 확인한다. Bilinear map을
`r=a+bξ+cη+dξη`로 쓰면 다음은 근사가 아닌 항등식이다.
\[
b=(-x_1+x_2+x_3-x_4)/4,\quad
c=(-x_1-x_2+x_3+x_4)/4,\quad
d=(x_1-x_2+x_3-x_4)/4,
\]
\[
w(\xi,\eta)=a_1\times a_2=w_0+u\xi+v\eta,\qquad
w_0=b\times c,\quad u=b\times d,\quad v=d\times c.
\]
따라서 `J_s²=||w||²`의 사각형 최소값은 convex quadratic minimization이다.
정확산술에서 다음 유한 후보 집합으로 최소점을 구할 수 있다.
1. 네 모서리 `(±1,±1)`.
2. `ξ=s=±1`의 edge: `v·v>0`이면
`η=clip(-v·(w0+su)/(v·v),[-1,1])`; `v=0`이면 모서리로 충분하다.
`η=s`의 edge는 `u,v`와 `ξ,η`를 교환한다. 이 clip은 constrained minimizer의
수학적 projection이며 geometry/Jacobian의 clamp가 아니다.
3. `D=(u·u)(v·v)-(u·v)²>0`이면
`[[u·u,u·v],[u·v,v·v]] [ξ,η]ᵀ=-[u·w0,v·w0]ᵀ`의 해가 사각형 안에 있을 때 추가한다.
`D=0`이면 minimum의 affine level set이 경계에 도달하므로 edge 후보로 충분하다.
최소점 `p*`와 `w*=w(p*)`에 대해 convex projection 성질은
`w*·w(p)≥||w*||²`이다. 그러므로 **어떤 비영 벡터 `s`에 대해서든 네 모서리 모두
`s·w_corner>0`를 증명하면 전체 patch의 `J_s>0`가 증명된다.** 역으로 유효 patch는
정확산술에서 `s=w*`가 이 certificate를 준다. 고정 normal에 대한 양의 projection을
가정하는 방법이 아니라 각 patch에서 certificate 방향을 찾는 방법이다.
자기 교차도 함께 배제된다. 임의의 두 자연좌표 `p,q`에 대해 bilinear 항등식
`r(p)-r(q)=Dr((p+q)/2)(p-q)`가 성립한다. 중점은 사각형 안에 있고 `J_s>0`이면
`Dr`의 rank는 2이므로 `r(p)=r(q)`는 `p=q`를 뜻한다. `g`는 Gram 행렬이므로
동시에 positive definite다. 별도 corner-only 또는 quadrature-only 판정은 충분하지 않다.
구현의 부동소수점 판정은 다음처럼 certificate와 후보 탐색을 분리한다.
- IEEE binary64 round-to-nearest, gradual underflow를 전제로 한다. Fast-math,
flush-to-zero 또는 검증하지 않은 연산 재결합을 사용하지 않는다. 좌표 차·곱·합을
powers-of-two로 재척도하여 overflow/underflow를 피한다. 길이 척도는 실제 좌표
차이에서 얻으며 `max(1,scale)` 같은 단위 의존 floor를 두지 않는다.
- 각 입력 double을 singleton interval로 보고 `+,-,×,÷` 각각의 계산 결과를
`nextafter(result,-∞)`와 `nextafter(result,+∞)`로 바깥으로 넓힌다.
곱/나눗셈은 endpoint 조합의 extrema를 사용하며 0을 포함하는 divisor는 사용하지
않는다. Exact powers-of-two 변환도 subnormal 손실이 있으면 enclosure에 포함한다.
이 규칙은 cancellation을 포함해 연산 횟수에 맞는 bound를 전파한다.
- 위 후보를 통상 double로 구해 `s=w(candidate)` 방향을 시도한다. 후보의 작은
오차나 사각형 바깥 위치는 certificate의 건전성을 해치지 않는다. 후보 탐색의 `D`
또는 division이 불확실하면 유한한 다른 후보를 시도한다. 각 `s` 자체는 고정 double
벡터로 취급하고, 원래 좌표로부터 interval로 다시 구한 네 `w_corner`에 대해
`lower(s·w_corner)>0`인지 검증한다. 최종 판정은 후보 최소값의 부호에 의존하지 않는다.
- 모두 양수이면 `J_s≥min_corner lower(s·w_corner)/upper(||s||)>0`라는
whole-patch lower bound가 생긴다. Certificate가 없으면 invalid/uncertain geometry로
fail closed하고 source element와 실패 bound를 진단한다. 임의 aspect/warp cutoff,
고정 epsilon margin 또는 arbitrary-precision dependency는 필요하지 않다.
수치적으로 증명 불가능한 경우를 정확한 기하학적 퇴화라고 단정하지 않는다.
평가점에서도 유한한 metric·gradient를 확인한다. Certificate 통과가 극단적 물성·두께나
metric 역행렬의 유한한 계산을 보장하지는 않으며, 그 실패를 normal 수정으로 숨기지 않는다.
### 2.2 Mesh topology와 source 접속
Source edge key는 정렬한 두 node identity이며 각 incidence에는 원래 directed pair를
보존한다. Incidence는 1 또는 2이고, 2이면 방향이 반대여야 한다. 3개 이상의 요소가
공유하는 edge, 같은 방향의 공유 edge, 중복 element patch는 거부한다. 자동 normal flip이나
좌표 기반 node welding은 하지 않는다. Disconnected components는 각각 유효하다.
각 vertex의 link graph를 검사한다. Link vertex는 해당 node에 닿는 source edge이고,
각 incident Q4는 그 두 edge 사이에 link edge 하나를 만든다. Link는 내부 vertex에서
하나의 cycle, boundary vertex에서 하나의 path여야 한다. 끊어진 link, 분기, 중복 연결은
vertex-only fan/nonmanifold junction이므로 거부한다. Fold의 각도 자체는 거부 이유가 아니다.
Source 공유 접합의 nonconforming T-edge 검사는 같은 instance에서 실제 source endpoint를
공유하는 edge 쌍으로 한정한다. Coarse edge `(p,q)`와 partial-edge 후보 `(p,v)` 또는
`(v,q)`가 존재하고 `p,q,v`가 서로 다른 source node일 때만 `v`를 검사한다. 후보 생성에서
좌표 일치는 source identity 공유를 대신하지 않는다. 각 edge를 coarse 후보로 고려한다.
후보의 `x_v=x_p+λ(x_q-x_p), 0<λ<1` 공선·열린 선분 포함 여부는 scale-aware interval
bound로 판정한다. 선분 밖임이 확정되면 제외하고, 포함이 확인되면 nonconforming T-edge로
거부하며, bound가 판정을 확정하지 못하면 uncertain 접속으로 진단한다. 이 검사는 source
endpoint를 공유하면서 부분적으로 겹치는 선분도 포함하지만 접속을 새로 생성하지 않는다.
모든 mesh node와 모든 edge 사이의 전역 좌표 scan은 수행하지 않는다. 후보 없는 서로 다른
component/instance의 우연 접촉과 같은 component 내 멀리 떨어진 부분의 좌표 접촉은 검사
대상이 아니다. 이는 [PRD STH-C01](PRD.md#sth-c01--모델-geometry와-자유도)의 source 접속
검증이며, patch 내부 injectivity 검사나 미지원 global contact 검출을 대신하지 않는다.
## 3. 20 DOF와 온도 보간
절점별 온도점은 `ζ_a=(-1,-1/2,0,1/2,1)`, `a=1,…,5`이다.
로컬 순서는 `p=5(i-1)+a`, 즉 node-major 20 DOF다. 외부 DOF는 `10+a`로 대응한다.
따라서 `11=bottom`, `13=midsurface`, `15=top`이다. 기하 절점 위치와 DOF11의 물리 위치를
혼동하지 않는다. 이 외부 위치 의미는 [Abaqus Heat Transfer Shell Elements](https://docs.software.vt.edu/abaqusv2025/English/SIMACAEELMRefMap/simaelm-c-shellelem.htm)의 공개 규약을 따른다.
혼동하지 않는다. 외부 위치 의미는 승인된 DS4 입력 계약을 따른다.
면내 shape와 미분은 각 절점의 부호 `(ξ_i,η_i)`를 사용한다.
@@ -81,7 +170,7 @@ N_{i,\xi}=\tfrac14\xi_i(1+\eta_i\eta),\quad
N_{i,\eta}=\tfrac14\eta_i(1+\xi_i\xi).
\]
두께는 `I−=[-1,0]`, `I+=[0,1]`의 **두 C0 quadratic 구간**으로 제안한다.
두께는 `I−=[-1,0]`, `I+=[0,1]`의 **두 C0 quadratic 구간**을 사용한다.
단일 quartic보다 각 구간의 저차 보간·국소 미분을 명료하게 유지하며, 중간 온도 1개를 공유한다.
구간 좌표는 `s=2ζ+1` 또는 `s=2ζ−1`; 두 경우 `s∈[-1,1]`, `ds/dζ=2`이다.
@@ -113,7 +202,7 @@ Fourier 법칙 `q=−k∇_*T`, 저장식 `ρc Tdot + div_* q=0`를 사용한다.
\[
\int_{V_*}\rho c\,v\dot T\,dV_*+
\int_{V_*}k\nabla_*v\cdot\nabla_*T\,dV_*
=\int_{\Gamma_{\rm edge,*}}v q_{\rm in}\,dA_*+\sum_j v_jP_j.
=\sum_j v_jP_j.
\]
지정 온도 외의 하중 없는 경계는 자연 단열이다. 집중 입력 `P_j`는 승인된 bottom/top
@@ -127,7 +216,7 @@ Fourier 법칙 `q=−k∇_*T`, 저장식 `ρc Tdot + div_* q=0`를 사용한다.
`C_e [J/K]`는 consistent 열용량, `K_e [W/K]`는 열전도 행렬이다.
`C_e` 대칭 양정치, `K_e` 대칭 양반정치와 유일한 요소 상수온도 null mode를 기대한다.
Reduced integration, hourglass stiffness, capacity lumping은 이 제안에 넣지 않는다.
Reduced integration, hourglass stiffness, capacity lumping, 인공 확산은 사용하지 않는다.
명시적인 component식과 분리 계산에 필요한 작은 행렬은 다음과 같다.
@@ -151,93 +240,255 @@ Z^1_{\rm seg}=\frac1{3}\begin{bmatrix}7&-8&1\\-8&16&-8\\1&-8&7\end{bmatrix}.
따라서 `C=ρc(h/2)(M⊗Z0)`, `K=k[(h/2)(S⊗Z0)+(2/h)(M⊗Z1)]`이다.
이는 해당 node-major 순서와 일정 두께·물성·product metric에서만 적용한다.
요소/DOF source identity와 stable scatter를 보존하여 full `C,K,F`를 조립한다.
요소/DOF source identity와 stable scatter를 보존하여 full `C,K,Q`를 조립한다.
## 5. 적분과 geometry 오차의 분리
Affine 평판(직사각형·평행사변형)의 일정 metric에서는 면내 `2×2` Gauss가 `M,S`에 정확하다.
두께 각 구간에서 `C`와 면내 전도는 4차 다항식이므로 **3점 Gauss**가 필요하다.
두께 전도는 2차이므로 2점으로 정확하지만 통일된 3점을 권장한다.
두께 전도는 2차이므로 2점으로 정확하지만 직접 검산은 통일된 3점을 사용한다.
이에 따라 affine의 직접 적분은 `2×2×(2 segments×3)=24`점으로 정확하다.
일반적인 평면 왜곡 Q4의 `S`는 역 Jacobian 때문에 다항식이 아니므로 이 exactness는 적용되지 않는다.
제안 기본은 모든 요소의 면내 `3×3` Gauss와 위 exact `Z`의 분리 계산이다.
승인된 production 계산은 모든 요소의 면내 `3×3` Gauss와 위 exact `Z`의 분리 계산이다.
직접 계산하면 `3×3×2×3=54`점이며 같은 행렬이어야 한다. 3점 규칙은
`s=(-√(3/5),0,√(3/5))`, `w=(5/9,8/9,5/9)`이고 각 두께 구간에서 `dζ=ds/2`다.
표면의 warped/curved metric에는 유한 고정 차수의 exactness를 주장하지 않는다.
`3×3 → 5×5 → 7×7` 적분 수렴, mesh 수렴, thin-shell 모델 오차를 별도 확인한다.
지원 geometry 범위에서 기본 차수가 불충분하면 검토를 통해 차수/범위를 고정하고 구현한다.
실행마다 허용오차에 맞춰 임의 adaptive quadrature를 선택하는 계약은 아니다.
3/5/7 차이는 진단용 evidence이며 3×3의 exactness 또는 전역 오차 bound가 아니다.
승인된 검증을 통과하지 못하면 조사·보고한다. 실패 후 지원 geometry 축소, 차수 변경,
adaptive quadrature 또는 tolerance 변경으로 통과시키지 않는다.
## 6. 가장자리·절점 열입력
## 6. 절점 열입력
Edge 순서는 `1:(1,2), 2:(2,3), 3:(3,4), 4:(4,1)`로 제안한다.
면내 자연좌표 edge curve를 `r_edge(u)`, `u∈[-1,1]`, `J_l=|dr_edge/du|`로 둔다.
두께 전체 균일 입력 `q_in [W/m²]`의 일관 하중은 다음과 같다.
\[
(F_e^{\rm edge})_{ia}=\int_{-1}^1\int_{-1}^1
N_i(\xi(u),\eta(u))L_a(\zeta)q_{\rm in}J_l\frac h2\,d\zeta du.
\]
Q4의 직선 edge, 일정 입력·두께에는 총 입력 `q_in h L_edge`가 정확히 보존된다.
각 endpoint의 5성분은 `q_in h L_edge (1,4,2,4,1)/24`이며 나머지 두 절점은 0이다.
이는 `∫L dζ=(1,4,2,4,1)/6`와 edge당 `∫N_i dl=L_edge/2`에서 나온다.
Edge 2점, 각 두께 구간 2점으로 이 하중은 정확하다. 내부 온도점에도 일관 하중이 생기는
것은 허용된 edge 입력의 결과이며, 내부 온도점 직접 집중 입력을 지원한다는 뜻은 아니다.
외부 syntax가 선열률 `λ [W/m]`를 쓰기로 결정되면 `q_in=λ/h`로 해석해야 한다.
두 단위를 혼용하지 않으며 외부 입력 단위·keyword는 I/O owner가 확정한다.
분포 edge 입력은 전부 연기되었다. 이전 초안의 edge 하중식·weight는 현재 kernel이나
검증 acceptance가 아니다. 입력 `Q`는 bottom/top의 signed inward nodal power만 포함한다.
Bottom/top 집중 열입력 `P [W]`는 지정 source node의 `a=1/5`에 그대로 더한다.
`h`, 면적 또는 적분 weight를 추가하지 않는다. 지정 온도와 같은 DOF의 입력도 full `F`에
`h`, 면적 또는 적분 weight를 추가하지 않는다. 지정 온도와 같은 DOF의 입력도 full `Q`에
남겨 반력 계산에 반영한다. 중복 application의 가산/충돌 의미는 I/O 계약을 따른다.
## 7. Backward Euler, 지정 온도와 최초 증분
고정 `Δt>0`에서 backward Euler를 제안한다. `t_n=nΔt`, `A=C/Δt+K`이면
각 실제 증분 `δ_n>0`에서 backward Euler를 사용한다. `A_n=C/δ_n+K`이면
\[
C\frac{T^{n+1}-T^n}{\Delta t}+KT^{n+1}=F^{n+1}+R^{n+1},\qquad
A_{ff}T_f^{n+1}=F_f^{n+1}+\frac{C_{ff}}{\Delta t}T_f^n
+\frac{C_{fc}}{\Delta t}T_c^n-A_{fc}g^{n+1}.
C\frac{T^{n+1}-T^n}{\delta_n}+KT^{n+1}=Q+R^{n+1},\qquad
(A_n)_{ff}T_f^{n+1}=Q_f+\frac{C_{ff}}{\delta_n}T_f^n
+\frac{C_{fc}}{\delta_n}T_c^n-(A_n)_{fc}g.
\]
여기서 `f/c`는 free/constrained, `T_c^{n+1}=g^{n+1}`이다. 동등한 RHS 표현은
`F_f−K_fc g^{n+1}+(C_ff/Δt)T_f^n+(C_fc/Δt)(T_c^n−g^{n+1})`다.
여기서 `f/c`는 free/constrained, `T_c^{n+1}=g`이다. 동등한 RHS 표현은
`Q_f−K_fc g+(C_ff/δ_n)T_f^n+(C_fc/δ_n)(T_c^n−g)`다.
**이전 constrained 온도 항은 첫 증분에서도 생략하지 않는다.** 이후 일정 `g`에서는
마지막 항이 사라진다. 잔차 tangent는 정확히 `A`; Newton 반복은 필요 없다.
Backward Euler의 mass/diffusion 행렬 구성은 [MFEM Example 16 원문](https://raw.githubusercontent.com/mfem/mfem/master/examples/ex16.cpp)에서도 확인할 수 있으며, 위 제약·초기화 정책은 별도의 FESA 제안이다.
마지막 항이 사라진다. 잔차 tangent는 정확히 `A_n`; Newton 반복은 필요 없다.
`C>0`이므로 단열·무구속 모델도 transient `A_ff`는 양정치다. 정상 전도 `K`의 상수온도
null mode를 과도 해석의 singular error로 해석하지 않는다. 일반 partition 대수에서 모든
DOF가 constrained이면 `0×0` free solve를 건너뛰고 prescribed state와 full reaction을 계산한다.
이 경우는 generic partition fixture로만 검증한다. 이번 공개 입력은 bottom/top만 구속할 수
있으므로 내부 3개 station까지 모두 구속하는 열 CLI 모델은 허용되지 않는다.
고정 geometry·물성·DOF 집합·Δt에서는 `A_ff`를 최초 한 번 factorize할 수 있다.
전체 시간이 `Δt`의 정수 배가 아니면 입력을 거부하는 방안을 권장한다. 짧은 마지막 증분은
고정 `Δt`와 factorization 재사용 의미를 바꾸므로 자동 도입하지 않는다. 시간 비교의 수치
판정과 종료 규칙은 사용자 검토와 I/O 계약에서 확정한다.
고정 geometry·물성·DOF 집합과 동일한 `δ_n`이면 `A_ff` factorization을 재사용한다.
마지막 실제 `δ_n`이 바뀌면 다시 factorize한다. 근접한 계수라는 이유로 재사용하지 않는다.
**초기 불일치 제안:** 입력 초기 온도를 5점에 복제한 전체 `T^0_raw`를 보존하며,
**승인된 초기 상태:** 절점/nset 초기 온도를 5점에 복제한 전체 `T^0_raw`를 보존하며,
시간 적분의 `T^0=T^0_raw`로 사용한다. `t=0`에는 지정 온도로 덮어쓰지 않는다.
지정 온도는 첫 완료 증분 `t_1`부터 정확히 적용한다. 입력 초기값과 `g`가 다른 경우
그 차이는 첫 증분에 들어가므로 `C_fc(T_c^0−g)/Δt`와 저장에너지가 모두 살아 있다.
그 차이는 첫 증분에 들어가므로 `C_fc(T_c^0−g)/δ_0`와 저장에너지가 모두 살아 있다.
이는 초기 불연속을 첫 유한 증분에서 처리하는 이산 정책으로, 실제 순간 impulse의
해상이나 `t=0+` 경계 만족을 주장하지 않는다. `t=0` 반력 열률은 시간 차분이 없어
미정의다. 0으로 만들어 ‘계산된 반력’으로 저장하지 않는다.
대안인 `t=0` 즉시 projection은 별도 에너지 jump `1ᵀC(T^0_projected−T^0_raw)`와
그 보존/보고 정책이 필요하다. 이 선택은 결과를 바꾸므로 사용자 결정 후 하나만 구현한다.
생략된 초기 온도의 값은 이 문서에서 결정하지 않는다. Reference는 이 default에 민감하다.
미정의다. `t=0`에는 temperature만 저장하며 flux/reaction dataset을 만들지 않는다.
생략된 초기 온도는 0이다. Model-level BC를 먼저 적용하고 step-level BC가 같은 target을
override한다. 같은 stage 안에서 상충하는 duplicate 값은 오류다. 이 단계별 해석을
raw 초기 상태의 projection으로 바꾸지 않는다.
각 증분은 solve → full 온도 복구 → residual/flux/유한성·수지 확인 → candidate commit
순서로 처리한다. 실패 candidate가 이전 완료 state나 최종 HDF5를 대체하지 않는다.
### 7.1 재현 가능한 종료 시점과 증분 수
입력 binary64 `d=nominal Δt`, `τ=T_end`는 모두 유한한 양수다. `INC`는 양의 uint32,
기본값 100이며 허용 증분 수의 상한이다. `d≥τ`이면 한 증분, `δ_0=τ`, 저장 시각 `τ`다.
그 외에는 **반복 덧셈 없이** 정수 `j`의 곱으로 정상 시각을 만든다.
정수비 근방의 동등성은 **입력 binary64로 반올림될 수 있는 실제 시간 구간**으로 판정한다.
입력 중심값을 오차 없는 물리량으로 간주하지 않는다. `fl`은 nearest-even이며 먼저
`q=fl(τ/d)`를 구한다. `INC+1`은 uint64 범위에서 계산한다.
`q`가 overflow하거나 `q>INC+1`이면 증분 수 초과다.
나머지는 `m=nearest_integer(q)`를 후보로 쓴다. Tie는 nearest-even으로 고정하며
`m≤2³²`는 정확히 표현된다. `INC`에 맞추려고 다른 정수 후보로 바꾸지 않는다.
원래 입력 `x>0`의 predecessor/successor를 `x−,x+`라 하면 rounding cell은
\[
I(x)=\left\langle (x^-+x)/2,\ (x+x^+)/2\right\rangle,
\qquad \text{endpoint-close}\iff m\ge1\ \land\ I(\tau)\cap mI(d)\ne\varnothing.
\]
꺾쇠의 두 끝점은 `x`의 binary64 significand 최하위 bit가 0(even)이면 포함하고,
1(odd)이면 제외한다. Subnormal도 같은 규칙이며 최소 양수의 predecessor는 0이다.
최대 finite의 successor는 부동소수점 `+∞` 대신 **수학적인 다음 lattice 값 `2^1024`**로
둔다. 최대 finite significand는 odd이므로 overflow midpoint를 cell에 포함하지 않는다.
Midpoint를 double로 반올림하지 않고 아래 정수·지수 표현으로 정확히 구성한다.
입력 변환은 correctly-rounded nearest-even을 전제로 하며 임의 입력 오차는 더하지 않는다.
`I(τ)=(a_L,a_U)`, `mI(d)=(b_L,b_U)`의 끝점 포함 flag를 보존하여
`L=max(a_L,b_L)`, `U=min(a_U,b_U)`를 구한다. `L<U`면 close, `L>U`면 nonclose다.
`L=U`면 그 접점이 **두 cell 모두에 포함될 때만** close다. 나눗셈·곱셈 오차와 입력 cell
폭을 중복 가산하는 epsilon은 없다. `q`는 정수 후보만 만들고, 겹침 판정에는 반올림된
`fl(md)`를 사용하지 않는다. 같은 endpoint를 허용하지 않는 cell은 잔여시간이 작아도 구별한다.
정확한 predicate는 임의정밀도 dependency 없이 다음 고정 폭 연산으로 구현할 수 있다.
- `e=frexp(τ).exponent`, `s=2^(−e)`로 두어 `sτ∈[1/2,1)`로 만든다. **원래 입력의
cell을 먼저 구성한 뒤** `s`를 곱한다. Normal화된 값의 prev/next로 cell을 다시 만들면
원래 subnormal 입력의 표현 범위를 줄이므로 금지한다. `s` 자체를 double로 만들 필요는 없다.
- `d<τ`, 유한 `q≤INC+1≤2³²`이면 정확한 `τ/d`는 `2³²+2^(−21)` 이하이므로
`sd>2^(−34)`이다. Scaled 중심값과 원래 cell의 midpoint는 모두 lattice `2^(−88)`의
정수배다. Normal 입력의 최소 반간격은 이 lattice보다 작지 않으며, subnormal은 더
성긴 원래 간격을 유지한다. Binade 경계의 비대칭 간격도 따로 보존한다.
- `X=τs2^88`, `Y=ds2^88`와 네 cell endpoint를 원래 binary64 significand/exponent에서
shift·덧셈으로 정확히 구성한다. Scaled cell 끝점은 2보다 작으므로 각 정수는 89 bit 이하다.
`m` 및 floor 교정용 `floor(q)+1`은 최대 `2³²+1`이므로 곱은 **122 bit 이하**다.
- 두 `uint64` 상당의 high/low word로 저장·비교한다. High word를 먼저 비교하고, 같으면
low word를 비교한다. 덧셈/뺄셈의 carry/borrow와 shift의 word 경계를 명시적으로 처리한다.
곱은 네 base-`2³²` digit와 정수 승수의 최대 두 digit로 계산할 수 있다. 각
`a_i*b_j+result_digit+carry≤2^64−1`이므로 중간값도 uint64에 들어가며 다음 digit으로
carry를 넘긴다. 최상위 128 bit 밖 carry가 없음을 검사한다. `__int128`, extended float,
compiler 연산 재결합 또는 arbitrary-precision library에 의존하지 않는다.
1. Endpoint-close이면 `N=m`. 아니면 중심값의 정확한 `N=ceil(τ/d)`를 쓴다.
`f=floor(q)`에서 정수 비교 `fY≤X<(f+1)Y`를 확인하고 필요하면 `f`를 한 번 ±1
교정한다. `q≤2³²`의 division 오차는 1보다 작으므로 한 번이면 충분하다.
`X=fY`이면 `N=f`, 아니면 `N=f+1`이다. Division 반올림을 증분 수에 그대로 전달하지 않는다.
`N`은 `INC`와 비교하기 전 uint32로 축소하지 않는다.
2. `N>INC`이면 solve 전에 실패한다. 정상 저장 시각은 `t_j=fl(j d)` (`j=1,…,N−1`),
마지막은 **입력 그대로** `t_N=τ`다. 모든 시각이 유한하고 엄격히 증가하며
`t_{N−1}<τ`인지 preflight에서 확인한다. `t_0=0`이다.
3. 정상 증분과 endpoint-close 마지막 증분의 연산 계수는 `δ=d`다. 실제로 짧아지는
마지막 증분만 `δ=fl(τ-(N−1)d)`로 만든다. 이 remainder는
`fma(-(N−1),d,τ)`처럼 곱-차를 한 번에 반올림한다. `δ>0`·유한성을 검사한다.
저장된 `fl(jd)`끼리 다시 빼서 정상 `δ`를 만들지 않는다. 개념적인 격자 `j d`와
저장 표현의 차이는 각 곱의 rounding cell로 제한된다. Endpoint-close 조정은 입력 cell
overlap으로 제한된다. 원래 단위의 `|τ−md|`는 해당 방향의 `τ` 반간격과 `m`배 `d`
반간격 합 이내이고, 계산된 FMA remainder에는 한 번의 최종 rounding만 남는다.
4. 해당 증분의 `A`, RHS, full residual, reaction, `ΔE` 수지, BE 소산은 모두 **같은
실제 `δ`**를 사용한다. Endpoint-close는 동일 계수를 재사용하며 실제 shortened
마지막 `δ`가 `d`와 다를 때 refactorize한다.
예를 들어 `(τ,d)=(10,3)`은 `N=4`, `δ=(3,3,3,1)`, `(τ,d)=(2,3)`은 `N=1,δ=2`다.
`(0.9,0.3)` 같은 binary64 정수비 근방은 불필요한 미소 증분 없이 3회가 되며
세 번 모두 `δ=0.3`, 마지막 저장 시각은 `0.9`다. 표현오차만으로 재분해하지 않는다.
`nextafter`로 endpoint-close 경계를 양쪽에서 이동한 사례, 아주 작고 큰 시간 척도,
`INC=N`/`N−1`, uint32 overflow 방지를 독립 검증한다. 이 격자는 adaptive time stepping이 아니다.
### 7.2 입력 cell 수정의 실제 검산 기록
2026-09-29 Python 3.11 `Fraction`을 독립 exact-dyadic oracle로 사용한 아래 scratch를
실행했다(exit 0). 이는 문서 수식 검산이며 solver build/test나 수치 correctness 승인이 아니다.
기존 `Bhat=nextafter(fl(d rad(q)+rad(fl(md))),+∞)`는 입력 cell을 누락하여 첫 사례를
nonclose로 오판했다. Scratch의 임의정밀도는 검산용이며 production dependency가 아니다.
```text
old-failure: m=31, exactD=8.326672684688674e-16 > Bhat=7.507374033534679e-16
corrected: close=True, N=31; INC=31 accept, INC=30 reject
tau=nextafter(.9,+inf), d=.3: close=False, N=4, delta=1.6653345369377348e-16
tau=.9, d=.3: close=True, N=3
tau=20*minsub, d=6*minsub: even/even midpoint included, close=True, N=3
tau=17*minsub, d=5*minsub: odd/odd midpoint excluded, close=False, N=4
power-of-two scales [-900,-400,0,400,900]: 5/5
maxfinite, original-subnormal, binade cases: 5/5
seed41 input-cell versus fixed-lattice predicate: 10000/10000
uint32 INC boundary: N=4294967295 accept; N=4294967296 reject
```
<details>
<summary>실행한 재현 scratch (Python 표준 라이브러리)</summary>
```python
import math, random, struct
from fractions import Fraction as F
def cell(x):
a=F(x); p=F(math.nextafter(x,0.0))
nx=math.nextafter(x,math.inf)
n=F(nx) if math.isfinite(nx) else a+(a-p)
even=(struct.unpack('>Q',struct.pack('>d',x))[0]&1)==0
return (a+p)/2,(a+n)/2,even
def overlap(a,b):
lo=max(a[0],b[0]); hi=min(a[1],b[1])
return lo<hi or (lo==hi and a[2] and b[2])
def predicate(t,d):
m=round(t/d); a=cell(t); b=cell(d)
return overlap(a,(m*b[0],m*b[1],b[2]))
def fixed(t,d):
e=math.frexp(t)[1]; s=F(2)**(-e+88)
a=cell(t); b=cell(d); m=round(t/d)
vals=[a[0]*s,a[1]*s,b[0]*s,b[1]*s,F(t)*s,F(d)*s]
assert all(v.denominator==1 for v in vals)
al,ah,bl,bh,x,y=map(int,vals)
assert max(al,ah,m*bl,m*bh,x,y).bit_length()<=122
return overlap((al,ah,a[2]),(m*bl,m*bh,b[2]))
def count(t,d):
if d>=t: return 1
if predicate(t,d): return round(t/d)
z=F(t)/F(d); return -(-z.numerator//z.denominator)
def count_fixed(t,d):
if d>=t: return 1
if fixed(t,d): return round(t/d)
e=math.frexp(t)[1]; s=F(2)**(-e+88)
x,y=int(F(t)*s),int(F(d)*s); f=math.floor(t/d)
if f*y>x: f-=1
elif (f+1)*y<=x: f+=1
assert f*y<=x<(f+1)*y
return f if x==f*y else f+1
def rad(x):
return max(F(x)-F(math.nextafter(x,-math.inf)),
F(math.nextafter(x,math.inf))-F(x))/2
t,d=5.351455173733616,.17262758624947144
m=round(t/d); q=t/d; p=m*d
old=math.nextafter(float(F(d)*rad(q)+rad(p)),math.inf)
distance=abs(F(t)-m*F(d))
assert distance>F(old) and predicate(t,d) and fixed(t,d)
print('old-failure:',m,float(distance),old,distance<=F(old))
print('corrected:',predicate(t,d),'N=',count(t,d),'INC31=',count(t,d)<=31,'INC30=',count(t,d)<=30)
t2=math.nextafter(.9,math.inf)
print('real-remainder:',predicate(t2,.3),'N=',count(t2,.3),'delta=',float(F(t2)-3*F(.3)))
print('decimal:',predicate(.9,.3),count(.9,.3))
tiny=math.ulp(0.0)
for a,b in [(20,6),(17,5)]:
print('midpoint:',a,b,cell(a*tiny)[2],cell(b*tiny)[2],predicate(a*tiny,b*tiny),count(a*tiny,b*tiny))
for e in [-900,-400,0,400,900]:
tt,dd=math.ldexp(t,e),math.ldexp(d,e)
assert predicate(tt,dd) and fixed(tt,dd) and count(tt,dd)==31
print('power-of-two-scales: 5/5')
maxf=float.fromhex('0x1.fffffffffffffp+1023')
for tt,dd in [(maxf,maxf/31),(20*tiny,6*tiny),(17*tiny,5*tiny),
(1.,math.nextafter(.25,0.)),(1.,math.nextafter(.25,math.inf))]:
assert predicate(tt,dd)==fixed(tt,dd)
assert count(tt,dd)==count_fixed(tt,dd)
print('extreme/binade cases: 5/5')
rng=random.Random(41)
for _ in range(10000):
dd=math.ldexp(rng.uniform(.5,1),rng.randint(-900,900))
mm=rng.randint(2,2**32)
tt=mm*dd
for _ in range(rng.randint(0,2)):
tt=math.nextafter(tt,rng.choice([0.,math.inf]))
assert predicate(tt,dd)==fixed(tt,dd)
print('seed41 cell-vs-fixed: 10000/10000')
assert count(float(2**32-1),1.)==2**32-1
assert count(float(2**32),1.)==2**32
print('uint32 INC limit: N=4294967295 accept; N=4294967296 reject')
```
</details>
## 8. 반력 열률, 에너지와 안정성
반력은 하중 제거 전 full residual에서 구한다.
\[
r^{n+1}=C(T^{n+1}-T^n)/\Delta t+KT^{n+1}-F^{n+1},\quad
r^{n+1}=C(T^{n+1}-T^n)/\delta_n+KT^{n+1}-Q,\quad
R_c=r_c,\quad r_f\approx0,\qquad
U(T)=\mathbf1^TC(T-T_{\rm ref}\mathbf1).
\]
@@ -247,16 +498,19 @@ U(T)=\mathbf1^TC(T-T_{\rm ref}\mathbf1).
`K1=0`과 대칭성에 의해 증분의 전역 수지는
\[
U^{n+1}-U^n=\Delta t\left(\mathbf1^TF^{n+1}+\mathbf1_c^TR_c^{n+1}
\Delta E=U^{n+1}-U^n=\delta_n\left(\mathbf1^TQ+\mathbf1_c^TR_c^{n+1}
+\mathbf1_f^Tr_f^{n+1}\right)
\]
이다. Free residual을 버리지 않고 수지 defect와 함께 확인한다. 정상상태에서는 총 입력과
이다. Physical balance defect `ΔE−δ_n(1ᵀQ+1_cᵀR_c)`와 `r_f`의 작음은 별도로
확인한다. 위 항등식에 `r_f`를 넣어 맞는다는 사실만으로 평형을 승인하지 않는다.
Full free residual과 energy는 내부 검증량이며 추가 공개 결과 quantity로 만들지 않는다.
정상상태에서는 총 입력과
반력 합이 0이다. 단열·무입력·무구속에서는 `U`가 보존되고 다음 이산 소산식이 성립한다.
\[
\tfrac12\big(\|T^{n+1}\|_C^2-\|T^n\|_C^2+\|T^{n+1}-T^n\|_C^2\big)
+\Delta t\|T^{n+1}\|_K^2=0.
+\delta_n\|T^{n+1}\|_K^2=0.
\]
이 식의 이차량은 안정성 척도이며 실제 열에너지 `U`와 다르다.
@@ -281,7 +535,7 @@ Global 성분은 Cartesian `qx,qy,qz`, local은 `q1,q2,q3`이며 모두 열유
`q3=−k(T_top−T_bottom)/h`; top이 뜨거우면 `q3<0`이다.
Bottom 외향 normal은 `−n`, top은 `+n`이므로 외향 성분은 각각 `−q3`, `+q3`다.
결정적인 필수 recovery 위치 제안은 면내 4점
결정적인 필수 recovery 위치는 면내 4점
`(-g,-g),(g,-g),(g,g),(-g,g)`, `g=1/√3`와 두께 trace
`(-1,+),(-1/2,interior),(0,−),(0,+),(1/2,interior),(1,−)`의 곱인 **24점/요소**다.
두께 trace의 부호는 어느 구간 내부에서 접근하는지를 뜻한다. Integration 점과 recovery
@@ -289,29 +543,77 @@ Bottom 외향 normal은 `−n`, top은 `+n`이므로 외향 성분은 각각 `
이 둘을 평균하여 원자료를 잃지 않는다. 연속해의 재료 계면을 뜻하는 것도 아니다.
Row는 최소 `(instance, element label, in-plane point, thickness station, trace side, time)`와
component/좌표계로 구별해야 한다. Global/local 중 영속 저장 형식은 I/O owner가 결정한다.
component/좌표계로 구별한다. 공식 flux는 global Cartesian `q`이며 직교 frame
`R=[e1,e2,n]`을 보존하여 후처리에서 `q_local=Rᵀq`를 구한다.
요소 사이 edge의 flux도 요소별 한쪽 값을 유지한다. 서로 다른 `k`, normal, thickness에서
벡터를 무조건 같게 만들지 않는다. 접합부의 보존량은 약형으로 조립되는 열률이다.
`−kBT` raw flux의 경계 적분은 이산 반력과 일반적으로 정확히 같지 않으며 반력을 대신하지 않는다.
절점 외삽, nodal averaging, flux smoothing은 이 제안의 필수 recovery가 아니다.
절점 외삽, nodal averaging, flux smoothing으로 원래 trace를 대체하지 않는다.
## 10. 검증 함의와 미결정 사항
## 10. 검증 함의와 독립 검토 인계
다음은 실행 결과가 아닌 검증 설계다. 수치 기준은 독립 reviewer가 정하며
외부 비교는 [공통 TOLERANCE](../TOLERANCE.md)를 그대로 따른다.
| 검증 | 관찰할 내용 |
| --- | --- |
| Shape/행렬 | Nodal interpolation, partition of unity, `ΣB=0`, explicit Z와 직접 적분 일치, `C` SPD, `K` PSD/rank 19, symmetry |
| Shape/행렬 | Nodal interpolation, partition of unity, `ΣB=0`, exact M/S/Z와 54점 직접 적분 일치, `C` SPD, `K` PSD/rank 19, symmetry, C·면내 K의 h scaling 및 두께 K의 1/h scaling |
| 평판 patch·회전 | Constant/면내 선형/두께 선형 온도에서 gradient·flux, 회전 시 T 불변과 q 회전, 법선 반전+station 재매핑 |
| 열입력 | Edge 총률 `q_in hL`·5점 weight, bottom/top 집중 입력, prescribed DOF와 겹친 입력의 반력 |
| 열입력 | Bottom/top 집중 입력 P 그대로, prescribed DOF와 겹친 P 보존과 반력 |
| 과도·초기 | 단열 균일온도 보존, 초기 불일치 첫 RHS, 공개 입력의 no-constrained, 에너지 수지·BE 소산; all-constrained는 generic partition fixture만 사용 |
| 면내 해석해 | 단열 직사각판 `T=T0+A cos(πx/L)exp[-α(π/L)²t]`, `α=k/(ρc)`; 면내 초기값은 5점 동일 |
| 두께 해석해 | 균일 초기온도 slab의 두 표면 온도 step에 대한 1D Fourier-series 해와 비교; 초기 jump와 첫 step 분리 |
| 두께 해석해 | 균일 초기온도 slab의 bottom 온도 step/top 단열에 대한 아래 series; 초기 jump와 첫 step 분리 |
| 재료·접속 | 두 물성/두께 strip의 직렬 전도는 아래 kernel/assembly oracle 범위로 구분; folded 등온/열률 수지와 junction 이상화의 별도 3D 한계 |
| 수렴·특이점 | `Δt` 1차 수렴, 면내 mesh/적분 차수 수렴, 집중 입력 근방 flux의 mesh 의존성과 곡률 오차 |
| 수렴·특이점 | 독립 semi-discrete eigenpair의 시간 1차 수렴, smooth affine mesh의 temperature L2 2차·gradient L2 1차, 집중 입력 근방 flux의 mesh 의존성과 곡률 오차 |
| Geometry·격자 | Whole-patch certificate, bow-tie/퇴화/uncertain metric, T-junction·비다양체·vertex-only fan·공유 edge 방향 오류, disconnected mesh, warped 3/5/7 적분 진단, 마지막 actual δ·endpoint-close·INC preflight |
| Reference | 기존 Step-1, `t=37`, NT11 bottom만 source identity로 비교; initial default·첫 증분·C 처리 차이 조사 |
단위 정사각형의 독립 exact surface oracle는 다음과 같다.
\[
M=\frac1{36}\begin{bmatrix}4&2&1&2\\2&4&2&1\\1&2&4&2\\2&1&2&4\end{bmatrix},\qquad
S=\frac1{6}\begin{bmatrix}4&-1&-2&-1\\-1&4&-1&-2\\-2&-1&4&-1\\-1&-2&-1&4\end{bmatrix}.
\]
`Z0·1=(1,4,2,4,1)/6`, `Z1·1=0`; 전체 thermal component마다 `K` null mode가
하나이며 모든 활성 DOF가 양의 체적 요소에 속하면 no-BC/disconnected 모델의
`C/δ+K`도 SPD다. 전체 K rank를 무조건 `N_dof−1`로 검사하지 않는다.
한쪽 trace fixture `(T1,…,T5)=(0,0,0,1,2)`는 `q3(0−)=0`, `q3(0+)=−4k/h`다.
최초 jump 독립 algebra fixture는 `C=[[2,1],[1,2]]`, `K=[[1,−1],[−1,1]]`,
`δ=1/2`, `Q=0`, raw `T0=(0,0)`, `g2=1`이다. 첫 두 증분의 free 온도는
`−1/5,1/25`, constrained reaction은 `24/5,36/25`이며 각각 `ΔE=δR_c`다.
이 기대값은 첫 `C_fc T_c_old/δ` 항을 검출한다. Generic all-constrained fixture와
함께 사용하되 공개 내부 station BC나 음수 clamp를 도입하지 않는다.
시간 수렴은 고정 mesh의 독립 generalized eigenpair `Kv=λCv`로
`T(t)=exp(−λt)v`와 BE의 `Π_j(1+λδ_j)⁻¹v`를 비교한다. 이 내부 초기 state를
공개 두께 비균일 초기 입력으로 확대하지 않는다. Temperature/gradient 공간 오차는
각각 `||T_h−T||_L2(V*)`, `||∇_*T_h−∇_*T||_L2(V*)`로 측정하고 시간 오차를 충분히
분리한다. 관찰 차수 `log(e_h/e_(h/2))/log(2)`의 수치 허용 구간은 reviewer가 고정한다.
Slab에서는 `x=z+h/2∈[0,h]`, bottom `Ts`, top 단열, 초기 `Ti`에 대해 `t>0`에서
\[
\theta=\frac{T-T_s}{T_i-T_s}
=\sum_{m=0}^{\infty}\frac4{(2m+1)\pi}
\sin\!\left(\frac{(2m+1)\pi x}{2h}\right)
\exp\!\left[-\alpha\left(\frac{(2m+1)\pi}{2h}\right)^2t\right],
\quad \alpha=\frac{k}{\rho c}.
\]
항별 미분으로 q를 구한다. 첫 생략 index가 `M`, `β=απ²t/(4h²)>0`이면 연속 항의
지수비는 최대 `r=exp[-8β(M+1)]<1`이다. 따라서 모든 `x`에서
\[
|\theta-\theta_M|\le\frac{4e^{-\beta(2M+1)^2}}{\pi(2M+1)(1-r)},\qquad
|q-q_M|\le\frac{2k|T_i-T_s|}{h}\frac{e^{-\beta(2M+1)^2}}{1-r}.
\]
`1-r`는 cancellation을 피하도록 `-expm1(-8β(M+1))`로 평가한다. 검증 시간·절단
`M`은 이 bound가 해당 oracle error budget보다 작도록 reviewer가 선택한다. 수치 적분한
slab energy는 temperature tail bound에 `ρc h·area·|Ti−Ts|`를 곱해 제한할 수 있다.
두 물성/두께 strip의 `R_th=Σ_e L_e/(k_e h_e w_e)` 정확해 검증은 양 끝 단면의 모든
station에 온도를 지정하는 kernel/assembly oracle이다. 내부 station 구속이 필요한 이
fixture는 공개 CLI acceptance 사례로 쓰지 않는다. 공개 입력으로 bottom/top만 구속한
@@ -323,10 +625,19 @@ Slab 해석해에 대한 두께 오차와 초기 `sqrt(αt)` 열침투 깊이를
Reference CSV 한 시점은 다른 4개 온도점, flux, reaction 또는 전체 이력의 evidence가 아니다.
Consistent C와 본 적분·초기화 선택이 reference 정확도에 미치는 영향은 아직 측정되지 않았다.
미결정 owner는 메인/사용자(물리·초기 정책), I/O owner(공개 의미), numerical reviewer(수치 허용 범위)다.
사용자 검토에는 thin-shell product geometry와 접속 이상화, C/BE, 최초 불일치·종료 시점,
24개 flux trace 제안을 포함한다. Geometry 판정 기준과 기본 적분의 지원 범위, 생략 초기값,
edge 외부 단위, t=0 결과의 validity 표현은 구현 전에 확정해야 한다.
현재 정식화 담당자의 확인은 문헌·요구사항 대조와 문서 정적 검토뿐이다. Solver build/test,
reference 수치 비교와 physics 검토는 미실행이다. 독립 numerical review의 결과는 별도
보고서를 참조하며 본 문서는 스스로 수치 타당성이나 결과 correctness를 승인하지 않는다.
Residual backward error는 `||r_f||/(||A_ff|| ||T_f||+||b_f||)`로, physical balance
defect는 `|ΔE−δ(1ᵀQ+1_cᵀR_c)|`를 저장·입력·반력 항 절댓값 합으로 정규화한다.
분모 0은 exact-zero case로 별도 처리한다. 행렬 오차도 해당 oracle norm과 연산 bound를
사용하며 단위 의존 `max(1,scale)`을 넣지 않는다. Geometry certificate·시간 predicate의
연산 bound와 convergence acceptance는 서로 다른 것이다.
Reviewer는 production kernel 착수 전에 roundoff/residual/energy 판정의 구체 bound,
mesh·시간열·series 절단·수렴차수 구간·warped quadrature 진단을 독립적으로 고정한다.
본 문서는 임의 tolerance 상수를 만들거나 공통 외부 5%를 내부 검증에 전용하지 않는다.
필수 NT11 비교는 `Step-1/t37/Part-1-1/source node/DOF11`의 729개 원본 row mapping과
[공통 TOLERANCE](../TOLERANCE.md)를 유지한다. 실패는 조사·상위 보고 대상이며
reference 보정·음수 clamp·방법 변경·tolerance 완화의 근거가 아니다.
물리·방법의 과거 미결정은 이번 승인 계약으로 대체되었다. 남은 독립 작업은 위 수치
acceptance 고정과 실제 구현 검증이다. Solver build/test, reference 수치 비교와 physics
검토는 본 정식화 작업에서 미실행이며, 스스로 수치 타당성이나 결과 correctness를 승인하지 않는다.
@@ -0,0 +1,504 @@
# 독립 numerical oracle scratch
이 파일은 numerical-review.md의 사전 계산을 재현하는 문서 부록이다. Solver 또는
production test가 아니며 외부 package·reference·출력 파일을 사용하지 않는다.
Python 표준 라이브러리만 사용한다. Fraction은 oracle에서 입력 double의 exact 산술을
확인하는 용도이며 production arbitrary-precision dependency를 요구하지 않는다.
T1 root에서 다음 명령은 각 Python block을 순서대로 메모리에서 실행한다.
```powershell
@'
from pathlib import Path
import re
p = Path('docs/shell-transient-heat-transfer/numerical-oracle-scratch.md')
blocks = re.findall(r'```python\n(.*?)```', p.read_text(encoding='utf-8-sig'), re.S)
assert len(blocks) == 5
scope = {}
for block in blocks:
exec(compile(block, str(p), 'exec'), scope)
print('ALL SCRATCH BLOCKS COMPLETED')
'@ | python -
```
실행 exit는 0이어야 한다. Assertion 범위는 exact algebra와 독립 oracle 자체에 한정된다.
수렴/이산 오차의 출력 숫자 자체를 production tolerance로 사용하지 않고 본 보고서의
roundoff·conditioning·오차 전파 기준을 적용한다. NumPy 부재 probe exit1은 실행 실패
이력으로 보고서에 남겼으며 아래 코드에는 NumPy를 쓰지 않는다.
## A. Slab spectral·quadrature·공간/시간 oracle
```python
import math,json
def zeros(n,m=None): return [[0.0]*(n if m is None else m) for _ in range(n)]
def tr(A):return list(map(list,zip(*A)))
def mm(A,B):return [[sum(a*b for a,b in zip(row,col)) for col in zip(*B)]for row in A]
def mv(A,x):return [sum(a*b for a,b in zip(row,x)) for row in A]
def solve(A,b):
a=[list(r)+[v]for r,v in zip(A,b)];n=len(b)
for k in range(n):
p=max(range(k,n),key=lambda i:abs(a[i][k]));a[k],a[p]=a[p],a[k]
d=a[k][k]
for j in range(k,n+1):a[k][j]/=d
for i in range(n):
if i!=k:
d=a[i][k]
for j in range(k,n+1):a[i][j]-=d*a[k][j]
return [r[-1]for r in a]
def chol(A):
n=len(A);L=zeros(n)
for i in range(n):
for j in range(i+1):
v=A[i][j]-sum(L[i][k]*L[j][k]for k in range(j))
L[i][j]=math.sqrt(v)if i==j else v/L[j][j]
return L
def eig(A):
a=[r[:]for r in A];n=len(a);v=[[float(i==j)for j in range(n)]for i in range(n)]
for z in range(100*n*n):
p,q=max(((i,j)for i in range(n)for j in range(i+1,n)),key=lambda ij:abs(a[ij[0]][ij[1]]))
if abs(a[p][q])<1e-15*sum(abs(a[i][i])for i in range(n)):break
phi=.5*math.atan2(2*a[p][q],a[q][q]-a[p][p]);c=math.cos(phi);s=math.sin(phi)
J=[[float(i==j)for j in range(n)]for i in range(n)]
J[p][p]=J[q][q]=c;J[p][q]=s;J[q][p]=-s
a=mm(mm(tr(J),a),J);v=mm(v,J)
ix=sorted(range(n),key=lambda i:a[i][i])
return [a[i][i]for i in ix],[[v[j][i]for i in ix]for j in range(n)]
def slabmat(ne):
n=2*ne+1;C=zeros(n);K=zeros(n);l=1/ne
for e in range(ne):
for i in range(3):
for j in range(3):
C[2*e+i][2*e+j]+=l/30*[[4,2,-1],[2,16,2],[-1,2,4]][i][j]
K[2*e+i][2*e+j]+=[[7,-8,1],[-8,16,-8],[1,-8,7]][i][j]/(3*l)
return C,K
C,K=slabmat(2);Cf=[r[1:]for r in C[1:]];Kf=[r[1:]for r in K[1:]]
L=chol(Cf);Li=tr([solve(L,[float(i==j)for i in range(4)])for j in range(4)])
lam,V=eig(mm(mm(Li,Kf),tr(Li)));W=mm(tr(Li),V)
yp=[-1-x for x in solve(Cf,[r[0]for r in C[1:]])]
coef=mv(tr(W),mv(Cf,yp))
def spectral(t,n=None):
amp=[a*(math.exp(-l*t)if n is None else math.exp(-n*math.log1p(l*t/n)))for a,l in zip(coef,lam)]
T=[1]+[1+x for x in mv(W,amp)]
if n is None:dT=[0]+mv(W,[-l*x for l,x in zip(lam,amp)])
else:
oldamp=[x*(1+l*t/n)for l,x in zip(lam,amp)]
old=[1]+[1+x for x in mv(W,oldamp)]
dT=[(a-b)/(t/n)for a,b in zip(T,old)]
E=sum(mv(C,T));R=mv(C,dT)[0]+mv(K,T)[0]
q=[]
for e,ss in [(0,-1),(0,0),(0,1),(1,-1),(1,0),(1,1)]:
q.append(-4*sum(a*b for a,b in zip([ss-.5,-2*ss,ss+.5],T[2*e:2*e+3])))
return T,q,E,R
def exact(t):
M=80
T=[1-sum(4/((2*m+1)*math.pi)*math.sin((2*m+1)*math.pi*x/2)*math.exp(-((2*m+1)*math.pi/2)**2*t)for m in range(M))for x in [0,.25,.5,.75,1]]
q=[sum(2*math.cos((2*m+1)*math.pi*x/2)*math.exp(-((2*m+1)*math.pi/2)**2*t)for m in range(M))for x in [0,.25,.5,.5,.75,1]]
E=1-sum(8/((2*m+1)*math.pi)**2*math.exp(-((2*m+1)*math.pi/2)**2*t)for m in range(M))
return T,q,E,q[0]
def errors(a,b):return [max(abs(x-y)for x,y in zip(a[0],b[0])),max(abs(x-y)for x,y in zip(a[1],b[1])),abs(a[2]-b[2]),abs(a[3]-b[3])]
out={"lambda":lam,"projected_initial_free":[1+x for x in yp],"slab":[]}
for t in [.001,.01,.05,.1,.5,1]:
s=spectral(t);x=exact(t)
out["slab"].append({"t":t,"penetration":math.sqrt(t),"spatial_error":errors(s,x),"BE4096_temporal":errors(spectral(t,4096),s),"BE8192_temporal":errors(spectral(t,8192),s),"spectral":s})
print(json.dumps(out,indent=2))
def gauss(n):
x=[];w=[]
for i in range(1,n+1):
z=math.cos(math.pi*(i-.25)/(n+.5))
for _ in range(30):
p0,p1=1.,z
for j in range(2,n+1):p0,p1=p1,((2*j-1)*z*p1-(j-1)*p0)/j
dp=n*(z*p1-p0)/(z*z-1);nz=z-p1/dp
if abs(nz-z)<2e-16: z=nz;break
z=nz
x.append(z);w.append(2/((1-z*z)*dp*dp))
return list(zip(x,w))
def dot(x,y):return sum(a*b for a,b in zip(x,y))
def cross(a,b):return [a[1]*b[2]-a[2]*b[1],a[2]*b[0]-a[0]*b[2],a[0]*b[1]-a[1]*b[0]]
def geometry(X,xi,eta):
signs=[(-1,-1),(1,-1),(1,1),(-1,1)]
N=[(1+s*xi)*(1+t*eta)/4 for s,t in signs]
D=[[s*(1+t*eta)/4,t*(1+s*xi)/4]for s,t in signs]
a=[sum(D[i][0]*X[i][j]for i in range(4))for j in range(3)]
b=[sum(D[i][1]*X[i][j]for i in range(4))for j in range(3)]
W=cross(a,b);J=math.sqrt(dot(W,W));norm=[z/J for z in W]
aa,ab,bb=dot(a,a),dot(a,b),dot(b,b);det=J*J
ra=[(bb*x-ab*y)/det for x,y in zip(a,b)];rb=[(aa*y-ab*x)/det for x,y in zip(a,b)]
G=[[D[i][0]*ra[j]+D[i][1]*rb[j]for j in range(3)]for i in range(4)]
return N,G,J,norm
def surface(X,n):
M=zeros(4);S=zeros(4)
for xi,wx in gauss(n):
for eta,wy in gauss(n):
N,G,J,no=geometry(X,xi,eta)
for i in range(4):
for j in range(4):
M[i][j]+=N[i]*N[j]*J*wx*wy;S[i][j]+=dot(G[i],G[j])*J*wx*wy
return M,S
Z0,Z1=slabmat(2)
# h=1: Z0 physical mass already; Z1 physical normal stiffness already.
def mats(X,n):
M,S=surface(X,n);C=zeros(20);K=zeros(20)
for i in range(4):
for j in range(4):
for a in range(5):
for b in range(5):
C[5*i+a][5*j+b]=M[i][j]*Z0[a][b]
K[5*i+a][5*j+b]=S[i][j]*Z0[a][b]+M[i][j]*Z1[a][b]
return C,K
def fn(A):return math.sqrt(sum(x*x for r in A for x in r))
def rel(A,B):return fn([[a-b for a,b in zip(r,s)]for r,s in zip(A,B)])/fn(B)
def response(X,C,K):
A=[[c/.1+k for c,k in zip(cr,kr)]for cr,kr in zip(C,K)]
T=solve(A,[float(i==0)for i in range(20)])
flux=[]
for xi,eta in [(-1/math.sqrt(3),-1/math.sqrt(3)),(1/math.sqrt(3),-1/math.sqrt(3)),(1/math.sqrt(3),1/math.sqrt(3)),(-1/math.sqrt(3),1/math.sqrt(3))]:
N,G,J,no=geometry(X,xi,eta)
for off,s in [(0,-1),(0,0),(0,1),(2,-1),(2,0),(2,1)]:
Ls=[s*(s-1)/2,1-s*s,s*(s+1)/2];Ds=[4*(s-.5),-8*s,4*(s+.5)]
flux.append([-sum((G[i][j]*Ls[a]+no[j]*N[i]*Ds[a])*T[5*i+off+a]for i in range(4)for a in range(3))for j in range(3)])
return T,flux
def rv(a,b):return math.sqrt(sum((x-y)**2 for x,y in zip(a,b)))/math.sqrt(sum(x*x for x in b))
cases={"affine":[[0,0,0],[1,0,0],[1,1,0],[0,1,0]],"distorted":[[0,0,0],[2,0,0],[1.5,1,0],[0,1,0]]}
for w in [.1,.5,2]:cases["saddle"+str(w)]=[[x,y,w*x*y]for x,y in [(-1,-1),(1,-1),(1,1),(-1,1)]]
ans={}
for name,X in cases.items():
ds=[]
for n in [3,5,7]:
C,K=mats(X,n);T,q=response(X,C,K);ds.append((C,K,T,q))
rows=[]
for a,b in [(ds[0],ds[1]),(ds[1],ds[2])]:
rows.append([rel(a[0],b[0]),rel(a[1],b[1]),rv(a[2],b[2]),max(math.sqrt(dot([x-y for x,y in zip(q,r)],[x-y for x,y in zip(q,r)]))for q,r in zip(a[3],b[3]))/max(math.sqrt(dot(r,r))for r in b[3])])
ans[name]=rows
print(json.dumps({"quadrature":ans},indent=2))
def direct54(X):
C=zeros(20);K=zeros(20)
for xi,wx in gauss(3):
for eta,wy in gauss(3):
N,G,J,no=geometry(X,xi,eta)
for off in [0,2]:
for s,ws in gauss(3):
ls=[s*(s-1)/2,1-s*s,s*(s+1)/2]
ds=[4*(s-.5),-8*s,4*(s+.5)]
H=[0.]*20;B=zeros(20,3)
for i in range(4):
for a in range(3):
p=5*i+off+a;H[p]=N[i]*ls[a]
B[p]=[G[i][j]*ls[a]+no[j]*N[i]*ds[a]for j in range(3)]
weight=wx*wy*ws*J/4
for i in range(20):
for j in range(20):
C[i][j]+=weight*H[i]*H[j];K[i][j]+=weight*dot(B[i],B[j])
return C,K
for name,X in cases.items():
C,K=mats(X,3);Cd,Kd=direct54(X)
print("54point",name,rel(C,Cd),rel(K,Kd))
assert rel(C,Cd)<1e-12 and rel(K,Kd)<1e-12
# Generalized eigen residual and C orthogonality, independently reassemble.
C,K=slabmat(2);Cf=[r[1:]for r in C[1:]];Kf=[r[1:]for r in K[1:]]
res=fn([[a-b*l for a,b,l in zip(kr,cr,lam)]for kr,cr in zip(mm(Kf,W),mm(Cf,W))])/fn(mm(Kf,W))
orth=fn([[v-float(i==j)for j,v in enumerate(row)]for i,row in enumerate(mm(mm(tr(W),Cf),W))])
print("slab eigen residual",res,"C orthogonality",orth)
assert res<1e-12 and orth<1e-12
def spatial_bounded(n,steps=32,t=.001,ng=12):
h=1/n;l=6*(1-math.cos(math.pi/n))/(h*h*(2+math.cos(math.pi/n)))
af=math.exp(-steps*math.log1p(l*t/steps));ase=math.exp(-l*t);ae=math.exp(-math.pi**2*t)
et=eg=etime=ese=0.
for i in range(n):
for s,w in gauss(ng):
x=(i+(s+1)/2)*h;measure=w*h/2
v=((1-s)*math.cos(math.pi*i/n)+(1+s)*math.cos(math.pi*(i+1)/n))/2
g=(math.cos(math.pi*(i+1)/n)-math.cos(math.pi*i/n))/h
et+=(af*v-ae*math.cos(math.pi*x))**2*measure
eg+=(af*g+math.pi*ae*math.sin(math.pi*x))**2*measure
etime+=((af-ase)*v)**2*measure
ese+=(ase*v-ae*math.cos(math.pi*x))**2*measure
return math.sqrt(et),math.sqrt(eg),math.sqrt(etime),math.sqrt(ese)
rows=[]
for n in [8,16,32,64]:
a=spatial_bounded(n);b=spatial_bounded(n,ng=16)
ratio=a[2]/a[3];quad=max(abs(x-y)/abs(x)for x,y in zip(a,b))
assert ratio<1/32 and quad<1e-8
rows.append(a)
fluxbytes=n*n*24*3*32*8
tempbytes=(n+1)**2*5*33*8
print("bounded space",n,"T,grad,time,semi",a,"time/space",ratio,"quad",quad,"historybytes",fluxbytes+tempbytes)
for a,b in zip(rows,rows[1:]):
pt=math.log(a[0]/b[0],2);pg=math.log(a[1]/b[1],2)
assert 1.95<=pt<=2.05 and .98<=pg<=1.03
print("bounded order",pt,pg)
time_orders=[]
for n in [16,32,64,128]:
e1=abs(math.exp(-n*math.log1p(1/n))-math.exp(-1))
e2=abs(math.exp(-2*n*math.log1p(1/(2*n)))-math.exp(-1))
order=math.log(e1/e2,2)
assert .97<=order<=1.02
time_orders.append(order)
print("fixed spatial time orders",time_orders)
```
## B. Exact Fraction integration·rank·jump
```python
from fractions import Fraction as F
def fmat(A,d):return [[F(v,d)for v in r]for r in A]
def mul(a,b):
c=[F(0)]*(len(a)+len(b)-1)
for i,x in enumerate(a):
for j,y in enumerate(b):c[i+j]+=x*y
return c
def der(a):return [i*a[i] for i in range(1,len(a))]
def integ(a):return sum(F(x,i+1)for i,x in enumerate(a))
sh=[[1,-3,2],[0,4,-4],[0,-1,2]]
assert [[integ(mul(a,b))for b in sh]for a in sh]==fmat([[4,2,-1],[2,16,2],[-1,2,4]],30)
assert [[integ(mul(der(a),der(b)))for b in sh]for a in sh]==fmat([[7,-8,1],[-8,16,-8],[1,-8,7]],3)
M=fmat([[4,2,1,2],[2,4,2,1],[1,2,4,2],[2,1,2,4]],36)
S=fmat([[4,-1,-2,-1],[-1,4,-1,-2],[-2,-1,4,-1],[-1,-2,-1,4]],6)
fac=[([1,-1],[1,-1]),([0,1],[1,-1]),([0,1],[0,1]),([1,-1],[0,1])]
assert [[integ(mul(a,c))*integ(mul(b,d))for c,d in fac]for a,b in fac]==M
assert [[integ(mul(der(a),der(c)))*integ(mul(b,d))+integ(mul(a,c))*integ(mul(der(b),der(d)))for c,d in fac]for a,b in fac]==S
Z0=[[F(0)for _ in range(5)]for _ in range(5)]
Z1=[[F(0)for _ in range(5)]for _ in range(5)]
for off in [0,2]:
for i in range(3):
for j in range(3):
Z0[i+off][j+off]+=F([[4,2,-1],[2,16,2],[-1,2,4]][i][j],30)
Z1[i+off][j+off]+=F([[7,-8,1],[-8,16,-8],[1,-8,7]][i][j],3)
C=[[M[i//5][j//5]*Z0[i%5][j%5]for j in range(20)]for i in range(20)]
K=[[S[i//5][j//5]*Z0[i%5][j%5]+M[i//5][j//5]*Z1[i%5][j%5]for j in range(20)]for i in range(20)]
def pivots(A):
a=[r[:]for r in A];p=[]
for k in range(len(a)):
d=a[k][k];p.append(d)
for i in range(k+1,len(a)):
for j in range(k+1,len(a)):a[i][j]-=a[i][k]*a[k][j]/d
return p
assert all(x>0 for x in pivots(C))
assert all(x>0 for x in pivots([r[:19]for r in K[:19]]))
assert all(sum(r)==0 for r in K)
assert sum(map(sum,C))==2
cc=[[F(2),F(1)],[F(1),F(2)]];kk=[[F(1),F(-1)],[F(-1),F(1)]]
old=[F(0),F(0)]
for expected,react in [(F(-1,5),F(24,5)),(F(1,25),F(36,25))]:
new=[(2*(2*old[0]+old[1])-1)/5,F(1)]
r=[2*sum(cc[i][j]*(new[j]-old[j])for j in range(2))+sum(kk[i][j]*new[j]for j in range(2))for i in range(2)]
assert new[0]==expected and r==[0,react]
assert 3*sum(new[j]-old[j]for j in range(2))==react/2
old=new
print("exact: C SPD 20/20; K rank19; K1=0; 1TC1=2; jump -1/5,1/25; Rc 24/5,36/25")
```
## C. Time predicate와 경계
```python
import math,random,struct
from fractions import Fraction as F
def input_cell(x):
v=F(x);prev=F(math.nextafter(x,-math.inf))
nxt=math.nextafter(x,math.inf)
hi=(v+F(nxt))/2 if math.isfinite(nxt) else v+(v-prev)/2
even=not(struct.unpack('>Q',struct.pack('>d',x))[0]&1)
return (prev+v)/2,hi,even
def cell_overlap(t,d,m):
a,b,ca=input_cell(t);c,e,cb=input_cell(d);c*=m;e*=m
lo=max(a,c);hi=min(b,e)
return lo<hi or (lo==hi and (ca or (a<lo<b)) and (cb or (c<lo<e)))
def dyadic_overlap(t,d,m):
_,ex=math.frexp(t)
factor=F(2)**(-ex)
aa,bb,ca=input_cell(t);cc,ee,cb=input_cell(d)
ints=[]
for x in [aa,bb,cc,ee]:
z=x*factor*2**88
assert z.denominator==1
ints.append(z.numerator)
a,b,c,e=ints;c*=m;e*=m
assert max(z.bit_length()for z in [a,b,c,e])<=122
lo=max(a,c);hi=min(b,e)
return lo<hi or (lo==hi and (ca or (a<lo<b)) and (cb or(c<lo<e)))
def fm(a,b,c):return float(F(a)*F(b)+F(c))
def grid_cell(t,d,inc=100):
if not(math.isfinite(t)and math.isfinite(d)and t>0 and d>0 and type(inc)is int and 1<=inc<=2**32-1):return "input-error"
if d>=t:return 1,False,t
if F(t)/F(d)>inc+1:return "INC"
_,e=math.frexp(t);q=math.ldexp(t,-e)/math.ldexp(d,-e);m=round(q)
close=cell_overlap(t,d,m)
assert close==dyadic_overlap(t,d,m)
if close:n=m
else:
factor=F(2)**(-e+88);X=F(t)*factor;Y=F(d)*factor
assert X.denominator==Y.denominator==1
X=int(X);Y=int(Y);f=math.floor(q)
if f*Y>X:f-=1
elif (f+1)*Y<=X:f+=1
assert f*Y<=X<(f+1)*Y
n=f if X==f*Y else f+1
if n>inc:return "INC"
if not close:assert n==-(-F(t)//F(d))
delta=d if close else fm(-(n-1),d,t)
assert 0<delta<=d and math.isfinite(delta)
assert (n-1)*d<t
return n,close,delta
pairs=[(10.,3.),(2.,3.),(.9,.3),(math.nextafter(.9,math.inf),.3),
(5.351455173733616,.17262758624947144),
(20*math.ulp(0.),6*math.ulp(0.)),(17*math.ulp(0.),5*math.ulp(0.)),
(math.nextafter(math.inf,0.),math.nextafter(math.inf,0.)/3)]
for t,d in pairs:print("input-cell grid",repr(t),repr(d),grid_cell(t,d))
rng=random.Random(901)
count=0
for i in range(10000):
m=rng.randint(2,2**32-1);d=math.ldexp(rng.uniform(.5,1),rng.randint(-1073,990))
t=m*d
if not(math.isfinite(t)and d>0):continue
t=math.nextafter(t,math.inf if i%2 else -math.inf)
if not math.isfinite(t):continue
q=t/d
if not math.isfinite(q)or q>2**32:continue
mm=round(q)
assert cell_overlap(t,d,mm)==dyadic_overlap(t,d,mm)
count+=1
print("input-cell/fixed-dyadic independent comparisons",count)
assert grid_cell(5.351455173733616,.17262758624947144)[0]==31
assert grid_cell(math.nextafter(.9,math.inf),.3)[0]==4
assert grid_cell(20*math.ulp(0.),6*math.ulp(0.))[0]==3
assert grid_cell(17*math.ulp(0.),5*math.ulp(0.))[0]==4
print("INPUT-CELL ASSERTIONS PASSED")
for t,d in [(1.,.1),(math.ldexp(.9,-1000),math.ldexp(.3,-1000)),(math.ldexp(.9,1000),math.ldexp(.3,1000)),(1e-323,5e-324)]:
print("input-cell scale",repr(t),repr(d),grid_cell(t,d))
assert grid_cell(10.,3.,3)=="INC" and grid_cell(10.,3.,4)[0]==4
assert grid_cell(float(2**32),1.,2**32-1)=="INC"
for bad in [0,-1,1.5,2**32]:assert grid_cell(1.,.1,bad)=="input-error"
for bad in [0.,-1.,math.inf,math.nan]:assert grid_cell(bad,.1)=="input-error"
print("CURRENT TIME GRID ASSERTIONS PASSED")
```
## D. 입력 rounding-cell의 독립 반증 검색
고정 seed의 검색은 폐기된 arithmetic-only predicate와 입력 rounding-cell의 차이를
검출한다. 이 출력은 predicate 실행 성공과 승인 시간 의미의 정합을 구별하기 위한
반증 evidence다. 상세 판정은 numerical-review.md §8을 따른다.
```python
import math,random
from fractions import Fraction as F
def rad(x):return max(x-math.nextafter(x,-math.inf),math.nextafter(x,math.inf)-x)/2
def pred(t,d,m):
q=t/d;p=m*d
B=math.nextafter(d*rad(q)+rad(p),math.inf)
return abs(float(F(t)-m*F(d)))<=B,B
rng=random.Random(41);found=[]
for z in range(20000):
d=rng.uniform(.125,.25);m=rng.randint(2,100);t=m*d
for dr in [-1,1]:
tt=math.nextafter(t,math.inf if dr>0 else -math.inf)
loT=(F(tt)+F(math.nextafter(tt,-math.inf)))/2
hiT=(F(tt)+F(math.nextafter(tt,math.inf)))/2
loD=(F(d)+F(math.nextafter(d,-math.inf)))/2
hiD=(F(d)+F(math.nextafter(d,math.inf)))/2
overlap=max(loT,m*loD)<=min(hiT,m*hiD)
ok,B=pred(tt,d,m)
if overlap and not ok:
found.append((tt,d,m,float(abs(F(tt)-m*F(d))),B));break
if found:break
print("input-cell overlap rejected by retired arithmetic predicate:",found)
t=math.nextafter(.9,math.inf);d=.3;m=3
print("specific .9 upper:",float(F(t)-m*F(d)),rad(t)+m*rad(d),"input-cell overlap",abs(F(t)-m*F(d))<=F(rad(t))+m*F(rad(d)))
```
## E. Library ULP 가정 없는 Fourier interval oracle
고정 dyadic 320-bit outward enclosure, Machin π, Taylor sin/cos/exp와 analytic
remainder를 사용한다. Python Fraction은 검증 generator 전용이다. 출력 interval의
outward binary64 hex endpoints를 C++ fixture에 고정할 수 있다. FESA solver와
C++ std::sin/std::exp의 ULP 보장을 가정하지 않는다.
```python
from fractions import Fraction as F
import math
BITS=320;SCALE=1<<BITS
def down(x):return F((x*SCALE).numerator//(x*SCALE).denominator,SCALE)
def up(x):return -down(-x)
def interval(a,b=None):return (down(F(a)),up(F(a if b is None else b)))
def add(a,b):return down(a[0]+b[0]),up(a[1]+b[1])
def neg(a):return -a[1],-a[0]
def sub(a,b):return add(a,neg(b))
def mul(a,b):
v=[x*y for x in a for y in b];return down(min(v)),up(max(v))
def div(a,b):
assert b[0]>0 or b[1]<0
return mul(a,(F(1)/b[1],F(1)/b[0]))
def scale(a,c):return mul(a,interval(c))
def atan_bounds(q):
a=sum((F(-1)**j)*q**(2*j+1)/(2*j+1) for j in range(64))
tail=q**129/129
return a,a+tail
a,b=atan_bounds(F(1,5)),atan_bounds(F(1,239))
PI=interval(16*a[0]-4*b[1],16*a[1]-4*b[0])
assert PI[1]-PI[0]<F(2)**-295
def exneg(x):
assert x[0]>=0
if x[0]>=512:return interval(0,F(2)**-512)
assert x[1]<=512
y=scale(x,F(1,1024));term=interval(1);s=term
for j in range(1,65):
term=scale(mul(term,neg(y)),F(1,j));s=add(s,term)
tail=F(1,2)**65/math.factorial(65)
s=add(s,interval(-tail,tail))
for _ in range(10):s=mul(s,s)
return s
def sincos_pi8(k):
k%=16;quad=(k+2)//4;rem=k-4*quad
x=scale(PI,F(rem,8));x2=mul(x,x)
st=x;ss=x;ct=interval(1);cs=ct
for j in range(1,32):
st=scale(mul(st,neg(x2)),F(1,2*j*(2*j+1)));ss=add(ss,st)
ct=scale(mul(ct,neg(x2)),F(1,(2*j-1)*2*j));cs=add(cs,ct)
ss=add(ss,interval(-F(1,math.factorial(65)),F(1,math.factorial(65))))
cs=add(cs,interval(-F(1,math.factorial(64)),F(1,math.factorial(64))))
q=quad%4
return [(ss,cs),(cs,neg(ss)),(neg(ss),neg(cs)),(neg(cs),ss)][q]
def series_interval(t):
Ts=[interval(1)for _ in range(5)];qs=[interval(0)for _ in range(5)];E=interval(1)
for m in range(80):
odd=2*m+1
x=scale(mul(PI,PI),F(t)*odd*odd/4);dec=exneg(x)
for j in range(5):
si,co=sincos_pi8(odd*j)
Ts[j]=sub(Ts[j],scale(div(mul(si,dec),PI),F(4,odd)))
qs[j]=add(qs[j],scale(mul(co,dec),2))
E=sub(E,scale(div(dec,mul(PI,PI)),F(8,odd*odd)))
assert max(hi-lo for lo,hi in Ts+qs+[E])<F(2)**-270
# Uniform analytic tail bound, using safe r upper at minimum time.
beta=scale(mul(PI,PI),F(t)/4)
first=exneg(scale(beta,161*161))
rr=exneg(scale(beta,8*81))
tailq=div(scale(first,2),sub(interval(1),rr))[1]
tailt=div(scale(first,4),scale(mul(PI,sub(interval(1),rr)),161))[1]
Ts=[add(v,interval(-tailt,tailt))for v in Ts]
qs=[add(v,interval(-tailq,tailq))for v in qs]
E=add(E,interval(-tailt,tailt))
return Ts,qs,E,tailt,tailq
for t in [.001,.01,.05,.1,.5,1.]:
Ts,qs,ee,tt,tq=series_interval(t)
assert tt<F(1,10**27) and tq<F(1,10**27)
# Certified brackets can be exported as outward binary64 hex endpoints for C++ fixtures.
width=max(float(hi-lo)for lo,hi in Ts+qs+[ee])
assert width<1e-26
print("Fourier certified",t,"max-width",width,"tailT",float(tt),"tailq",float(tq),
"E bracket",math.nextafter(float(ee[0]),-math.inf).hex(),math.nextafter(float(ee[1]),math.inf).hex())
print("FOURIER INTERVAL ASSERTIONS PASSED")
```
@@ -1,121 +1,589 @@
# Shell 과도 열전달 독립 수치 검토
- 검토일: 2026-09-28. 상태: **Proposed 계획의 독립 검토 완료 / 구현 검증 미실행**.
- 판정: 제안된 약형·행렬·BE partition·반력·한쪽 heat flux는 수학적으로 일관된다.
아래 미결정 사항 때문에 구현 착수 계약이나 제품 정확도가 승인되었다는 뜻은 아니다.
- 범위: [PRD](PRD.md)의 FEM-REQ-SHELL-THERMAL-001–018,
[정식화](formulation.md) §§1–10, [상세 계획](design-plan.md) §§1–12.
- 기준: [프로젝트 지침](../../AGENTS.md), [프로젝트 PRD](../PRD.md),
[ARCHITECTURE](../ARCHITECTURE.md), [ADR-028](../ADR.md#adr-028),
[공통 TOLERANCE](../TOLERANCE.md). 상위 계약·소스·테스트·reference는 수정하지 않았다.
- 검토일: 2026-09-29. **정식화 일관성 및 구현 전 수치 검증 설계 PASS /
solver 구현·정확도·reference·physics 검증 미실행**.
- 독립 검토 대상: formulation commit e052b6e, source 접속 범위 보완 efa5403,
입력 rounding-cell 보완 dd409b5의 실제 diff 및 §7.1–7.2.
[PRD STH-C07](PRD.md#sth-c07--독립-검증과-수치-계약의-선행-조건)의 위임 범위에서
아래 식·상수·fixture·판정 기준을 production kernel 구현 전에 고정한다.
- 상태: **DONE_WITH_CONCERNS**. 고정 5점의 초기 침투층 및 warped 3×3 적분 오차를
정량 기록했다. 지원 방법·geometry·외부 tolerance 변경을 제안하거나 승인하지 않는다.
- 기준: [정식화](formulation.md), [승인 ADR](ADR.md), [공통 수치 기준](../TOLERANCE.md),
[개발 규칙](../standards.md), [수치 주의점](../engineering-notes.md).
재현 계산은 [scratch appendix](numerical-oracle-scratch.md)에 있다.
Reviewer는 정식화 작성자가 아니며 production·test·reference를 수정하지 않았다.
## 1. 영향이 큰 판단과 필요한 결정
## 1. 판정과 이전 검토의 정리
| ID / 판정 | 근거·영향 | 필요한 조치 / owner |
| --- | --- | --- |
| NR-01 / 물리 결정 대기 | Midsurface metric을 두께 전체에 쓰는 product-domain은 실제 offset volume의 곡률 항을 생략한다. 같은 normalized ζ를 연결한 두께 단차·fold는 실제 3D 접합부가 아니다. | 메인/사용자: 이 이상화를 승인하고 지원 한계를 명시. Numerical/I/O: 곡률·왜곡·퇴화·orientation 판정과 허용 geometry 범위를 구현 전에 확정. 높은 적분 차수로 모델링 오차가 제거된다고 주장하지 않음. |
| NR-02 / 초기 상태 계약 대기 | Raw 초기 온도를 첫 BE 증분의 이전 상태로 쓰는 식은 맞지만 t=0 즉시 projection과 첫 RHS·에너지가 다르다. Reference 입력에는 step 전 `*Boundary / Set-10,11,11`과 step 안의 `Set-10,11,11,1000.`이 모두 있다. | 메인/사용자: 첫 증분 jump 의미 결정. I/O: 생략 초기값, 생략 boundary 값, pre-step boundary와 step 재정의의 의미를 명시. 단일 transient step 지원을 이유로 pre-step 데이터를 무시하지 않음. |
| NR-03 / 공개 입력 결정 대기 | Edge kernel의 면적당 입력과 `q_in hL`은 일관되지만 DS4 `*DFLUX, S1…S4` 대응은 확인되지 않았다. 기존 reference에는 edge 하중이 없다. | I/O/메인: 공개 syntax·단위·대상 edge·중복 입력 의미 확인. 입증되지 않은 keyword 지원을 선언하거나 line-rate에 h를 중복 곱하지 않음. |
| NR-04 / 구현 전 수치 계약 부족 | 3×3 표면 적분은 일반 warped Q4에서 exact가 아니다. Consistent C와 P2 두께 보간은 maximum principle을 보장하지 않는다. 고정 5점은 초기 두께 경계층을 임의로 해상할 수 없다. | Numerical owner: 아래 geometry·적분·수렴·내부 residual 판정을 고정하고 독립 사례로 검증. 실패 뒤 clamp/lumping/인공 확산 또는 tolerance 완화 금지. 방법 변경은 새 검토 대상. |
| NR-05 / 문서 정합 보완 반영 | 최초 계획의 24개 integration Gauss 위치와 정식화의 54점 직접 적분·24개 trace recovery가 달랐다. 최신 계획에서 구분됨. | 유지할 계약: production 3×3 surface×exact Z, 직접 검산 54점, affine 검산 24 Gauss점, 공식 recovery 24 trace점. 숫자 24가 같다는 이유로 위치를 대체하지 않음. |
| NR-06 / 검증 범위 보완 반영 | 직렬 strip 저항식은 단면 전체 end 온도 조건을 전제로 한다. 내부 station 지정 온도와 all-constrained는 현재 공개 입력 범위 밖이다. | 최신 계획처럼 exact strip 및 all-constrained를 assembly/algebra fixture로 분류. CLI bottom/top 지정 사례는 해당 경계층을 포함한 oracle로 검증. |
2026-09-28 Proposed 검토는 당시 evidence다. NR-01의 product-domain 및 normalized
station 접속, NR-02의 raw 초기값과 BC 상속/override는 현재 승인 계약으로 해결되었다.
NR-03의 edge 입력 제안은 **011 전체 연기**로 대체된다. 과거 edge weight는 이번
구현·검증 요구가 아니다. 종료 시간의 정수배 제한 제안도 폐기되었다.
NR-01의 물리 지원 제한, NR-02의 시작 의미, 종료 시간이 Δt의 정수 배가 아니면 거부하는
제안은 사용자 요구에 영향을 주는 결정이다. Z 계수·scatter 순서·구간 Jacobian·동일
물리량의 검증 metric 고정은 합의된 물리를 구현하기 위한 상세 계약이다. 영속 flux
좌표계·trace identity·t=0 validity는 I/O 계약이며 이 보고서가 대신 승인하지 않는다.
## 2. 독립 수식·정확 산술 검산
Python 표준 라이브러리 `fractions.Fraction`을 명령 내 scratch로 사용했다. 별도 계산 파일,
production/test 변경 또는 solver 실행은 없었다. 아래 입력과 결과로 검산을 재현할 수 있다.
| 항목 | 검산과 결과 |
| 항목 | 독립 판단 / 현재 경계 |
| --- | --- |
| DOF/shape | Node-major `p=5(i−1)+a`, a=1…5. 각 길이 1인 두께 구간의 좌표 u∈[0,1]에서 shape `[1−3u+2u², 4u−4u², −u+2u²]`를 곱하여 항별 정확 적분. |
| 두께 질량 | `Z0_seg=(1/30)[[4,2,−1],[2,16,2],[−1,2,4]]` 재현. 조립 뒤 `Z0·1=(1,4,2,4,1)/6`, 합 2. |
| 두께 전도 | 위 다항식을 u로 미분하여 적분하면 `Z1_seg=(1/3)[[7,−8,1],[−8,16,−8],[1,−8,7]]`. 조립 뒤 `Z1·1=0` 정확 성립. |
| 체적/계수 | `dz=(h/2)dζ`, `∂z=(2/h)∂ζ`이므로 `C=ρc(h/2)M⊗Z0`, `K=k[(h/2)S⊗Z0+(2/h)M⊗Z1]`. 접선과 normal의 직교성으로 교차항 0. C는 energy/temperature, K는 power/temperature. |
| 평판의 rank/에너지 | 단위 정사각형, `h=2`, `k=ρc=1`의 수치 단위 fixture. `M=(1/36)[[4,2,1,2],[2,4,2,1],[1,2,4,2],[2,1,2,4]]`, `S=(1/6)[[4,−1,−2,−1],[−1,4,−1,−2],[−2,−1,4,−1],[−1,−2,−1,4]]`. C의 정확 LDL pivot 20개 모두 양수, K의 leading 19×19 pivot 모두 양수, `K·1=0`, `1ᵀC1=2` 정확 성립. 즉 C SPD, K PSD/rank 19. |
| Edge weight | 한 endpoint당 `q_in hL(1,4,2,4,1)/24`; 두 endpoint의 합은 `q_in hL`. 두께 station에 1/5씩 나누면 잘못된 일관 하중이다. |
| 적분 차수 | Affine Q4 면내 M/S는 2×2 exact. 두께 mass·면내 전도는 4차여서 구간별 3점, 두께 전도는 2차. Production 3×3×2×3=54점은 exact Z 분리식과 같은 quadrature를 표현한다. 실제 54점 코드의 동등성은 미검증. |
| 약형·분리 행렬·부호 | Q4×두 C0 P2, 실제 midsurface metric, exact Z, consistent C와 BE는 일관된다. 단위·node-major 순서·법선 방향을 확인했다. |
| Whole patch | affine normal image의 separating vector certificate 및 bilinear midpoint identity의 injectivity 유도가 성립한다. 적분점만 검사하는 것은 불충분하다. |
| Source 접속 | 전역 node-edge scan은 disconnected/contact 범위를 침범할 수 있었다. efa5403의 source endpoint를 공유하는 coarse/partial edge 후보 제한을 독립 확인했다. |
| 초기 jump | raw t=0 보존 및 첫 RHS의 C_fc T_c_old 항이 필요하다. Slab의 δ→0 자유 초기 극한은 분석용 oracle이며 raw 초기 dataset을 바꾸는 projection 절차가 아니다. |
| 시간 | 표현오차 수준 endpoint-close와 실제 shortened δ를 구별하는 predicate를 검산했다. 저장 시각의 subtraction과 연산 δ를 혼동하지 않는다. |
| Recovery | 4 면내점×6 두께 trace의 global q와 R=[e1,n×e1,n], 0−/0+ 분리, full residual 반력 부호가 일관된다. t=0 flux/reaction은 없다. |
| 남는 위험 | consistent C의 초기 undershoot, 고정 두께 해상도, metric conditioning, 3×3 유한 적분 오차. Clamp/lumping/차수 변경/지원 축소로 숨길 수 없다. |
유효 geometry와 양의 적분 weight에서는 `vᵀCv>0`, `vᵀKv≥0`의 에너지 형태가 유지된다.
전역 K의 상수 null mode는 **연결된 thermal component마다 하나**이며, 분리된 mesh를
항상 rank N−1로 검사하면 안 된다. 모든 활성 DOF가 양의 체적 요소에 연결되면
`A=C/Δt+K`는 구속 없이도 SPD다. 극단적인 h 또는 Δt는 conditioning을 악화시킬 수
있으므로 factorization 성공만으로 정확도를 판단하지 않는다.
PASS는 **설계의 수학적 일관성과 독립 사전 oracle**에 한정한다.
실제 MSVC solver·HDF5·외부 reference 결과가 통과했다는 뜻은 아니다.
첫 증분 partition의 독립 반례 fixture는 다음과 같다. 일관된 단위로
`C=[[2,1],[1,2]]`, `K=[[1,−1],[−1,1]]`, `Δt=1/2`, `F=0`,
`T_raw⁰=(0,0)`, 두 번째 DOF의 지정 온도 `g=1`을 둔다.
## 2. Geometry와 mesh 판정
- 올바른 `A=[[5,1],[1,5]]`; `T_f¹=−1/5`, `T_f²=1/25`.
- `R_c¹=24/5`, `R_c²=36/25`; `E=1ᵀCT`에서 두 증분 모두 `ΔE=Δt R_c` 정확 성립.
- 최초 RHS를 일정 g의 이후 증분 식으로 바꾸면 `T_f¹=+1/5`가 되어 검출된다.
이 fixture는 signed 온도와 consistent C의 undershoot가 에너지 안정성과 모순되지
않음을 보이며 실제 20 DOF shell 해나 외부 reference 통과 증거는 아니다.
### 2.1 Whole-patch certificate와 uncertainty
`r=CΔT/Δt+KT−F`의 constrained 성분은 모델에 공급되는 양의 반력 열률이다.
전역 수지는 `ΔE=Δt(1ᵀF+1_cᵀR_c+1_fᵀr_f)`이며, free residual을 포함한 대수적
동일성 확인과 free residual 자체가 작은지의 확인을 분리한다. `q=−kBT` 경계 적분은
일반적으로 이 반력과 같지 않다. 입력 없는 단열계의 열에너지 보존과 BE의 C-norm
소산도 서로 다른 검증이다. 초기 t=0 반력은 미정의이고 material history update는 없다.
정식화 §2.1의 r=a+bξ+cη+dξη, w=w0+uξ+vη를 사용한다.
Convex set W={w(ξ,η)}의 원점에 가장 가까운 w*는 w*·w≥||w*||²를 만족한다.
따라서 네 corner에 대한 **같은** 고정 double 벡터 s의 양의 interval dot product가
전체 patch를 증명한다. 후보 최소점 계산이 부정확해도 최종 certificate가 양수이면
안전하다. 후보 탐색 실패를 수학적 퇴화의 증명으로 보고하면 안 된다.
r(p)−r(q)=Dr((p+q)/2)(p−q)이므로 전체 patch에서 rank(Dr)=2이면 단사다.
이 증명은 bilinear map 전용이며 일반 고차 curved element로 확대하지 않는다.
## 3. 구현 후 반증 가능한 검증 계획
**고정 규칙:** binary64, u=2^-53, epsilon=2^-52, round-to-nearest,
gradual underflow; 입력 double은 singleton interval이다. 각 사칙연산 endpoint를
양쪽 nextafter로 한 번씩 넓힌다. 0을 포함한 divisor는 사용하지 않는다.
sqrt는 비음수 endpoint를 계산한 뒤 같은 방식으로 넓힌다. 재척도는 exact
power-of-two이며 손실된 subnormal도 enclosure에 포함한다. Fast-math/FTZ와
증명되지 않은 재결합은 이 bound의 전제가 아니다.
| 사례·수준 | 독립 oracle와 판정 의도 |
모든 lower(s·w_corner)>0이면
J_lower=min lower(s·w_corner)/upper(sqrt(s·s))>0를 보존한다.
증명되지 않으면 다른 유한 candidate를 시도하고 모두 실패하면 uncertain geometry로
fail closed한다. 정확히 중복된 좌표·source node, 비유한 값, 증명된 w=0은 invalid다.
고정 Jacobian epsilon·warp/aspect cutoff는 **0개**다. Certificate 뒤에도 각 사용점의
metric inverse/gradient가 유한하고 interval에서 0 divisor가 배제되는지 확인한다.
κ(g)는 진단과 오차 전파에 쓰며 지원 cutoff가 아니다.
### 2.2 Source topology와 T-edge uncertainty
Edge incidence는 1 또는 반대 순서의 2; 3 이상·같은 방향·중복 patch는 거부한다.
각 vertex link는 하나의 path 또는 cycle이며 끊어진 fan·branch·중복 link는 거부한다.
T-edge 후보는 **같은 instance에서 source endpoint를 공유한**
(p,q)와 (p,v) 또는 (v,q)이고 p,q,v는 서로 다른 node다. Coarse 방향을 양쪽으로
시도한다. a=q−p, b=v−p의 interval에서 |a_j| lower가 가장 큰 성분 j를 택하고
λ=b_j/a_j, 잔여 성분 b_i−λa_i를 계산한다.
- 잔여 성분 interval 하나가 0을 제외하거나 λ가 [0,1] 밖이면 비접합이다.
- 공선 및 0<λ<1이 증명되면 nonconforming T-edge다.
- 공선/열린 구간 포함이 불확실하면 uncertain 접속으로 실패한다. 정확한 끝점 일치는
별도 중복 좌표 진단이며 열린 구간으로 clamp하지 않는다.
- Interval이 0을 포함한다고 exact 공선을 주장하지 않는다. Exact double-rational
검산은 선택 사항이며 새 arbitrary-precision runtime dependency를 요구하지 않는다.
Interval fail-closed만으로 계약을 충족한다.
모든 node 대 모든 edge scan, 좌표 welding, source 후보가 없는 component 접촉 검출은
하지 않는다. 같은 component의 멀리 떨어진 부분끼리의 접촉도 대상이 아니다.
이는 global contact 검증을 완료했다는 뜻이 아니다.
### 2.3 고정 geometry fixture군
| 군 | 좌표/변환과 필수 판정 |
| --- | --- |
| Shape/행렬, kernel | 위 exact M/S/Z를 기대값으로 사용. Node-major Kronecker 순서, symmetry, 상수 null, SPD, 양의 에너지, 직접 54점/분리식 동등성 확인. h를 바꿔 C 및 면내 K의 h scaling과 두께 K의 1/h scaling을 각각 검사. |
| 온도 patch/방향, kernel | 평판에서 `T=A+b·r+c z`를 주고 `q=−k(b+c n)` 확인(b는 평판 접선). 임의 rigid rotation에서 T 불변, q 회전. Connectivity 반전 시 station a↔6−a와 면 조건도 함께 반전. |
| 한쪽 trace, recovery | 면내 상수, 두께 값 `(0,0,0,1,2)`의 임의 kernel state는 `q3(0−)=0`, `q3(0+)=−4k/h`. 24개 trace 중 0−/0+가 별도 identity로 저장되는지 확인. 이 임의 state는 공개 초기조건 fixture가 아님. |
| 열용량/하중/수지, CLI+assembly | 단열·균일 초기 온도는 불변. Edge만 주면 총 입력 q_in hL; 집중 입력은 P 그대로. Prescribed DOF의 P도 full F에 남아 reaction을 바꾸는지 확인. 반력 부호, free residual, ΔE와 BE 소산을 각각 확인. |
| BE, algebra | 위 두 증분 fixture, no-constrained 및 generic all-constrained. 이전 constrained state를 raw 초기값으로 유지하고 실패 candidate가 이전 state/output을 덮어쓰지 않는지 확인. |
| 면내 해석해, CLI | 단열 직사각판 `T=T0+A cos(πx/L)exp[−α(π/L)²t]`, `α=k/(ρc)`. 초기값은 station 전체 동일. 두께 균일해가 유지되는지와 면내 temperature/flux 오차 확인. |
| 시간 수렴, 독립 semi-discrete oracle | 같은 고정 mesh의 generalized eigenpair `Kv=λCv`로 `T(t)=exp(−λt)v`, BE 증폭률 `(1+λΔt)^−n` 대조. dt, dt/2, dt/4에서 1차 오차를 확인해 공간오차와 분리. Eigenvector 초기값은 algebra fixture로 한정. |
| 두께 slab, CLI | x=z+h/2, bottom Ts, top 단열, 초기 Ti. `θ=(T−Ts)/(Ti−Ts)=Σ[m≥0]4/[(2m+1)π] sin[(2m+1)πx/(2h)] exp[−α((2m+1)π/(2h))²t]`. 양의 시간에서 5 station T와 미분한 q, 전체 energy 및 reaction 비교. Series 절단 오차를 별도 제한. |
| Mesh/두께 해상도 | Smooth affine mesh에서 면내 temperature L2의 기대 2차, gradient L2의 기대 1차를 확인. Slab은 고정 5점이므로 dt→0의 두께 공간오차 plateau와 열침투 깊이 √(αt)를 보고. 면내 refinement로 두께 오차가 사라진다고 판정하지 않음. |
| 재료·fold·두께 접속 | 직렬 저항 ΣL/(khw)는 full-end prescribed assembly fixture. CLI에서는 bottom/top 조건의 실제 oracle을 별도로 사용. Fold 등온은 q=0; 비균일 경우 global q의 동일성을 강요하지 않고 공유 DOF에 조립되는 열률 수지를 검사. |
| Warped/왜곡 geometry | 허용될 geometry군에 대해 3×3→5×5→7×7에서 K/C 및 대표 T/q의 변화 확인. Geometry mesh 오차, quadrature 오차, offset-volume 생략 오차를 분리. Degenerate·자기 교차·orientation 불일치는 음성 사례. |
| Affine | (0,0,0),(1,0,0),(1,1,0),(0,1,0); 평행사변형 (0,0,0),(2,0,0),(3,1,0),(1,1,0). Certificate PASS. |
| 평면 왜곡 | (0,0,0),(2,0,0),(1.5,1,0),(0,1,0). J>0; M exact, S의 3×3 exactness 주장 금지. |
| Saddle | r=(ξ,η,aξη), a=.1,.5,2; 자연 corner를 node로 사용. w=(-aη,-aξ,1), s=(0,0,1), J_lower=1의 exact 근거. 모두 지원. |
| 긴 aspect | rectangle 1×2^-40; 유효 certificate 필요. 조건수와 metric 계산을 별도 기록하며 aspect만으로 거부 금지. |
| Scale/이동/회전 | 유효군의 2^-200,1,2^200배; unit patch에 (2^20,-2^20,2^20) 이동. Signed axis permutation 및 z축 3-4-5 회전. 입력 반올림을 oracle interval에 포함. |
| 적분점 검사의 반례 | r=(ξ,(ξ−1/4)η,0). 3×3 모든 점의 J는 비영이지만 ξ=1/4에서 J=0. Whole-patch 거부. |
| Invalid | bow-tie 순서, 반복 source/좌표, collinear 네 node, NaN/Inf. 자동 재순서/normal flip 금지. |
| Uncertain | affine b=(1,1,0),c=(1,1+2^-52,0): cancellation interval이 0을 포함할 수 있다. 증명된 certificate 또는 uncertain 실패만 허용; 임의 epsilon/형상 보정 금지. |
| 접속 양성 | 반대 순서 edge 공유 Q4 두 개, 90도 fold, h=1/2 두께 단차, k=1/10 재료 접속; disconnected 두 patch. |
| 음성/대조 | edge 3 incidence, 동일 directed edge, node 하나만 공유한 두 fan; coarse(p,q)/partial(p,v),v=(p+q)/2. 각 음성의 별도 source identity/instance 좌표 접촉 대조군은 자동접속·T-edge 거부가 없어야 함. |
내부 roundoff 판정은 exact oracle를 기준으로 행렬 norm·연산 수·machine epsilon과
conditioning을 반영하여 구현 전에 고정한다. 예를 들어 free residual은
`||r_f||/(||A_ff|| ||T_f||+||b_f||)`, 수지 defect는 저장·입력·반력 항의 절댓값 합으로
정규화하되 분모 0은 별도 처리한다. Zero floor `max(1,scale)`는 추가하지 않는다.
Roundoff 한계, geometry margin, quadrature 변화 한계와 관찰 order의 허용 구간은 아직
확정되지 않았다. 위 기대 order 자체를 임의 5% 통과 기준으로 바꾸거나 공통 외부
tolerance를 이 검사들에 전용하지 않는다. 이 미결정 때문에 검증 실행 준비 완료 판정은 보류한다.
## 3. 내부 roundoff, 행렬 및 solve 기준
## 4. 기존 reference의 비교 계약과 관찰
### 3.1 공통 오차 언어
Read-only 원본: `reference/shell-thermal/shell-thermal.inp`,
`reference/shell-thermal/shell-thermal NodalTemperature.csv`.
이번 검토에서 실제로 읽은 위치는 `C:/git/FESADev/reference/shell-thermal`이며
현 worktree로 복사하지 않았다. 별도 `reference-model.md` 없이 이 절에 필요한 계약을 둔다.
외부 5%/1%와 무관하게 γ_m=mu/(1−mu),u=2^-53이며 mu<1을 요구한다.
행렬/shape는 Frobenius norm, solve는 ∞ norm, 물리 integral은 아래 L2 norm이다.
분모가 0이면 numerator도 정확히 0일 때만 통과한다. max(1,scale)은 금지다.
- 제안 family: `(shell-thermal, Step-1, t=37, bottom temperature, NT11, blocking)`.
온도 단위는 원본 입력의 일관된 단위이며 scalar이므로 좌표 회전과 무관하다.
- Row identity: `(instance=Part-1-1, source node label, source DOF=11)`와 step/time.
CSV에 instance 열이 없으므로 입력의 단일 instance와 source label 집합을 검증한 뒤 복원.
HDF5는 모든 station을 보존하되 비교 view는 bottom만 명시적으로 선택한다.
- Snapshot 정적 검사: 729행·729개 고유 node label, 모두 finite, 단일 Step-1/t37 frame.
min `−15.5742216`, max `1000`; `S=1000`, near-zero band `B=10`.
음수 전부가 near-zero는 아니다. `−15.5742216`은 상대오차 5% 대상이다.
- [공통 정책](../TOLERANCE.md): `|r|≤B`이면 `|f−r|≤B`, 나머지는 상대오차≤0.05,
모든 행의 relative RMS≤0.01을 함께 만족해야 한다. S=0은 exact-zero 규칙 적용.
Reference의 음수·미소값을 clamp하거나 단위/온도 영점을 이동하여 scale을 바꾸지 않는다.
- Missing/extra/duplicate/nonfinite, frame/instance/station mismatch를 tolerance 전에 거부.
현재 수행한 CSV 자체 검사는 FESA HDF5와의 일대일 대응 검증을 대신하지 않는다.
P_ij는 계산한 signed 합의 절댓값이 아니라 **각 quadrature/분리 항의 절댓값 합**이다.
Cancellation 때문에 P와 |X|를 바꿔 쓰지 않는다. 단위 affine fixture의 closed-form
rational oracle에 대한 고정 roundoff bound는 다음과 같다.
원본 NT11 한 frame은 내부 4 station, flux, reaction, edge load, 다중 재료/두께 또는
전체 시간 이력을 검증하지 못한다. 초기값·pre-step 경계·jump·consistent C·시간 적분
차이가 reference 비교에 미치는 영향은 아직 측정되지 않았다. 비교 실패는 이 차이들을
분리 조사할 근거이며 FESA가 Abaqus 내부 알고리즘과 같아야 한다는 요구가 아니다.
| 검사 | 고정 bound / 연산 수 근거 |
| --- | --- |
| M,S,Z 및 분리 C,K | \|\|X_hat−X_exact\|\|F ≤γ_4096\|\|P_X\|\|F. 점당 shape·metric·gradient·weight·entry scalar path128 이하×9+분리/scatter64 이하를 포괄한다. |
| 직접54점 대 분리 | \|\|X_54−X_sep\|\|F ≤γ_8192\|\|P_54\|\|F+γ_4096\|\|P_sep\|\|F. 54×128+reduction<8192. |
| Affine gradient·H·R | component 절댓값 전파 scale P에 error norm≤γ_2048\|\|P\|\|. Partition sum은 P=Σ\|H\|, gradient sum은 P=Σ\|\|B_p\|\|; 비영 patch는 exact gradient 크기도 포함. |
| Symmetry | \|\|X−Xᵀ\|\|F≤2E_X. E_X는 위 bound. 강제 대칭화로 검사를 통과시키지 않는다. |
| 상수 null | \|\|K1\|\|2≤sqrt(20)E_K; \|1ᵀC1−2\|≤20E_C. Fraction 검산은 exact equality. |
| h scaling | h=1/2,1,2, k/rho/c=1; C와 면내 K는 h, normal K는 1/h. 두 계산의 E를 해당 배율로 더한다. |
## 5. 수행 범위와 잔여 위험
Operation budget은 단순 loop의 최장 path 상한이다. 실제 경로가 이보다 길면 결과에 맞춰
상수를 늘리지 않고 경로/interval을 다시 검토한다. 일반 warped/translated geometry에는
affine bound를 그대로 쓰지 않는다. §2 interval을 같은 quadrature에 전파하여 X^I를
얻고 E_X=||max(|lower(X^I)−X_mid|,|upper(X^I)−X_mid|)||F를 사용한다.
각 point의 κ2(g)와 κ2(g)γ_64를 기록한다. Inverse perturbation 상대 증폭은
κ2(g)γ_64/(1−κ2(g)γ_64)이고 분모가 양수가 아니면 이 상계는 불확실하다.
직접 interval inverse가 유한 enclosure를 증명하면 사용하고 아니면 uncertain이다.
κ를 geometry 지원 조건으로 바꾸지 않는다. Roundoff와 quadrature 절단 차이는 별개다.
수행: 관련 계약의 정적 대조, 위 exact 유리수 행렬/두 증분 검산, 원본 입력·CSV의
read-only 정적 확인, 문서 정합 확인. 미수행: solver 구현·build/CTest, FESA HDF5 생성,
외부 reference 수치 비교, 곡면/접합부/수렴 실험, 실제 결과의 physics evaluation.
초기·geometry·edge syntax·내부 수치 판정·영속 결과 계약을 확정하고 위 검증을 실행해야
한다. 본 보고서의 판정은 **정식화 초안의 독립 일관성 검토**로 한정한다.
### 3.2 Exact fixture와 rank
단위 square,h=2,k=rho=c=1:
M=1/36 [[4,2,1,2],[2,4,2,1],[1,2,4,2],[2,1,2,4]],
S=1/6 [[4,-1,-2,-1],[-1,4,-1,-2],[-2,-1,4,-1],[-1,-2,-1,4]].
길이1 ζ 구간의 Z0=1/30 [[4,2,-1],[2,16,2],[-1,2,4]],
Z1=1/3 [[7,-8,1],[-8,16,-8],[1,-8,7]]을 station(1,2,3),(3,4,5)에
조립한다. C=M⊗Z0,K=S⊗Z0+M⊗Z1, node-major 순서다.
정확 유리수 LDL에서 C의20 pivot>0,K의 leading19 pivot>0,K1=0,1ᵀC1=2를
확인했다. C SPD,K PSD/rank19가 독립 증명된다.
Floating 결과는 matrix bound 외에 C의 λ_min>E_C, K를 1⊥에 제한한
λ_min>E_K와 null eigenvalue 크기≤E_K를 확인한다. Eigen 연산 자체의 backward
error enclosure도 더한다. Separation이 불가능하면 rank는 inconclusive이며
clamp로 PASS를 만들지 않는다. 전역 K는 component마다 null mode 하나다.
Gradient fixture는 단위 affine 평판의 접선 b에 대해 T=A+b·r+cz,
A=2,b=(3,-2,0),c=4로 q=−k(b+cn)이다. Warped patch의 임의 global b에
같은 식을 강제하지 않으며 그 경우 실제 reciprocal tangent의 projection을 사용한다.
Constant T=2, 면내/두께선형도 각각 검사한다. 회전 후 q는 같은 회전,
normal 반전 후 station a↔6−a와 조건까지 재매핑한다. RᵀR=I,detR=+1.
Trace state(0,0,0,1,2)는 q3(0−)=0,q3(0+)=−4k/h이며 평균하지 않는다.
### 3.3 Free solve, certified residual과 forward 전파
n=free DOF 수,A=A_ff,b=같은 δ와 full old state로 만든 RHS다.
Full residual의 free 성분 및 독립 A T_f−b를 확인한다.
β_limit=nγ_(n+1)+γ_(2n+4)+γ_(4n+8).
첫 항은 Cholesky의 |L||L|ᵀ bound와 || |L||L|ᵀ ||∞≤n||A||∞에서,
나머지는 triangular solves 및 residual 평가의 dot-product bound에서 나온 보수적인
**허용 상한**이다. Backend가 전제를 충족하지 못하면 더 엄밀한 backward-error
evidence가 필요하며 이 상한을 늘리지 않는다. n=0은 free 검사 대신 full reaction·수지를
검사한다. 큰 n의 상한을 실제 오차로 가정하여 forward budget에 넣지 않는다.
실제 sparse 행 i의 nonzero 수 ν_i로 독립 r_tilde=fl(A T_f−b)를 평가하고
E_r,i=γ_(2ν_i+2)(Σ_j|A_ij||T_f,j|+|b_i|),
R_cert=upper(|| |r_tilde|+E_r ||∞),
D_lower=lower(||A||∞||T_f||∞+||b||∞),
β_cert=upper(R_cert/D_lower).
모든 scale sum/norm은 outward interval로 계산하며 underflow 항도 §2 enclosure로
보존한다. D_lower=0은 exact-zero numerator/denominator를 별도 확인한다.
**β_cert≤β_limit**가 residual 판정이다. Full C dT/δ+KT−Q 평가에도 자신의 실제
row count·signed term absolute scale로 같은 의미를 검산한다. 구체적으로
E_r,full,i=γ_(2ν_C,i+2ν_K,i+6)(Σ|C_ij||dT_j|/δ+Σ|K_ij||Tnew_j|+|Q_i|)
+γ_1Σ|C_ij|(|Tnew_j|+|Told_j|)/δ로 둔다.
관측 residual에 roundoff enclosure를 더하므로 우연히 작은 computed residual을 그대로
정확한 잔차라고 가정하지 않는다.
n≤100의 소형 fixture는 독립 dense inverse로 Z≈A^-1를 구하고
θ=upper(||I−ZA||∞)<1을 확인한다. 그러면
upper(||A^-1||∞)≤upper(||Z||∞)/(1−θ)가 rigorous inverse-norm bound다.
κ∞(A)β_cert<1일 때 relative forward bound는
2κ∞(A)β_cert/(1−κ∞(A)β_cert).
정확 coefficient의 E_A/E_b도 numerator에 E_A||T_f||∞+E_b로 더한다.
이때 stored A의 inverse 상계를 J라 하면 J E_A<1을 확인하고 정확 A*의 inverse는
J/(1−J E_A)로 넓힌다. 따라서 absolute forward bound는
[J/(1−J E_A)](R_cert+E_A||T_f||∞+E_b)다(E_A/E_b는 ∞ norm bound).
이 denominator의 조건을 증명하지 못하면 inconclusive이며 clamp/상한 확대는 금지다.
큰 **구조화된 convergence fixture**는 inverse나 최악 β_limit 없이 다음 C-norm
전파를 사용한다. 정확 oracle C*,K*에 대한 증분 residual의 component enclosure를
R_n=|r_tilde|+E_r+E_C|dT|/δ+E_K|Tnew|+E_Q
로 둔다(E_C/E_K는 성분별 절댓값 오차 상계). 단열 homogeneous error의
(C*+δK*)^-1C*는 C* norm에서 contraction이므로
B_T^n=B_T^(n−1)+δ||R_n||2/sqrt(λ_lower(C*))
가 온도 C-norm error의 사전 고정된 누적 상계다. 초기 nodal cos의 interval 오차를
B_T^0에 포함한다. Normalized k=rho=c=h=1 fixture에서는 이 C-norm이 실제
product-volume temperature L2다. Sampling/오차 norm 적분의 interval도 더한다.
n×n unit square의 tensor mass는 Mx⊗My⊗Cz다.
Uniform Q1 1D mass의 Gershgorin lower는 1/(6n), 두 길이1/2 P2 segment의
Cz lower는 1/60이므로 **λ_lower(C*)=1/(2160n²)**.
K generalized upper는 면내 12n²+12n² 및 thickness
||Kz||∞/λ_lower(Cz)=(64/3)/(1/60)=1280으로
**λ_upper(K*,C*)=24n²+1280**.
따라서 gradient L2 solver budget은 sqrt(24n²+1280) B_T^n이다.
이 analytic bounds는 §6의 unit affine tensor fixture 전용이며 일반 mesh의 condition
보증으로 확대하지 않는다. 큰 fixture에 production sparse inverse를 요구하지 않는다.
§4 에너지 roundoff 상계의 β_limit는 보수적인 acceptance로 유지할 수 있지만,
§6–7의 oracle 차이/관찰 order budget에는 β_cert 또는 위 R_n 전파를 사용한다.
## 4. 에너지, 반력과 BE 소산
C=[[2,1],[1,2]],K=[[1,-1],[-1,1]],δ=1/2,Q=0,rawT0=(0,0),g2=1에서
T_f=(-1/5,1/25),Rc=(24/5,36/25)를 exact 확인했다. Raw constrained 항을 빼거나
먼저 g로 projection하면 실패한다. Signed undershoot를 삭제하면 안 된다.
이 fixture와 all-constrained 검사는 algebra다. Full-end strip ΣL/(khw)도 assembly
전용이며 공개 내부 station BC 지원을 뜻하지 않는다.
dT=Tnew−Told,ΔE=1ᵀC dT,Rc=(C dT/δ+K Tnew−Q)_c.
온도 영점을 뺀 두 큰 U를 차분하는 대신 ΔE를 직접 계산해도 같은 물리량이다.
**물리 수지 defect**
D_E=|ΔE−δ(ΣQ+ΣRc)|,
S_E=Σ_ij|C_ij dT_j|+δΣ_i|Q_i|+δΣ_c|Rc_i|.
m=전체 DOF 수,D_f=||A_ff||∞||T_f||∞+||b_f||∞라 하면
B_E=δnβ_limit D_f+δ||1ᵀK||1||Tnew||∞+γ_(4m+12)S_E
+δγ_(2m+6)(|| |C||dT| ||1/δ+|| |K||Tnew| ||1+||Q||1)
+γ_1Σ_ij|C_ij|(|Tnew_j|+|Told_j|).
두 번째 줄은 full residual 평가의 dot-product/division/subtraction 오차이며,
마지막 항은 온도 subtraction 오차다. D_E≤B_E와 free β≤β_limit를 **각각**
만족해야 한다. K 상수 mode 항은 §3 null bound도 만족해야 하므로 K 오류를 면제하지
않는다. 보고 metric D_E/S_E에서 S_E=0이면 exact-zero 분기다. Actual free residual을
수지에 다시 넣어 0을 만드는 대수 항등식은 물리 평형 판정을 대신하지 않는다.
단열·무입력·무구속의 별도 C-norm 소산 defect는
D_C=|0.5(TnewᵀCTnew−ToldᵀCTold+dTᵀCdT)+δTnewᵀKTnew|.
S_C는 앞의 세 C quadratic form의 **성분별 절댓값 합**에 각각 0.5를 곱하고,
마지막 K quadratic form의 성분별 절댓값 합에 δ를 곱하여 모두 더한 값이다.
B_C=δ||Tnew||1β_limit D_f+γ_(6m+20)S_C
+0.5||dT||2||C−Cᵀ||2||Told||2
+γ_1(|Tnew|+|Told|)ᵀ|C||dT|.
Residual work, 합/subtraction, symmetry defect, dT 오차 전파를 구분한 식이다.
D_C≤B_C,||Tnew||C²≤||Told||C²+2B_C를 검사하고 S_C=0은 exact-zero다.
이 제곱 norm은 실제 열에너지와 단위가 다르다.
P/−P,constrained DOF의 P,disconnected component를 별도 실행한다. 집중 입력은
signed power 그대로 Q에 가산한다. 단열 균일 T=2 보존,no-BC transient,component별
수지,정상상태 ΣQ+ΣRc=0을 검사한다. Raw flux 경계 적분과 FE 반력의 exact equality를
강제하지 않는다.
## 5. Warped quadrature의 사전 진단
§2 affine/distorted/saddle a=.1,.5,2에서 k=rho=c=h=1.
각 3/5/7 surface rule과 exact Z로 C_r,K_r를 독립 계산한다.
대표 응답은 old=0,Q=e_(node1,station11),δ=.1,noBC의 한 BE 증분이다.
d_X(r,s)=||X_r−X_s||F/||X_s||F(X=C,K),
d_T=||T_r−T_s||2/||T_s||2,
d_q=max_trace||q_r−q_s||2/max_trace||q_s||2.
| fixture | d_C(3,5)/(5,7) | d_K(3,5)/(5,7) | d_T(3,5)/(5,7) | d_q(3,5)/(5,7) |
| --- | --- | --- | --- | --- |
| distorted | roundoff/roundoff | 2.189065e-7/5.970737e-12 | 2.385149e-6/6.505495e-11 | 1.070060e-6/2.918586e-11 |
| saddle .1 | 3.944257e-7/7.758135e-13 | 3.862043e-7/7.500898e-13 | 7.916236e-7/1.286406e-12 | 9.283663e-7/1.601583e-12 |
| saddle .5 | 1.667355e-4/1.431780e-7 | 1.643476e-4/1.408553e-7 | 3.301794e-4/2.150240e-7 | 3.703751e-4/2.649276e-7 |
| saddle 2 | 5.631484e-3/2.213607e-4 | 5.618988e-3/2.217451e-4 | 6.690109e-3/1.033767e-4 | 5.483294e-3/8.641481e-5 |
표의 반올림 값을 tolerance로 쓰지 않는다. 판정은 (1) production3점이 독립3점
oracle와 §3 bound로 일치,(2) 각 fixture d(5,7)≤d(3,5)+양쪽 roundoff bound,
(3) affine 모든 차이는 roundoff 이내다. Appendix에서 재계산한 metric과 interval을
함께 보고한다. d가 roundoff 이하이면 ratio/order를 만들지 않는다. 독립54점/분리
C,K의 scratch 상대 norm 차이는 최대7.73e-16이었다.
이는 **선정 fixture의 적분 수렴 진단**이다. 7점을 exact해로 보거나 모든 허용
geometry에 같은 d 상한을 요구하지 않는다. 일반 모델의 진단이 비단조이면 조사·보고하며
지원 geometry cutoff,production 차수 변경/adaptive quadrature로 자동 교체하지 않는다.
## 6. 시간·면내 수렴
### 6.1 고정 spatial eigen oracle의 시간1차
Affine 단열 unit square,균일8×8 mesh,h=1,k=rho=c=1,
station 동일 v_ij=cos(πi/8)를 쓴다. Consistent Q1의 정확한 discrete eigenvector이며
λ_h=6(1−cos(π/8))/[(1/8)²(2+cos(π/8))].
독립 assembly로 Kv−λCv 및 C normalization을 §3 bound로 확인한다.
τ=1/λ_h,δ=τ/N,N=16,32,64,128,256,INC≥256.
T_exact=e^-1v,T_BE=(1+1/N)^-Nv,
e_N=||T_BE−T_exact||C/||v||C,p=log2(e_N/e_2N).
모든 pair의 p 허용구간은 **[.97,1.02]**.
독립 기대값은 .9817124,.9907338,.9953356,.9976599.
log(1+x)의 전개에서 e_N=e^-1/(2N)+O(N^-2)이며 N≥16의 독립 값으로 구간을
고정했다. Roundoff propagation이 e_N/1000을 넘으면 order는 inconclusive다.
공간오차와 시간오차를 분리하는 고정 spatial oracle다.
### 6.2 Fourier rectangle의 면내2차/gradient1차
Unit square,h=1,k=rho=c=1,noBC/noQ,T0=cos(πx)를 모든 station에 동일하게 둔다.
**τ=.001,N=32,δ=τ/32,INC=32**로 고정한다. 연속해는
T=cos(πx)exp(−π²t),∇T=(-πsin(πx)exp(−π²t),0,0).
Affine n×n,n=8,16,32,64를 모두 유지한다. Norm은 product-volume temperature와
global gradient L2다. 평가 적분은 cell별12×12 Gauss×exact thickness이며
12→16점 oracle 변화≤해당 오차의1e-8을 확인한다.
독립 exp(-λ_hτ) spatial oracle와 finite-N modal BE oracle를 동시에 사용한다.
BE-exponential 온도 L2 차이는 continuum semidiscrete 오차의 **1/32 이하**여야 한다.
작은 시간과 bounded N은 이 시간오차 분리를 유지하면서 모든 accepted frame을 저장한다.
별도 시간1차 검사는 §6.1에 그대로 남긴다.
| n | Semidiscrete T L2 | BE T L2 | BE gradient L2 | Time/spatial ratio |
| ---: | ---: | ---: | ---: | ---: |
| 8 | 9.903629946e-3 | 9.902647628e-3 | 2.487150504e-1 | .0001089543 |
| 16 | 2.482303232e-3 | 2.481327028e-3 | 1.245973612e-1 | .0004305809 |
| 32 | 6.209759342e-4 | 6.200013808e-4 | 6.232870148e-2 | .0017171164 |
| 64 | 1.552689976e-4 | 1.542952827e-4 | 3.116810477e-2 | .0068632659 |
모든 pair의 기존 허용구간을 유지한다: T p **[1.95,2.05]**,gradient p **[.98,1.03]**.
독립 finite-N 기대값은 T=(1.9967024,2.0007686,2.0065775),
gradient=(.9972203,.9993050,.9998262).
Q1 smooth interpolation O(h²)/O(h),λ_h=π²+O(h²)와 위 독립 계산에 근거한다.
12/16점 변화 최대1.30e-13이며 criterion1e-8보다 작다.
§3의 실제 residual 기반 B_T 및 gradient budget이 측정 오차의1/1000 이하여야
order 판정을 사용한다. 집중 source의 특이 해에 이 smooth order를 요구하지 않는다.
**실행 예산:** n=64에서도 sparse solve32회다. Mandatory flux와 initial 포함
temperature raw doubles는
4096×24×3×32×8 + 4225×5×33×8 = **81,074,472 bytes(약77.3MiB)**.
모든4개 mesh의 같은 payload 합은107,772,960 bytes(약102.8MiB), solve총128회다.
Model의 local frames는 시간마다 중복하지 않는 승인 schema이므로 largest mesh에서
추가4096×4×9×8=1,179,648 bytes다. Identity/HDF5 overhead와 실제 peak memory를
evidence에 별도 기록한다. 이 예산은 history/flux 생략이나 제품 출력 변경을 전제하지 않는다.
초기 문서의 τ=.1,N=16n²는 n64에서 flux만 약154.6GB인 과도한 fixture였으므로
**kernel 구현 전** 독립 계산으로 위 설계로 대체했다. Production 결과에 맞춘 tolerance
완화가 아니며 기존 공간 order·시간오차 분리·모든 frame 계약을 그대로 유지한다.
## 7. 고정5점 slab: 허용 이산 오차와 plateau
Unit area,h=1,k=rho=c=1,Ti=0,Ts=1,bottom step/top 단열이며 면내 node 조건은 같다.
x=(0,1/4,1/2,3/4,1) T,six trace x=(0,1/4,1/2−,1/2+,3/4,1)의 q3,
E=1ᵀCT,bottom 전체 Rc를 비교한다. Unit square1요소 및2×2 mesh다.
t=0은 rawTi이며 flux/reaction을 생성하지 않는다.
무차원 시간 t*=αt/h²,T/ΔT,q/(kΔT/h),E/(ρc area hΔT),R/(k areaΔT/h).
이 fixture의 척도는 모두1이다. T/flux는 위 위치의 max absolute error,
E/R은 absolute error다. Station 상대오차나 near-zero floor를 만들지 않는다.
### 7.1 Fourier 및 독립 P2 spectral oracle
정식화 §10 series에서 M=80,최소 t*=.001.
β=π²t*/4,r=exp[-8β(M+1)]이면 온도 tail≤
4exp[-β(2M+1)²]/[π(2M+1)(1-r)],flux/reaction tail≤
2exp[-β(2M+1)²]/(1-r),energy tail≤temperature tail.
최소 시간에서도 모든 무차원 tail은 **1e-27 미만**이다.
1−r는 -expm1 표현과 수학적으로 동등한 exact interval subtraction으로 평가한다.
**Oracle 산술은 다음 방식으로 고정**하며 C++17 libm의 막연한 1ULP 보장을 가정하지 않는다.
Appendix E가 실제 실행 가능한 표준 라이브러리 generator다.
1. π=16 atan(1/5)−4 atan(1/239), 각 alternating arctan을 j=0…63까지 exact
유리수로 합한다. 첫 생략항 q^129/129로 양쪽 enclosure를 만든다.
2. 이후 모든 interval endpoint는 격자 2^-320으로 floor/ceil하여 outward 보존한다.
원래 π enclosure 폭은2^-295 미만이다.
3. Slab station의 phase는 정수 kπ/8이다. 정수 quadrant reduction으로 |r|≤π/4에
옮긴 뒤 sin/cos Taylor 각각32항을 쓴다. Remainder는 보수적으로1/65!,1/64!다.
4. exp(-x),0≤x≤512는 y=x/1024의64차 Taylor,
remainder≤(1/2)^65/65!를 더하고 interval squaring10회로 복구한다.
x≥512이면 0≤exp(-x)≤2^-512(e>2)을 사용한다.
5. M=80 series를 위 interval로 합한다. Tail을 더하기 전 각 quantity enclosure 폭은
**2^-270 미만**이어야 한다. 여기에 앞의 analytic Fourier tail을 더한다.
필요한 binary64 fixture는 lower를 float 변환 후 nextafter(-∞),upper를 변환 후
nextafter(+∞)하여 export한다. 이 endpoint 차이가 B_oracle이며 별도 floor가 없다.
실제 generator에서 최소 시간의 tailT=1.658661e-30,tailq=4.194724e-28,
최종 enclosure 최대 폭8.389448e-28이었다. 나머지 시각의 폭은8.93e-89 이하였다.
이 방식은 테스트 oracle generator 전용이며 production arbitrary-precision dependency가
아니다. Spectral oracle의 exp도 동일 enclosure를 적용하고 eigenpair/coefficients의
backward-error interval을 함께 전파한다. Appendix A의 통상 math 계산은 표 수치의
간단한 재현용이며 interval E가 transcendental correctness의 근거다.
T10은 generator의 T/q/E/R interval을 outward binary64 fixture로 고정하여 사용할 수
있다. 또는 binary64 outward Taylor가 동일한 최종 B_oracle 및 §6–7의 오차 budget을
증명하면 그 경로도 가능하다. 2^-270은 고정 고정밀 generator의 내부 검사이며 production
정밀도 요구가 아니다. 무차원 interval error를 실제 단위로 옮길 때 T는 |ΔT|,
q는 k|ΔT|/h,E는 rho c area h|ΔT|,R은 k area|ΔT|/h를 곱한다.
온도 기준 이동은 interval 폭을 바꾸지 않는다. 선형 sampling/trace/energy/reaction row L은
성분별 |L| times input error로 전파하며 signed 값 자체를 절댓값으로 바꾸지 않는다.
독립1D P2 FE는 길이1/2인 두 segment를 조립한다.
C_seg=(1/2)/30 [[4,2,-1],[2,16,2],[-1,2,4]],
K_seg=2/3 [[7,-8,1],[-8,16,-8],[1,-8,7]].
Bottom elimination 뒤4×4 generalized eigenproblem의 Cholesky+Jacobi 계산은
3D product kernel을 호출하지 않는다.
Raw 초기 정책에서 δ→0 자유 초기 극한은
T_f(0+)=-C_ff^-1C_fc=(-5/34,3/17,-1/34,1/17).
y0=T_f(0+)−1,C-orthonormal v_j의 a_j=v_jᵀC_ff y0로
T_f(t)=1+Σa_jv_j exp(-λ_jt).
이는 첫 finite-step RHS의 극한이며 FESA t=0 state/dataset을 변경하지 않는다.
λ=(2.4686647564,22.9461660098,77.0631371712,198.6985026508).
Scratch eigen residual3.994e-16,C-orthogonality8.069e-16.
Solver 허용치는 이 숫자 자체가 아니라 §3의 coefficient/eigen/residual enclosure다.
BE oracle는 exp(-λt)를 (1+λt/N)^-N으로 대체한다.
N>1의 마지막 반력은 modal derivative -λ×final amplitude로 구할 수 있어
작은 δ의 subtraction cancellation을 피한다.
### 7.2 사전 고정 물리 오차와 허용치
각 t*마다 별도 N=4096,8192,δ=t*/N,INC=N run.
아래 첫 값은 **δ→0 P2 공간오차**, 괄호는 위로 반올림한 **D_space**다.
초기의 큰 허용치는 method의 제한이며 제품 전 범위 정확도 보증이 아니다.
| t*/sqrt(t*) | T∞ | q∞ | E | Rc |
| --- | --- | --- | --- | --- |
| .001/.03162 | .1399996(.141) | 10.64818(10.7) | .03031345(.0304) | 10.86927(10.9) |
| .01/.10000 | .01524952(.0153) | 1.274464(1.28) | .00467600(.00468) | .823687(.824) |
| .05/.22361 | .00597851(.00598) | .252470(.253) | .00129690(.00130) | .00631488(.00632) |
| .1/.31623 | .001387325(.00139) | .135427(.136) | .000665986(.000666) | .009524106(.00953) |
| .5/.70711 | .000141091(.000142) | .0283969(.0284) | .0000944866(.0000945) | .000525742(.000526) |
| 1/1 | .0000413942(.0000414) | .00815759(.00816) | .0000158234(.0000159) | .0000478014(.0000479) |
Continuum 오차 허용치는 **D_space+D_time(N)+B_oracle+B_solver**다.
D_time은 독립 modal 차이의 max/absolute norm이다. 온도/flux sampling L,
energy row1ᵀC,reaction row K_cf−C_cf C_ff^-1K_ff 각각에 대해
LΣa_jv_j[(1+λ_jt/N)^-N−exp(-λ_jt)]를 계산한다.
Cancellation이 불확실하면 Σ|La_jv_j|×|amplitude차이|로 enclosure한다.
B_oracle은 series tail+summation/eigen interval,B_solver는 §3 행렬/solve 오차를
L에 전파한 값이다. **B_oracle+B_solver>D_space/1000이면 plateau 판정은
inconclusive**다. 이 조건은 geometry 지원 cutoff가 아니다.
N=4096의 독립 D_time은 아래와 같다. N=8192에서는 대략 절반이며 appendix가 재현한다.
| t* | T∞ | q∞ | E | Rc |
| ---: | ---: | ---: | ---: | ---: |
| .001 | 7.5820e-7 | 8.7770e-6 | 4.6839e-8 | 6.0181e-6 |
| .01 | 1.6648e-5 | 1.9198e-4 | 1.8013e-6 | 1.5383e-4 |
| .05 | 2.9817e-5 | 2.2173e-4 | 7.3403e-6 | 1.9484e-4 |
| .1 | 3.0985e-5 | 1.9790e-4 | 1.1177e-5 | 1.6441e-4 |
| .5 | 6.8903e-5 | 1.1406e-4 | 4.3931e-5 | 1.0879e-4 |
| 1 | 8.0305e-5 | 1.3224e-4 | 5.1125e-5 | 1.2621e-4 |
물리 오차 허용치 외에 **동시에** FESA와 finite-N BE oracle가
B_oracle+B_solver 내에서 일치해야 한다. 큰 초기 공간오차가 잘못 구현된 방법까지
허용하지 않도록 분리했다. Table rounding은 production 결과와 무관하게 고정했다.
### 7.3 Plateau 정량 보고
E_space=||oracle_P2−Fourier||와 FESA error e_N에 대해
|e_N−E_space|≤D_time+B_oracle+B_solver를 확인한다.
면내1×1과2×2의 slab 차이는 §3 roundoff 이내여야 한다.
두께 지배 해가 면내 refinement로 좋아진다고 판정하면 실패다.
관측 plateau의 기준은 D_time+B_oracle+B_solver≤E_space/10.
N=8192에서 .001/.01/.05/.1의 네 quantity는 충족한다. 뒤 시각에서 미충족하는
T/E/R은 **미도달**로 보고하고 δ→0 spectral 공간오차를 따로 제시한다.
예컨대 t*=1의 E는 시간오차가 공간오차보다 커서 두 N만으로 plateau가 관측되었다고
할 수 없다. 침투 깊이 sqrt(αt),station 간격 h/4,segment 길이 h/2,E_space,D_time을
함께 보고한다. 초기 큰 flux/reaction 오차는 음수 온도 제거로 고칠 수 없다.
## 8. 시간격자 predicate, count와 재분해
### 8.1 독립 반증과 수정 범위
Arithmetic-only D≤d rad(q)+rad(p)는 입력 표현오차를 모두 포함하지 않았다.
τ=5.351455173733616,d=.17262758624947144,m=31에서
exact |τ−md|=8.326672684688674e-16 > 기존 Bhat=7.507374033534679e-16인데,
입력 rounding-cell I(τ)와 mI(d)는 겹친다. 이 반례를 blocking으로 보고했고
승인된 입력 표현오차 의미에 맞춰 **input-cell overlap**으로 보완했다.
이전 predicate의 실행 성공을 의미 정합 승인으로 대체하지 않았다.
Appendix D는 이 실패 이력을 보존하며 현재 판정은 C를 따른다.
### 8.2 현재 고정 predicate
입력 double x의 원래 rounding-cell은
L_x=(prev(x)+x)/2,U_x=(x+next(x))/2다. Significand LSB가 짝수면 두 endpoint
포함, 홀수면 두 endpoint 제외(nearest-even). 최대 finite의 U_x는 next=Inf를
평균하지 않고 x+(x−prev(x))/2라는 유한 real dyadic threshold다.
최소 subnormal의 half-ULP도 binary64로 먼저 반올림하지 않는다.
τ,d를 공통 power-of-two로 τ∈[1/2,1)에 옮겨 q=fl(τ/d),
m=nearest-even integer(q)를 만든다. **원래 입력의 cell도 동일 scale로 옮긴다.**
Normal화된 값의 이웃으로 cell을 다시 만들면 subnormal 입력의 표현오차를 축소하므로
잘못이다. Endpoint-close는 I(τ)와 mI(d)의 exact 교집합이 비어 있지 않은 경우다.
두 구간이 한 점에서만 닿으면 그 점이 양쪽에 모두 포함되는지도 확인한다.
INC≤2^32−1,q≤INC+1 및 위 scale에서 cell endpoint는 고정 dyadic 격자 2^-88에
정확히 놓인다. Endpoint 정수와 m곱의 보수적 최대122bit는 두 uint64 또는 네 uint32
limb 안에 들어간다. Carry/shift를 검사하는 bounded 정수 산술이며 새 arbitrary-precision
runtime dependency가 아니다. Reviewer의 Fraction은 이 bounded 비교의 독립 oracle다.
나눗셈·곱셈 roundoff는 후보/연산 bound 설명에 쓰고 cell을 더 넓히는 별도 epsilon으로
중복 가산하지 않는다.
q overflow 또는 INC+1 초과는 preflight 실패; close면 N=m이다. Nonclose에서는
X=τs2^88,Y=ds2^88,f=floor(q)로 fY≤X<(f+1)Y를 exact 비교하고 필요하면 f를
한 번 ±1 교정한다. X=fY이면 N=f, 아니면 N=f+1이다. Division 오차를 count에
그대로 전달하지 않는다. d≥τ는 한 증분 δ=τ의 선행 분기다.
INC 기본100,양의 uint32,N≤INC. 정상 stored time=fl(jd),마지막=입력τ,
strict increase/finite를 검사한다. Close의 마지막 δ는 **d 그대로**다.
실제 shortened만 fma(-(N−1),d,τ)를 쓰며 A/RHS/reaction/energy/dissipation은 같은 δ다.
| fixture | 필수 결과 |
| --- | --- |
| τ=10,d=3 | N=4,stored3,6,9,10;δ=3,3,3,1;마지막 refactor |
| τ=2,d=3 및 d=τ | 한 증분 δ=τ |
| τ=.9,d=.3 | N=3,stored 마지막 .9;세 δ=.3,factorization 재사용 |
| τ=1,d=.1 | N=10;stored .30000000000000004 등을 보존,δ=.1 유지 |
| nextafter(.9,-∞),d=.3 | .8999999999999999:close,N=3,δ=.3 |
| 두 ULP 아래 τ,d=.3 | .8999999999999998:nonclose,N=3,finalδ=.2999999999999998 |
| nextafter(.9,+∞),d=.3 | .9000000000000001:nonclose,N=4,finalδ=1.6653345369377348e-16 |
| 기존 반례 τ=5.351455173733616,d=.17262758624947144 | input cells overlap,N=31,모든 δ=d;가짜 tiny 증분 제거 |
| Subnormal midpoint | a=최소 positive subnormal:τ=20a,d=6a는 even/even 닫힌 접촉,N=3;τ=17a,d=5a는 odd/odd 열린 접촉,N=4,마지막δ=2a |
| Scale/최대값 | (.9,.3)×2^-1000 및×2^1000은 같은 count/close;τ=1e-323,d=5e-324는 N=2;τ=maxfinite,d=fl(maxfinite/3)는 N=3 |
| INC | N/N−1 경계,생략100;0,-1,1.5,4294967296 실패;τ=2^32,d=1,INC=2^32−1은 사전 count 실패 |
| 비유한/비표현 | Inf/NaN/0/음수,binary64에서0으로 underflow하는 token 거부;거대한 τ/d는 cast/allocation 전 실패 |
| 저장 가능성 | finite-positive/uint32 count 제한 아래 정상 간격 d는 ulp(jd)보다 충분히 크다. Generator의 finite/strict-increase 보호 분기는 중복/Inf 시간을 주입한 내부 fixture로 검사;불가능한 정상 공개 입력을 만들지 않음. |
특히 nextafter(.9,+∞)/.3의 tiny remainder를 유지하는 근거는 단순 predicate 결과가
아니다. exact |τ−3d|=1.6653345369377348e-16이고 두 입력 cell 반폭 합은
1.3877787807814457e-16이므로 **실제 rounding cells가 분리된다**. 이를 가짜 표현오차
증분으로 흡수하면 입력 표현오차 범위를 넘는다. 반면 위 m=31 반례는 cells가 겹치므로
추가 tiny 증분을 만들면 안 된다. 작은 실제 δ의 나쁜 conditioning은 §3과 candidate
보존으로 처리하며 nominal δ로 바꾸지 않는다.
현재 predicate의 independent seed901 10,000쌍은 원래 cell Fraction overlap과
2^-88 정수 비교가 모두 일치했고 최대122bit를 넘지 않았다. 원래 subnormal,maxfinite,
열린/닫힌 접촉,반례 및 scale 경계를 따로 검산했다. FMA oracle는 exact Fraction
product-sum을 한 번 float로 반올림한다. Production dependency 요구가 아니다.
## 9. 외부 reference와 evidence 경계
필수 원본은 reference/shell-thermal/shell-thermal.inp와
reference/shell-thermal/shell-thermal NodalTemperature.csv다.
이 검토는 원본 복사·수정·생성 또는 Abaqus 실행을 하지 않았다.
2026-09-28 read-only snapshot(729row,min=-15.5742216,max=1000)은 과거 관찰이며
이번 solver comparison을 대신하지 않는다.
Family=(shell-thermal,source Step-1,exact t=37,NT11 bottom,scalar,
사용자 일관 온도 단위,blocking),identity=(Part-1-1,source node label,DOF11).
CSV instance 생략은 입력 단일 instance와 label 집합 검증으로만 복원한다.
Missing/extra/duplicate/nonfinite/frame/instance/station 오류는 tolerance 전에 실패한다.
공통 TOLERANCE 그대로 S=max|reference|=1000,B=.01S=10;
|r|≤B이면 |f−r|≤B,그 밖에는 |f−r|/|r|≤.05,
모든729row RMS/S≤.01을 동시에 만족한다. S=0이면 모든 값 exact0이다.
−15.5742216은 relative 대상이며 clamp/온도 영점 이동 금지.
추가 Abaqus CSV는 필수가 아니다. 한 frame으로 내부 station/flux/reaction/시간 이력을
검증했다고 할 수 없다.
## 10. 실행 기록과 완료 경계
git rev-parse --show-toplevel로 T1 worktree를 확인했다(exit0).
NumPy/SciPy import 환경 probe는 NumPy 부재로 exit1;설치 없이 표준 라이브러리
scratch로 전환했다. [Appendix](numerical-oracle-scratch.md)의 명령은 Python blocks를
지정 순서로 메모리에서 실행한다. 실제 scratch exit0의 핵심 출력:
- Fraction:C positive pivot20,K positive leading pivot19/K1=0,1ᵀC1=2;
jump−1/5,1/25 및 reaction24/5,36/25.
- 직접54점/분리식:다섯 geometry 상대 norm 차이 최대7.73e-16.
- Slab eigen residual3.994e-16/C-orthogonality8.069e-16,§7 공간/시간 오차.
- §5 quadrature,§6 예상 차수,§8 time predicate/INC assertions.
최종 재현 명령은 appendix 첫 PowerShell block이며 exit 0,
마지막 출력은 ALL SCRATCH BLOCKS COMPLETED였다. Input-cell 독립 비교는
10,000/10,000 및 CURRENT TIME GRID ASSERTIONS PASSED를 출력했다.
실행 예산 보완 뒤 five-block 재검산도 exit 0이었다. Bounded-space order/time ratio,
FOURIER INTERVAL ASSERTIONS PASSED 및 최종 완료 출력을 확인했다.
상대 링크 파일 존재/Markdown 표 열 수 정적 검사도 exit 0(8 links,2 files),
두 소유 파일만 대상으로 한 git diff --check는 exit 0이었다.
수학 scratch는 production correctness를 증명하지 않는다. **미실행:**MSVC build/CTest,
구현 단위·통합 검사,FESA 결과 생성,NT11 comparison,실제 solver convergence/energy/
physics 검토. 필요한 수치 기준은 고정되었으며 구현 결과 뒤에 선택할 미정 tolerance는
없다. 필수 실패 시 구현·conditioning·oracle 전제를 조사하고 physics/method/acceptance
충돌은 메인에 보고한다.