FESADev/docs/mitc4-shell-feast-eigen/formulation-eigen.md

제약 후 질량 공간, 유한 spectrum과 FEAST

• Feature ID: mitc4-shell-feast-eigen
• Owner: formulation-agent
• Input contract revision: r2
• 기본 보고서: formulation.md

질량 유도 M1~M10, research evidence의 S03~S08, P01, R03/R04/R07를 적용한다. 아래 reduction은 r2 scale로 개선된 원형상 stabilized stiffness와 derived lumped mass의 pencil에 정확히 대응한다. 기존 nonlinear trial tangent와는 별개다.

E1. 조립, 단위 scaling과 실제 global 제약

Element gather A_e를 stable source node/component 순서로 정의하면 K=sum A_e^T K_e A_e, M=sum A_e^T M_L,e A_e이다. 이하 M은 항상 lumped이다. Element contribution→stable COO key→고정 순서 reduction으로 두 행렬을 조립한다. Lumped nodal block이므로 M의 node 간 block은 0이다. Worker는 global storage를 동시에 갱신하지 않는다. Source identity와 equation index는 다른 번호다.

혼합 단위를 가진 6 DOF를 그대로 Euclidean norm으로 합하지 않는다. 활성 shell의 기존 최대 edge 길이들 중 최대값 ell>0 [L]을 characteristic length로 선택하고 D=blockdiag(ell I3,I3)를 node마다 둔다. q=D y에서 y는 무차원, Khat=D^T K D [ML²/T²], Mhat=D^T M D [ML²]이다. 이는 congruence이며 physics와 lambda를 변경하지 않는다.

원 global homogeneous BC의 free injection을 C라 한다. 고정 DOF의 행은 0, free DOF의 행은 하나의 1이며 free 순서는 DofManager의 stable 순서다.

[ y=C y_f,\qquad K_f=C^TKhat C,\qquad M_f=C^TMhat C. \tag{E1} ]

Nonzero prescribed base state는 이 procedure의 입력 오류다. 먼저 local 회전을 임의 삭제한 뒤 global BC를 적용하지 않는다. 특히 global UR3와 drilling은 같지 않다.

E2. Null/positive basis를 node마다 결정

각 active node의 원 질량 null은 M8로부터 z_i=[0;d_i] 한 방향이다. 따라서 node별 실제 free 좌표에서 null이 존재하는 필요충분조건은 constrained rotational component마다 해당 d_i 성분이 정확히 0인 것이다. 그때 Z_i=C_i^T z_i는 unit null vector이며 otherwise null dimension은 0이다. Translation constraint는 이 판정에 영향을 주지 않는다. 모든 회전이 고정이면 null도 없다.

이 판정은 analytic nullspace와 실제 입력 BC의 교집합이다. 작은 nonzero constrained director 성분을 임계값으로 0으로 바꾸지 않는다. 거의 정렬된 경우 작은 양의 physical mass가 남을 수 있으며 이를 보존한다. Numerical conditioning 검증을 통과하지 못하면 명시적 실패로 처리하고 free DOF를 추가 삭제하지 않는다.

각 node에서 Z_i의 Euclidean 직교 여공간 basis를 P_i^f라 한다. 회전 3개가 free이면 기존 tangent basis (a_i,b_i)를 사용할 수 있다. Free 회전 2개 안에 null이 있으면 그 평면에서 d_i에 수직인 단위벡터 하나, 회전 1개가 null이면 positive 회전 basis는 없다. Free translation axis는 그대로 유지한다. Null이 없으면 P_i^f=I다. Source component 순서, basis sign/tie는 고정한다. Householder 등 동등한 작은 orthogonal construction을 쓸 수 있으나 physical rank를 수치 eigenvalue cutoff로 정하지 않는다.

[ P=\operatorname{blockdiag}(P_i^f),\quad Z=\operatorname{blockdiag}(Z_i),\quad [P\ Z]^T[P\ Z]=I,\quad M_p=P^TM_fP\succ0,\quad M_fZ=0. \tag{E2} ]

여기서 P_i^f와 M2의 tangent projector는 서로 다른 기호다. 실제 positive dimension n=cols(P), null dimension nz=cols(Z)이다. Reference의 full clamped flat mesh에서는 625개 interior node, nf=3750, nz=625, n=3125가 예상된다. 실제 assembly와 rank evidence로 확인한다.

