FESADev/docs/mitc4-shell-feast-eigen/numerical-review-controls.md

MITC4 eigen 내부 수치 기준

• Feature ID: mitc4-shell-feast-eigen
• Owner: numerical-review-agent
• Input contract revision: r2
• 기본 보고서: numerical-review.md
• 적용: IEEE binary64, MSVC x64 Debug, 실제 oneMKL FEAST adapter 및 같은 산술의 내부 검증.

V1~V6, E1~E8, S1~S8의 식을 사용한다. 아래 값은 구현 전에 고정한 내부 수치 계약이며 reference fitting 값이 아니다. 외부 판정의 공통 tolerance, 기존 nonlinear Newton 종료 기준과 별개다. 실패 결과를 보고 완화하지 않는다.

C1. Scale, 공통 zero branch와 근거

ell, D, Khat/Mhat, Kf/Mf, A, s, B=A/s는 E1/E6대로 사용한다. 같은 비교의 두 행렬에는 같은 ell을 적용한다. Matrix norms는 Frobenius, vector norms는 Euclidean이며 별도 inf 표기만 infinity norm이다. Mixed force/moment 또는 translation/rotation 원값을 직접 더하지 않는다. Scaled sum-of-squares와 stable ordered sum을 쓰고 finite-first로 검사한다. 모든 ratio는 분모가 양수일 때 계산한다. numerator=denominator=0은 metric 0, denominator=0,numerator!=0은 named failure다. Zero eigenvector, nonpositive modal mass 및 invalid scale은 이 zero-pass의 대상이 아니다. 어떤 branch에도 max(1,physical_scale), mass/eigenvalue floor나 clamp를 추가하지 않는다.

공통 기본값은 tau_alg=1e-12, tau_modal=1e-10이다. 기존 linear NR-D02의 algebraic/rigid 기준을 유지하며 mass polynomial과 trace의 같은 단위 항등식에도 적용한다. Unit roundoff u=2^-53에서 gamma_k=ku/(1-ku)이다. 4096개의 단순 연산을 가정한 gamma_4096≈4.55e-13에 비해 1e-12는 작은 dense block·27점 적분의 정상 산술 오차를 구분할 여유가 있다. 이는 모든 geometry의 forward-error theorem이 아니다. Rank·conditioning·원 residual·독립 oracle를 따로 요구한다. 현재 rational oracle 및 r2 covariance 오차는 이 기준보다 충분히 작지만 기준을 그 결과에 맞춰 고른 것은 아니다.

tau_modal은 1e-12 내부 FEAST stopping보다 100배 큰 후검증 한계다. 외부 1% 판정을 내부 residual 허용값으로 쓰지 않는다. 같은 계산의 재실행 구조/row/ diagnostic 순서는 exact이고, eigenvector의 repeated-basis 비유일성은 C4를 따른다.

C2. Matrix, mass, stabilization와 reduction

