Files
FESADev/docs/hex8-structural-analysis/formulation.md
T

8.8 KiB
Raw Blame History

Hex8 정식화와 검증 불변량

적용: PRD의 FEM-REQ-HEX8-002–007,012–018,021. 물리 선택은 ADR에 있다. 이 문서는 승인된 식의 구현용 유도이며 외부 결과에 맞추어 physics를 바꾸는 권한이 아니다. 수치 controls는 별도 원본, 출력 tensor/identity는 I/O를 따른다.

Shape·기준 geometry·적분

절점 signs는 ---,+--,++-,-+-,--+,+-+,+++,-++다. Point signs는 ---,+--,-+-,++-,--+,+-+,-++,+++, coordinate는 각 ±1/√3이고 weight1이다.

[ N_a={1\over8}(1+\xi_a\xi)(1+\eta_a\eta)(1+\zeta_a\zeta),\quad J_{0,iJ}=\sum_aX_{ai}N_{a,J},\quad G_a=J_0^{-T}\nabla_\xi N_a. ]

J0 column은 자연 ξ/η/ζ 미분이다. Gauss reference 위치는 ΣNaXa, volume weight는 det J0다. Translation 차이를 먼저 사용해 large coordinate offset의 취소 오차를 줄일 수 있지만 기준 좌표와 source identity는 변경하지 않는다. Shape의 partition of unity·linear reproduction와 gradient 합0을 검증한다. Invalid det/inverse는 abs/clamp나 corner 재배열로 숨기지 않는다.

Cell 전체의 Jacobian 인증

Trilinear map의 det J는 각 자연 좌표에 대해 최대2차인 tensor polynomial이다. s=(ξ+1)/2,t=(η+1)/2,r=(ζ+1)/2에서 degree(2,2,2) Bernstein basis로 표현할 수 있다. 각 축 -1/0/+1의27개 determinant 값에서 다음1D 변환을 세 축에 적용한다.

[ b_0=p(-1),\quad b_2=p(1),\quad b_1=2p(0)-(p(-1)+p(1))/2. ]

Bernstein basis는 nonnegative이고 합1이므로 모든 coefficient의 양의 lower bound는 cell 전체 det>0를 보증한다. 모든 sample가 양수인 것만으로는 충분하지 않다. Coefficient bound가 불확실하면 de Casteljau midpoint subdivision 후 child cell을 검사한다. Sample의 nonpositive upper bound 또는 완전히 nonpositive cell은 invalid이며, depth/resource 한계에서도 positive를 인증하지 못하면 uncertain geometry failure다.

인증용 연산은 outward-rounded interval arithmetic으로 수행한다. Addition/subtraction/ multiplication·변환·subdivision·current Xa+ua까지 실제 FP 오차를 포함하며 NaN/Inf/overflow를 거부한다. 단순 coefficient estimate에 임의 epsilon을 더해 positive라고 선언하지 않는다. Current geometry에도 같은 인증을 적용하고 point det F>0를 확인한다. Trial failure가 reference cache나 accepted state를 변경하지 않는다. 세부 종료 controls는 numeric-controls를 따른다.

Linear elasticity와 K0

(\mu=E/[2(1+\nu)],\lambda=E\nu/[(1+\nu)(1-2\nu)]). (E>0,-1<\nu<.5)에서 bulk/shear modulus가 양수다. Engineering strain vector는 [ε11,ε22,ε33,2ε12,2ε13,2ε23]; stress vector는 실제 shear stress다. Normal C block은 diagonal λ+2μ/off-diagonal λ, shear diagonal μ다.

[ K_0=\sum_g B_g^TCB_gw_g,\qquad f=K_0u,\qquad U={1\over2}u^TK_0u. ]

동일 point의 block 형태는 (\lambda G_{ai}G_{bk}+\mu\delta_{ik}G_a\cdot G_b+ \mu G_{ak}G_{bi})다. B/Voigt 코드와 독립 tensor/polynomial oracle를 구분한다. Unconstrained K0는 symmetry/PSD 및6 rigid modes를 가지며 BC 적용 후 mechanism/singularity를 solver가 명시적으로 처리한다. Full integration 잠김을 안정화로 바꾸지 않는다.

Total Lagrangian energy·내력·접선

[ F=I+\sum_a u_a\otimes G_a,\quad E_G=(F^TF-I)/2, \quad W={\lambda\over2}(\operatorname{tr}E_G)^2+\mu E_G:E_G. ]

(S=\partial W/\partial E_G=\lambda\operatorname{tr}(E_G)I+2\mu E_G), (P=FS). (\delta E_G=\operatorname{sym}(F^T\delta F))이므로 (\delta W=P:\delta F), (\delta F=\sum\delta u_a\otimes G_a)다.

[ U=\sum W_gw_g,\quad f_{ai}=\sum P_{iJ}G_{aJ}w_g. ]

(\delta P=\delta F S+F\delta S)를 분리하면

[ \mathbb C_{IJLN}=\lambda\delta_{IJ}\delta_{LN} +\mu(\delta_{IL}\delta_{JN}+\delta_{IN}\delta_{JL}), ] [ (K_{geo}){ai,bk}=\sum\delta{ik}G_{aJ}S_{LJ}G_{bL}w_g, \quad (K_{mat}){ai,bk}=\sum G{aJ}F_{iI}\mathbb C_{IJLN}F_{kN}G_{bL}w_g. ]