E3. 개선된 stiffness에서 decoupling의 정확한 증명

공통 node director는 shell geometry가 한 번 정의하고 ElementFactory::CreateShell이 node identity로 각 element에 전달한다. MITC4 kernel의 PhysicalTransformation20, DrillingTransformation4, Stiffness에서 M4의 구조를 확인했다. 옛 계수만 S1/S2로 대체하며 모든 incident element의 k_d,e는 finite positive다. Physical kernel과 공통 director·nodal drill map은 그대로다.

한 node에만 [0;d_i]가 있는 full vector z_i에 대해 T_p A_e z_i=0이다. Drilling 부분은 해당 node에만 작용하므로

[ Khat z_i=\kappa_i z_i,\qquad \kappa_i=\sum_{e\ni i}k_{d,e}>0,\qquad K_f Z=Z\operatorname{diag}(\kappa_i). \tag{E3} ]

마지막 등식은 E2에서 z_i가 실제로 free space에 있을 때 성립한다. 따라서

[ K_{pz}=P^TK_fZ=0,\quad K_{zz}=Z^TK_fZ=\operatorname{diag}(\kappa_i)\succ0, \quad K_p=P^TK_fP,\qquad K_p p=\lambda M_p p. \tag{E4} ]

Finite mode의 massless 성분은 z=0이다. 이 projection은 arbitrary small coupling을 버리는 근사가 아니라 승인된 physical/drilling separation의 항등식이다. Null이 없는 부분 제약 node에서는 P=I이고 개선된 drilling contribution이 K_p에 그대로 남는다. 이를 빼면 원 constrained stabilized pencil을 바꾸므로 금지한다. 여기서 원 pencil은 r2의 실제 BC 적용 후 K/M이다. 옛 K와 새 K의 spectrum 동일성은 요구하지 않는다. S3의 complete K 공변성은 이 partial-BC contribution까지 포함한다.

실제 floating-point matrix와 basis의 M_f Z, K_fZ-Z diag(kappa), P^T Z를 별도 evidence로 남긴다. 증명된 analytic null만 구조적으로 제거하고 계산된 tiny coupling을 근거 없이 threshold로 버리지 않는다. 이 항등식의 roundoff 검증 및 복원 후 full free residual이 NR 기준을 넘으면 assembly/director/conditioning 실패로 중단한다. 추가 massless mode를 찾아 제거하거나 개선된 K를 재조정하지 않는다. S5의 nonlinear Hessian에는 corotation coupling이 있으므로 E3/E4를 모든 nonlinear 상태, preload 또는 prestressed eigenproblem에 적용하지 않는다.

E4. 일반 Schur 식과 비정상 rank의 의미

E3의 구조가 없는 일반 symmetric pencil을 비교하면, Kzz가 가역일 때만

[ z=-K_{zz}^{-1}K_{zp}p,\quad (K_{pp}-K_{pz}K_{zz}^{-1}K_{zp})p=\lambda M_p p. \tag{E5} ]

이는 massless row를 푼 정확한 elimination이며 massive DOF의 근사 Guyan reduction이 아니다. 예를 들어 K=[[6,2],[2,4]], M=diag(2,0)이면 lambda=5/2, z=-p/2이다. 첫 행만 남기는 lambda=3은 잘못이다. E5는 independent small-pencil 검산과 잘못된 projection을 잡는 기준이다. 이 Feature의 production shell 경로에는 E3이 적용되므로 일반 dense Schur 구현을 요구하지 않는다. E3과 불일치하는 input을 E5로 자동 fallback하여 승인 geometry 문제를 숨기지 않는다.

• nz=0: null 처리 없이 실제 constrained Kf,Mf의 문제다.
• Local positive rank 0: 모든 free DOF가 null이면 해당 node의 finite 좌표는 없고 양의 Kzz에 의해 그 free coordinate는 0으로 복원한다. Fully constrained node도 finite 좌표 0개다. 다른 node의 유효 finite 좌표는 계속 보존한다.
• Global n=0 또는 N>n: 요청 모드가 존재할 수 없으므로 명시적 model/solve 실패다. 정적 해석의 유효 0×0 Kff 계약과 이 양의 mode-count 요구를 혼동하지 않는다.
• 예상 Kzz가 singular/nonpositive이거나 활성 node에 mass contribution이 없음: 유효 shell의 E3 전제가 깨진 것이다. Common K/M null이면 임의 lambda에 대한 indeterminate direction이며 ordinary finite eigenmode가 아니다. 제거하거나 pseudoinverse, stiffness/mass epsilon으로 정규화하지 않고 진단한다.
• 일반 indefinite Kzz의 singular/nonzero Kpz case는 descriptor-pencil의 추가 constraint를 요구한다. 이번 positive stabilized shell과 다른 문제이며 자동 해석하지 않는다.