ID	Metric / exact predicate	한계·zero 처리
C2-01	V1 eta_sym(Q); mass·physical K·linear complete K·whitened A	<=1e-12; Q=0이면 C1. 대칭 copy 자체만 독립 evidence로 사용하지 않음
C2-02	V1 ||R-S||F/(||R||F+||S||F); congruence, moment→matrix, 직접 적분, Cholesky reconstruction	<=1e-12; 단위·같은 scale 필수
C2-03	Scalar moment entry oracle: |Ck_ij-Ck*_ij|/[m_e(t/2)^k]	k=0,1,2 각각 <=1e-12; positive m_e=integral rho dV. Exact C1=0도 이 차원 일치 scale로 검사
C2-04	Mass total |m-m*|/(|m|+|m*|); uniform translation energy	<=1e-12; 각 방향을 따로 검사, trace(Muu)를 총질량으로 쓰지 않음
C2-05	Shifted spatial first moment difference divided by m_e ell	<=1e-12; 한 source node를 공통 Xref로 사용. hi di 포함
C2-06	M9 difference vs independent pairwise integral / (q^T ML q+q^T MC q)	absolute mismatch <=1e-12; raw Delta/scale >=-1e-12. Zero motion/drill-only는 C1. Negative roundoff를 값을 0으로 바꾸는 데 사용하지 않음
C2-07	V2 mass-null action, MZ, KZ	각각 <=1e-12; nz=0은 N/A. 정확한 analytic null만 대상
C2-08	||P^T Z||F, ||[P Z]^T[P Z]-I||F/sqrt(nf)	<=1e-12; 각 node basis의 해당 norm도 <=1e-12
C2-09	mi>0,ji>0,ji-hi²/mi>0; consistent physical20 및 모든 nonempty Mp nodal block Cholesky	strict, finite. Zero/negative pivot·variance는 실패. Consistent rank20/lumped rank5는 analytic expected rank로 검증
C2-10	|(ji-hi²/mi)-integral rho Ni(z-hi/mi)²dV|/(ji+hi²/mi)	<=1e-12; vi/ell², |hi|/sqrt(mi ji)도 기록하며 별도 rank cutoff로 쓰지 않음
C2-11	S2 여덟 초기 rotation diagonals	모두 finite >=0, 하나 이상 >0; zero 포함 분모8. kd>0 finite. Negative를 허용오차로 0 처리하지 않음
C2-12	|kd'-kd|/(|kd'|+|kd|) 및 직접 trace oracle	<=1e-12; global 회전/독립 tangent 재생성 후 계산. Trial별 frozen kd는 exact unchanged
C2-13	E3 Kzz=diag(kappa_i), kappa_i=sum incident kd,e	각 kappa finite >0, expected diagonal relation C2-02; Kpz 검사와 원 free residual을 함께 요구
C2-14	Scaled constrained component norm V3	0 exact; 복원한 source constrained entries도 exact zero

Positive/null dimension은 E2의 constrained director component exact-zero 교집합으로만 정한다. Constrained director의 tiny nonzero 성분을 축정렬이라고 반올림하지 않는다. 각 nonempty dimension≤6 scaled nodal Mp에

rcond_inf = 1/(||Mp||inf ||Mp^-1||inf) >= 2^-32 = 2.3283064365386963e-10

을 요구한다. 작은 block의 Cholesky triangular solves로 inverse action norm을 구하며 global inverse는 만들지 않는다. Product overflow를 피하는 동치 ratio로 계산한다. 이 조건은 physical rank 변경이 아니라 binary64 whitening의 numerical admissibility다. 64 연산의 gamma_64/rcond를 약 3.06e-5 이하로 제한하는 최소 guard이며 모든 eigenvalue의 상대 정확도를 보장하지 않는다. Cholesky reconstruction과 full modal residual을 통과해야 한다. 한계 미만·unresolvable estimate는 conditioning failure이고 좌표 삭제, 더 큰 mass, 다른 solver로 우회하지 않는다. 빈 nodal block은 skip하지만 global n=0 또는 요청 N>n은 명시적 failure다.

C3. Eigenpair 및 물리 내부 판정

계산된 c개 전체 eigenpair를 검사한 뒤 첫 N개만 공개한다. V3/V4의 actual M은 ML다.

ID	Metric	한계
C3-01	V1 full original free residual eta_r; V3 whitened eta_w	각각 <=1e-10
C3-02	추가 operator-scale residual rho_B=||Bx-mu x||2/||x||2	<=1e-10; dimensionless. Frobenius norm의 크기로 오류가 희석되는 것을 방지
C3-03	V3 eta_norm,i, pairwise eta_orth,ij, eta_G	각각 <=1e-10; positive mass를 먼저 확인. Gram full rank는 이 bound와 독립 column 수로 확인
C3-04	V4 energy/Rayleigh consistency eta_energy,i	<=1e-10; zero lambda도 C1/V4 분모 사용
C3-05	V5 rigid K action	각 6개 field <=1e-10; physical rigid mass Gram는 Cholesky 성공, vector 수6. Drill과 혼동하지 않음
C3-06	|c_rho lambda_i'-lambda_i|/s 또는 coordinate |lambda_i'-lambda_i|/s	<=1e-10; repeated group 내부 pairing은 C4. 동일 물리 BC·loads 사용
C3-07	Exact known-spectrum test eigenvalue error / known s	<=1e-10이고, 알려진 인접 distinct gap이 있으면 absolute error도 <=gap/8
C3-08	Simple/반복군 C4 projector error V5	<=1e-6; 테스트의 해당 군과 나머지 spectrum 사이 known gap은 >=1e-3 s로 선택