Kt=Kmat+Kgeo는 energy Hessian이며 load-independent CLOAD tangent는0이다. S는 symmetric이고 total Kt도 symmetric이나 stressed trial에서 SPD라고 가정하지 않는다. Zero state에서 F=I/S=0/geo=0/total=K0다. 유한 rigid R에서 F=R, Green/W/f는0이며 strained 상태를 global R로 회전하면 energy invariant·force/tangent covariance를 만족해야 한다. Affine deformation의 energy/stress 및 multi-element constant-stress patch가 독립 oracle다.

공식 비선형 stress는 (\sigma=FSF^T/\det F), Green strain은 기준 global tensor다. Linear 공식 strain은 infinitesimal tensor다. 저장 strain의 shear는 실제 tensor 성분이고 B 계산용 engineering shear의 절반이다.8점 값을 그대로 저장한다.

Consistent mass와 generalized eigen

[ M_{ai,bj}=\delta_{ij}\sum_g\rho_0N_aN_bw_g, \quad T={1\over2}v^TMv. ]

M은 reference-state 상수다. Positive weights와8점 interpolation matrix의 full rank에서 SPD가 따른다. Axis-aligned brick의 scalar block은 (\rho V/216)에 각 축 sign이 같으면2/다르면1인 세 factor의 곱이다. Direction별 총합은 rho V, row sum은 rho V/8다. Production에서는 lump하지 않는다.

8점 shape evaluation의 minimum singular value 제곱은1/27이다. Exact local quadrature에서 (m_{e,lower}=\rho\min_gw_g/27)는 local M의 positive lower bound다. Global scatter·principal free restriction에도 active-node incidence를 고려한 lower bound가 보존된다. 예를 들어 모든 active node에 적어도 하나의 element가 있는 경우 min_e bound는 conservative global lower bound다. 실제 FP matrix의 roundoff와 SPD evidence를 포함해야 하며 이 수학적 식만으로 설치 backend의 실행을 검증했다고 주장하지 않는다.

(K_{ff}\phi=\lambda M_{ff}\phi)를 풀며 constrained component0으로 복구하고 (\phi^TM\phi=1), (f=\sqrt\lambda/(2\pi))를 구분한다. Proven mass lower bound가 있다면 (\lambda_{max}\le|K_{ff}|\infty/m{lower})는 search upper-bound 근거가 된다. Diagonal ratio classification scale 자체는 upper bound가 아니다. 검색/count의 실제 증거와 모든 backend pair의 original K/M 잔차·mass norm/Gram을 검증한 뒤 최저 N을 공개한다. Repeated modes는 cluster/subspace로 판단하며 arbitrary rotation의 개별 vector를 외부 비교하지 않는다.

Static partition과 accepted 상태

Linear은 Kff uf=Ff−Kfc uc다. Nonlinear free residual은 fint−F이고 Ktff correction을 풀되 prescribed 값을 factor별 absolute target으로 복구한다. 외력과 prescribed target을 같은 factor로 증가시키며 이전 accepted state로부터 trial을 시작한다. Empty free solve는 유효한 recovery 경로다. Constrained reaction/full residual을 구분하고 모든 required recovery가 끝나야 commit한다.

HHT constant-mass 식

Public alpha=−.05, beta=(1−alpha)²/4=.275625, gamma=.5−alpha=.55다. Algorithmic A와 physical endpoint a를 구분한다.

[ u_p=u_n+hv_n+h^2(.5-\beta)A_n,\quad v_p=v_n+h(1-\gamma)A_n, ] [ A_{n+1}=(u_{n+1}-u_p)/(\beta h^2),\quad v_{n+1}=v_p+\gamma hA_{n+1},\quad a_{n+1}=(A_{n+1}+\alpha a_n)/(1+\alpha). ] [ R=MA_{n+1}+(1+\alpha)(f_{n+1}-F_{n+1})-\alpha(f_n-F_n), \quad K_{eff}=M/(\beta h^2)+(1+\alpha)K_t. ]

이전 free endpoint가 (Ma_n+f_n-F_n=0)이면 R=(1+alpha)(Ma_next+f_next−F_next)다. Endpoint 표현을 쓰면 tangent도1/(1+alpha)로 함께 환산해야 한다. Incomplete previous endpoint를 동등하다고 가정하거나 weighted force에서 reaction을 계산하지 않는다. Physical reaction은 Ma+f−F의 constrained projection이다. Initial u/v0에서 Mff a0=Ff(0+)−ff(0), A0=a0다. 곡선은 endpoints에서 실제 평가하며 prescribed dynamic q/v/a/A는 exact0다.

Physical momenta는 interpolated current x/v를 사용해 (p=\int\rho v,dV_0), (L=\int x\times\rho v,dV_0)다.8점 consistent quadrature와 nodal mass 표현이 일치해야 한다. Angular-momentum rate는 (\int x\times\rho a,dV_0)이며 v×v 항은0이다. External/support moment는 같은 current 위치로 구한다. HHT 및 force trapezoid work는 frame별 exact energy conservation을 의미하지 않는다. Wext는 endpoint force·du trapezoid 누적, error=T+U−T0−U0−Wext다. Alpha0 linear unforced benchmark에서는 독립 보존을 확인한다.

정식화 검증 경계

Analytic brick/affine energy, independent tensor polynomial integrals,6 rigid modes, scaling/ covariance, gradient/tangent directional step sweep와 small-load limit를 결합한다. Warped geometry는 independent higher-order quadrature로 discretization error를 보고한다. Generalized eigen은 작은 independent dense oracle와 original sparse-pencil metrics로 검증한다. HHT는 analytic oscillator/비선형 constant-mass 제조해에서 u/v/physical a·h/h2/h4 수렴을 검증한다. 외부 comparison의5%/1%로 이 내부 검사를 대신하지 않는다.