E5. SPD whitening, 복원과 sparsity

M_p는 node별 최대 6×6 SPD block이다. 각 block을 Cholesky하여 M_p=L L^T를 얻고 다음을 정의한다.

[ A=L^{-1}K_pL^{-T},\quad p=L^{-T}x,\quad A x=\lambda I x,\quad \phi=D C P L^{-T}x. \tag{E6} ]

Inverse를 dense global matrix로 만들지 않고 작은 block triangular solve로 congruence를 구성한다. phi는 source node의 global 6 DOF로 복원되고 constrained entries는 정확히 0이다. phi^T M phi=x^T x다. Identity mass를 전달하는 E6도 실제 lumped mass에서 얻은 generalized SPD 문제다. Consistent mass 또는 인공 unit physical mass로 바꾼 것이 아니다.

Projection/whitening은 node 간 topology를 늘리지 않는다. 기존 node coupling마다 작은 dense block만 생기며 전체 O(n²) Schur matrix와 dense inverse를 피한다. Sparse shifted FEAST factorization 자체의 fill은 별도 존재한다. Full M_C는 element별 24×24와 세 4×4 moment로 검증할 수 있고 global dense consistent matrix는 불필요하다. Allocation/nnz/factorization memory 실패는 불완전 모드의 성공이 아니라 명시적 failure다.

E6. 유한 search interval과 요청 개수

요청 N은 .inp *FREQUENCY의 양의 정수이며 현재 reference는 10이다. Reference Value/Freq/행수는 interval, scale, count 또는 parameter 선정에 사용하지 않는다.

먼저 s=||A||inf [T^-2]를 안정된 합으로 계산한다. s>0이면 B=A/s로 두고, 반환값은 lambda=s mu다. A=0이면 s=min_e E_e/(rho_e ell²)>0라는 모델 유래 time scale을 사용하며 B=0이다. 이 zero-matrix case도 N개의 0 모드를 FEAST에서 계산·검증한다. Invalid/underflow/overflow scale을 임의 1로 바꾸지 않는다.

각 행의 off-diagonal absolute sum을 r_i라 하면 Gershgorin에 의해 min_i(Bii-ri)<=mu<=max_i(Bii+ri)다. 안정된 directed-outward accumulation 또는 동등하게 증명된 roundoff bound를 포함하여 finite enclosure [gL,gU]를 만든다. B=0에서는 알고리즘용 dimensionless enclosure [-1,1]을 사용한다. Nonzero B의 endpoint margin은 NR가 정할 dimensionless coefficient와 ||B||inf 및 endpoint ULP로 결정하고, outward nextafter를 적용해 a0<gL<=gU<b0를 확보한다. 이는 물리 mass/eigenvalue floor나 zero clamp가 아니다.

CSR count c(a,b)는 endpoint가 분리된 interval에 대해 Intel fpm[49]=1이 반환한 count와 status다. c(a0,b0)=n을 먼저 확인한다. 다음은 finite search의 결정적 규칙이다.