C3-02는 symmetric B에서 computed mu가 실제 spectrum으로부터 residue norm 이내에 있다는 직접 spectral-decomposition bound를 사용한다. 가장 낮은 mode인지, 올바른 multiplicity인지는 증명하지 않으므로 C5 count/endpoint evidence가 별도로 필수다.

FEAST의 res는 backend diagnostic이고 C3의 독립 residual을 대신하지 않는다. Zero eigenvalue에서 backend relative-residual의 분모가 0이면 해당 diagnostic은 unavailable로 기록하고 이유를 남긴다. Raw e/x, 자체 metric 및 최종 scalar 결과는 finite여야 한다. Undefined backend diagnostic을 physical residual=0으로 바꾸거나 mandatory finite 결과에 NaN을 저장하지 않는다. Trace epsout은 strict interval에서 finite여야 한다.

Raw mu=lambda/s에서 tau_zero=1e-9, tau_cluster=1e-8로 고정한다. mu<-tau_zero는 유의한 negative spectrum으로 procedure failure, |mu|<=tau_zero는 near-zero classification, 그 외 양수는 positive다. Boundary equality는 near-zero에 포함한다. Returned lambda는 부호와 raw 수치를 보존한다. 양수에서만 frequency를 파생하고 exact zero는 f=0, 작은 음수에서는 실수 frequency를 만들지 않는다. C3-02보다 10배 큰 near-zero band와 100배 큰 cluster band는 수치 분류용이다. 양의 작은 mode를 제거하거나 N을 세는 방식에 쓰지 않는다. s/λ의 nonfinite 또는 nonzero mu의 λ 환산 underflow-to-zero는 실패하고 physical scale을 1로 바꾸지 않는다.

C4. Repeated subspace와 결정론

Raw lambda 오름차순의 stable sort, exact equal tie는 original solver column 순서다. 인접 |lambda_i-lambda_(i-1)|/s<=1e-8의 connected group은 diagnostic cluster일 뿐, distinct 값을 merge/평균하거나 vector를 섞을 권한이 아니다. 알려진 exact repeated 군의 density/rotation tests에서 V5의 mass-whitened projector norm을 사용한다. Known gap≥1e-3 s와 C3-02이면 first-order invariant-subspace sensitivity는 약1e-7 수준이므로 C3-08의 1e-6은 방향 비유일성을 허용하면서 잘못된 공간을 검출한다. 이 조건 밖 실제 cluster는 raw eigenpair·Gram/residual을 계속 검사하고 특정 개별 vector 방향의 정확도를 새로 약속하지 않는다. Mode shape 외부 comparison은 없다.

Source DOF 순서로 largest absolute scaled global component를 찾고 sign을 양수로 정한다. 동률은 첫 source DOF다. Zero vector는 실패다. 동일 binary/runtime/thread/입력의 세 반복에서는 sparse structure·COO reduction·row/diagnostic identity가 exact 동일해야 한다. Scalar 재현 오차는 C3-06, repeated vector는 C4를 적용한다. CPU/runtime를 바꾼 bitwise eigenvalue/vector equality를 요구하지 않는다.

C5. FEAST controls와 bounded interval 정책

Actual call은 dfeast_scsrgv, (B,I), uplo='U', sorted duplicate-free 1-based CSR다. 각 call은 feastinit하고 아래 C-index fpm을 설정한다. 명시하지 않은 값은 해당 설치의 초기값을 유지·기록한다. Private backend를 넘어 vendor type을 노출하지 않는다.

fpm	값	근거·역할
[0]	0	Status를 deterministic diagnostic으로 수집
[1]	16	공식 허용 contour 값; 기본8보다 적분 resolution에 여유, 모든 attempt 고정
[2]	12	1e-12 double stopping; 독립 C3 residual 필수
[3]	20	유한 refinement 한계; 미수렴 성공 금지
[4]	1	아래 고정 initial subspace 사용
[5]	1	Endpoint-scaled trace stopping; zero mode에서 C3-02로 독립 검증
[13]	0	Eigenpairs가 아닌 subspace-only 결과 금지
[26], [27]	1, 1	CSR validity와 second-matrix SPD check
[49]	count 1, solve 0	Total interval count와 실제 eigenpair solve 분리
[63]	0	초기 PARDISO 정책 사용; 임의 iparm tuning 금지

근거는 Intel parameter table과 E7의 API/status 출처다. Count는 documented total이며 stochastic estimate로 부르지 않는다. 공식 count/status의 floating-point 한계가 있으므로 아래 redundant evidence와 실제 known spectra가 필수다.

Gershgorin 각 off-diagonal absolute sum은 매 덧셈 뒤 nextafter(+infinity)로 outward 누적한다. Bii-ri/Bii+ri도 각각 바깥으로 한 ULP 이동한다. Overflow는 실패다. beta=||B||inf를 같은 outward 방식으로 구한다. Nonzero B에서는 beta>0; 정확한 B=0은 E6의 algorithmic enclosure [-1,1]에 대응하여 beta=1로 둔다. 이는 physical mass/λ floor가 아니다. ulp(v)는 v의 두 finite adjacent spacing 중 큰 값이며

delta(v)=1e-8 beta + 64 ulp(v)

로 정의한다. 초기 endpoints는 a0=nextafter(gL-16delta(gL),-infinity), b0=nextafter(gU+16delta(gU),+infinity)다. Finite strict enclosure와 total count=n을 확인한 뒤 E6의 upper-prefix search와 empty-prefix lower 축소를 한다.

Interior endpoint v를 bracket에 반영하려면 (a0,v-delta(v)), (a0,v), (a0,v+delta(v))의 count가 모두 info=0이고 같아야 한다. 이 guard count는 eigenvalue 근처의 endpoint를 피하는 추가 evidence다. Guard가 enclosure 밖이면 그 부분은 알려진 outside bound를 이용하며 실제 interval이 strict하지 않으면 평가하지 않는다. Zero count는 성공 count=0일 때만 믿고 info=1을 empty proof로 바꾸지 않는다. 이미 알려진 prefix count와 monotonicity가 맞지 않으면 전체 실패다.

제어	유한 값 / 정확한 처리
Total count API calls	전체 analysis 최대512. 성공·실패·guard·최종확인·retry 모두 포함
Search optimization	Upper/lower 각각 최대32 midpoint 단계; optimization count calls 최대384로 남은128 예약
Endpoint collision	같은 midpoint에서 offsets 0,+1,-1,+2,-2,+4,-4 times delta, 이 순서 최대7후보. 현재 bracket 내부 후보만 사용
Search 종료	E6의 count target N<=c<=min(n,2N+8) 또는 폭 <=8delta(mid)/새 representable midpoint 없음/위 budget. 마지막 certified c>=N bracket 사용
최종 확인	Enclosure 안의 a에 prefix0, b에 prefix c>=N, 두 endpoint guard, actual (a,b) count=c. 불일치/예약 budget 부족이면 실패
Initial solve capacity	E6 그대로 m0=min(n,max(c+8,ceil(3c/2))); c=n이면 n
Actual solve calls	최초 포함 최대3; 각 call loops≤20, control/seed 동일
허용 solver retry	info1/2/3, count/returned-m mismatch 또는 endpoint separation failure만. 저장된 더 넓은 certified bracket으로 확대·재count, m0=min(n,max(new formula,2*previous m0))
Retry bracket 순서	Upper search에 저장된 현재 b보다 큰 성공 b를 오름차순으로 선택, lower는 이전의 더 작은 certified prefix0 또는 a0; 최종 fallback은 원 enclosure
금지 retry	Nonfinite eigenpair, wrong Gram/residual, negative spectrum, invalid/SPD/index/allocation/inner-solver 오류는 즉시 실패. Seed·mass·solver·threshold 변경 금지

새 유효 bracket이 없으면 같은 실패 call을 무조건 반복하지 않는다. Optimization 한계에서는 이미 certified enclosure를 사용할 수 있지만 실제 solve·최종 count가 실패하면 부분 결과를 공개하지 않는다. Returned m=c, interval 내부와 mu-a>=delta(a), b-mu>=delta(b), 모든 c columns C3 통과를 확인한다. N번째가 cluster 내부여도 전체 c를 검사하고 출력은 N개다. Bounds·prefix/guard/status/count/m0/m/loop/epsout·chosen ranks 및 모든 실패를 보존한다. 후속 tolerance 완화나 reference Value를 bounds 설정에 사용하는 것은 금지한다.