1. 모든 count는 원 enclosing lower bound a0부터의 prefix count로 평가한다. lo=a0의 알려진 prefix count는 0, hi=b0의 count는 n이다. Upper search에서 midpoint v의 count가 N 미만이면 lo=v, 아니면 hi=v로 둔다. 성공한 hi와 count를 항상 보존한다. N<=count<=min(n,2N+8)인 upper endpoint를 얻으면 이 최적화를 끝낼 수 있다. 이는 workspace 목표이며 N의 변경이나 mode cap이 아니다.
2. Cluster 때문에 목표 count로 줄일 수 없어도 count>=N인 마지막 safe upper를 쓴다. Representable midpoint가 없거나 NR endpoint separation 폭에 도달하면 subdivision을 끝낸다. 경계에 있는 반복 고유치의 일부만 count하도록 tolerance를 조정하지 않는다.
3. 불필요하게 먼 lower bound를 줄이기 위해 a=a0, 위에서 얻은 b=hi로 시작한다. v=(a+b)/2의 prefix count가 0이면 a=v로 이동하고, 양수이면 이 하한 축소를 끝낸다. 따라서 (a0,a)의 empty-prefix evidence를 잃지 않는다. 아주 작은 고립 cluster에서 넓은 구간이 검색 실패를 유발하는 Intel의 경고를 이 과정과 후검증으로 다룬다.
4. 최종 (a,b)에서 다시 count하여 c>=N이고 알려진 prefix와 일치함을 확인한다. Capacity는 m0=min(n,max(c+8,ceil(3c/2)))로 한다. c=n이면 m0=n이다. Checked integer arithmetic을 쓰며 큰 multiplicity도 전체 count를 수용한다.
5. Solve 성공 후 모든 c개를 독립 검증하고 정렬하여 첫 N개를 결과 candidate에 넣는다. 같은 eigenvalue의 multiplicity를 지우거나 평균내지 않는다. N번째가 repeated/cluster의 내부에 있으면 내부 계산 group은 보존하고 출력 row는 정확히 N개다.

Midpoint·endpoint 이동은 overflow-safe하게 수행한다. NR는 count-call budget, endpoint separation coefficient와 perturbation 횟수를 구현 전에 고정한다. Budget은 각 search 및 endpoint 재시도까지 포함하는 유한 정수다. -5 등의 endpoint 실패 시 기록된 bracket 안에서 NR가 정한 간격으로 양옆의 endpoint 후보를 순서대로 평가하며, 성공한 prefix/count evidence 없이는 bracket을 갱신하지 않는다. 마지막 valid bracket이 없거나 count monotonicity/total count가 맞지 않으면 실패다. Search 최적화 budget에 도달하면 이미 검증된 enclosing bracket을 쓸 수 있으나 실패한 count를 성공으로 바꾸지 않는다.

Postsolve에서는 returned m=c, finite paired vectors/values, Gram rank, interval 내부 및 endpoints와의 separation을 확인한다. Endpoint에 너무 가깝거나 info=1/2/3, count/solve mismatch이면 NR가 정한 유한 회수 안에서 safe bracket 확대·재count·capacity 확대만 허용한다. 매번 count/solve history를 보존한다. 그 뒤에도 불일치하면 실패하고 부분 모드 또는 다른 eigensolver로 fallback하지 않는다. Formulation은 이 수학적 알고리즘을 NR에 인계한다. 구현 진입 bundle은 다음 stage가 수치값·budget을 고정한 뒤 완성된다.

Full enclosing bound, empty-prefix, interval count와 solve-count/Gram 검증은 낮은 N개의 누락을 탐지하는 서로 다른 evidence다. Intel의 documented total count를 모든 conditioning, endpoint 및 floating-point 상황의 completeness theorem으로 확대하지 않는다. 독립 known-spectrum의 zero/negative/cluster/repeated/size1/request=n case와 실제 설치 runtime 검증이 필수 연결점이며 선택 및 수치 기준은 NR가 소유한다.

E7. FEAST 경계와 deterministic 실행

Actual eigenpairs는 dfeast_scsrgv로 E6의 (A,I) 또는 E6에서 동일하게 scale한 (B,I)를 풀어 구한다. Real symmetric, SPD second matrix, n>0, 0<m0<=n, strict finite emin<emax를 만족한다. uplo=U로 A와 identity mass의 각 CSR을 전달하며 row/column 증가 순서, duplicate-free, 1-based row pointer, terminal nnz+1을 검증한다. 두 pattern이 같을 필요는 없다. MKL_INT, n, nnz+1, column index와 n*m0*sizeof(double)의 변환/곱은 checked다. FEAST vendor type과 scratch resource는 private backend에만 있다.