C6. Initial subspace, thread와 storage

Count scratch는 E7의 conservative n×n, identity initial x다. 실제 solve에서 m0=n도 identity다. m0<n이면 uint64 recurrence를 column-major로 한 entry마다 한 번 실행한다. Seed 0x4D495443344D4645, z <- 6364136223846793005*z+1442695040888963407 mod 2^64, x=2*(double(z>>11)*2^-53)-1이다. Overflow wrap은 이 generator에만 의도적이며 size/index 산술에는 허용하지 않는다. Seed를 바꾸며 성공을 찾지 않는다.

Column ascending의 nonpivoted Householder QR를 사용한다. Scaled norm, sign(0)=+1, v1=x1+sign(x1)||x||로 cancellation을 피하고 최종 R diagonal이 양수인 Q를 선택한다. Original column norm에 대한 |Rjj|/||input_col_j||2<=1e-12 또는 zero는 initial-space rank failure다. Post-QR ||Q^TQ-I||F/sqrt(m0)<=1e-12를 요구한다. 이것은 모든 eigenspace overlap의 보장이 아니며 C5 count와 C3 Gram/residual이 계속 필요하다.

Backend call 동안 oneMKL thread-local count=1을 설정하고 모든 exit에서 이전 값을 복원한다. 공식 mkl_set_num_threads_local의 반환 previous setting을 사용한다. 전역 설정을 바꾸거나 vendor workspace를 동시 공유하지 않는다. 설치 runtime/ABI/유효 thread count·seed·controls를 기록한다. Element contribution 병렬 생성은 허용하되 stable reduction 순서는 유지한다.

Mode cap을 임의의 상수로 만들지 않는다. 정량 representation 경계는 1<=N<=n<=MKL_INT_MAX, nnz+1<=MKL_INT_MAX, 모든 row/column in-range, n+1, 6*nodes, n*m0, n*n, 각 count*sizeof(T)가 internal index 및 min(SIZE_MAX,PTRDIFF_MAX) 안에 있고 해당 container max_size를 넘지 않는 것이다. 곱/합은 overflow 전에 검사하며 double의 정수 근삿값으로 capacity를 계산하지 않는다. Conservative double n×n scratch는 Win64에서 표현상 n≤1073741823이고 실제 사용 가능 메모리 한계는 훨씬 작다. Allocation/factorization 실패는 명시적으로 종료하며 N을 줄이지 않는다. Reference n=3125의 count x만 78,125,000 bytes다. 이를 전체 peak memory라고 부르지 않고 CSR, FEAST reduced storage 및 complex shifted sparse fill을 따로 기록한다.

C7. 승인된 정적 회귀 기준

새 scale을 적용하되 기존 linear NR의 symmetry/frame/energy 1e-12, rigid/full residual/global equilibrium 1e-10을 유지한다. NR-MFE-001 static response는 C2 congruence 1e-12와 complete full residual 1e-10이다. Static moment admissibility는 기존 |d·m|/||m||<=1e-12, zero moment는 별도 branch다.

기존 nonlinear NR §6는 그대로다: rotation/frame 및 no-applied-moment tangent symmetry 1e-12, directional tangent error 1e-6(zero analytical direction 1e-10), finite-rigid objectivity 1e-10, factor-one global force/moment balance 각 1e-8. 기존 finite-difference schedule [1e-3,3e-4,1e-4,3e-5,1e-5,3e-6,1e-6,3e-7] 및 interior index3 기준을 유지한다. 근거는 등록된 tangent test와 기존 승인 numerical contract다. 새 scale의 Dkd=D²kd=0 및 gamma derivative를 검산한다. Nonlinear Hessian 전체의 SPD 또는 linear Tdrill과의 entry equality를 추가하지 않는다. NSNR-069/070의 네 채널·typed zero-scale·100 iteration·no cutback/retry 종료 정책도 기존 계약을 보존한다. 새 eigen C1 zero branch로 해당 승인 typed absolute branch를 덮어쓰지 않는다.