Call마다 feastinit한 128-entry parameter 배열에서 이번 call의 설정을 채운다. CSR check, SPD check를 켜고 subspace-only mode는 끈다. Zero mode에도 적용할 수 있도록 trace stopping을 사용하고 독립 residual을 반드시 후검증한다. epsout만으로 pass하지 않는다. Contour points, stopping exponent, finite loop bound는 NR가 정한다. User initial subspace를 선택하고 고정 seed/명시된 integer generator 및 고정 QR 순서로 full-column-rank 초기 공간을 만든다. m0=n에서는 identity를 쓴다. Initial guess가 모든 eigenspace와 충분히 overlap한다는 보장은 없으므로 count/solve/Gram 검증으로 부족한 projection을 검출한다. Seed를 바꾸며 통과할 때까지 반복하지 않는다. Backend의 thread/재현 설정과 seed는 실행 evidence에 기록한다. 동일 source 입력의 sparse structure, reduction, row 및 diagnostic 순서는 고정하며 repeated mode의 개별 벡터 방향은 수학적으로 유일하지 않다는 제한을 유지한다.

Solve는 e,res의 m0 entries와 column-major x[n*m0]를 준비하며 첫 m columns만 소비한다. Count call의 m0는 output count로 overwrite되므로 solve 요청수/실제 allocation capacity와 분리한다. 설치 API가 count-only workspace의 작은 용량을 명시하지 않은 경우 보수적으로 count scratch는 m0=n, x[n*n]까지 준비하고 count 단계 사이에 재사용한다. Solve 전 해제할 수 있다. Reference n=3125에서 이 scratch는 약 78.1 MB (74.5 MiB)이고 10모드의 실제 solve subspace는 이보다 작다. Count용 보수적 allocation은 dense Schur와 별개지만 O(n²) memory 비용이므로 Planning/Implementation이 peak 및 shifted sparse factorization fill을 측정한다. Allocation 실패 시 N을 줄이지 않는다.

Intel 근거: API/storage/capacity, C-index parameter와 count, status/endpoint/고립 cluster 주의. 2026-09-06 확인. info=0 이외의 빈 interval, 미수렴, subspace 부족, count/inner solver, allocation 또는 SPD 실패는 status와 stage를 보존하여 진단한다. Count success는 eigenpair success가 아니며, caller의 count authority는 끝까지 .inp다.

E8. Mode identity, normalization과 상태

모든 returned eigenpair를 lambda 오름차순으로 함께 정렬한다. Equal lambda의 tie는 original solver column order로 안정 정렬한다. Near-equal 값을 평균하거나 merge하지 않는다. 양의 modal mass로 column별 normalize하고, 가장 큰 absolute scaled global component의 sign이 양수가 되도록 한다(source DOF 순서로 tie 결정). 서로 다른 eigenvalues의 벡터를 Gram 보정 목적으로 임의 혼합하지 않는다. Failed orthogonality는 solve/validation 실패다. 정확한 repeated eigenspace의 다른 orthonormal basis는 같은 해이고 mode-shape external matching은 없다. Density scaling/좌표 테스트는 repeated group의 projector로 비교한다.

Normalized phi^T M phi=1의 수치 convention에서 phi_u의 차원은 M^-1/2, phi_theta는 M^-1/2 L^-1이고 modal amplitude가 M^1/2 L를 가진다. Mode vector는 실제 응답 변위가 아니다. Scalar lambda [T^-2]와 source step/mode identity를 분리하고 internal node/equation number를 source label로 노출하지 않는다.

원형상·density 검증→consistent/moment 검증→lumped 검증→K/M 조립→실제 BC→E2E6 admissibility→E6E7 FEAST→global 복원/독립 검증→candidate 결과 순서다. 반환 lambda는 부호를 보존한다. NR 기준상 유의한 음수는 불안정/invalid spectrum 실패, near-zero는 분류와 raw value를 함께 유지하며 0으로 clamp하지 않는다. 양수에서만 f=sqrt(lambda)/(2pi)를 파생한다. 정확한 0은 f=0이고 작은 음수에는 실수 frequency를 만들지 않는다. 이 validity의 저장 방식은 I/O owner가 정한다.

N개 검증이 끝나기 전에 accepted state를 바꾸지 않는다. 검증된 candidate HDF5를 완성·검증·close한 뒤 final replacement를 수행하고 성공 뒤 준비된 state를 교환한다. Replacement 이후 state 교환은 실패하지 않도록 준비한다. 제어된 solve/validation/write/ finalization failure에서 이전 accepted state와 기존 final HDF5를 보존한다. Crash/power loss까지 원자성을 새로 주장하지 않는다. Live log lifecycle은 별도다. Modal residual의 constrained 성분을 static reaction으로 출력하지 않는다.