Files
FESADev/docs/hht-implicit-dynamics/formulation-historical-docs-cleanup-20261004.md
T

16 KiB
Raw Blame History

HHT implicit dynamics 정식화

  • 상태: independently-reviewed-candidate; physics-and-validation-pending.
  • Owner: formulation-agent, 2026-09-22. 문서 gate 산출물이며 구현 승인이나 수치 검증 통과가 아니다.
  • 기준: 프로젝트 지침, 제품 PRD, 공통 ADR, Feature PRD, 요구사항, 연구.
  • 관성의 전개식·초기 일관성·정확도 근거는 formulation-inertia가 소유한다. 영속 schema, C++ 저장 구조와 tolerance는 정의하지 않는다.

1. 결론과 적용 범위

추천 후보는 유한회전 kinetic energy에서 유도한 관성력과, 물리적 가속도 a / 알고리즘 가속도 A를 분리하는 HHT 확장이다. 일정 질량의 표준 HHT force weighting으로 정확히 환원한다. 단순히 M(q)를 기존 일정 질량 잔차에 대입하는 식은 추천하지 않는다. 이 차이는 상세 §5에 유도한다.

사용자가 shell 속도 보간의 운동에너지를 그대로 적분하는 consistent mass를 선택했다. 따라서 HHT 기본값과 이번 구현 범위는 consistent mass다. 같은 운동학의 scalar-moment nodal lumped 식은 비교용 이론 및 기존 eigen 회귀 참고이며 이번 HHT 구현 요구가 아니다. 선행 eigen 결과나 동적 reference와의 적합도로 선택한 것이 아니다. 유한회전 지원 범위를 작은 회전으로 축소하지 않는다.

아래 HHT 확장과 0+ 초기화는 FESA 채택 제안이다. 이번 consistent mass 선택이 다른 physics·알고리즘·공개 output 계약 전체의 승인을 의미하지 않는다.

구분 정식화 경계
물리 초기 형상 기준 Total Lagrangian FESA-MITC4, 큰 변위·유한회전, 작은 탄성 변형률, 균질 등방 선형 탄성, 고정 두께
시간 단일 step, DIRECT, h=0.1, T=10, t_n=n h, 초기 frame 및 100 increments. INC=1000은 increment 상한
하중/구속 t=0+부터 node 2의 global 1–3 force 각각 −100000; 외곽 지정 DOF 1–6은 모든 시점에 0
초기값 q_0=v_0=0; a_0는 동적 초기 평형·대수 일관성으로 구하며 A_0=a_0
제외 물리적 damping, contact, 소성·손상, follower/분포하중, nodal moment 입력 확대, 임의 amplitude, 시간변화 구속, adaptive/cutback, multi-step/restart, 다른 요소 dynamics

원본 INP의 section-point 입력 5를 새 질량 또는 내력 quadrature 지정으로 해석하지 않는다. Abaqus 내부 요소·질량·회전·초기화와의 동등성은 목표가 아니다.

2. 좌표·운동학·재료와 내력

절점 좌표 순서는 정확히

[ q_i=[u_x,u_y,u_z,\phi_X,\phi_Y,\phi_Z]_i^T,\quad v=\dot q,\quad a=\ddot q. ]

u는 global translation, phi는 초기 형상과 고정 global basis에 대한 total rotation vector다. 각 절점에서

[ R_i=\exp([\phi_i]_\times),\quad Q_i=R_iQ_i^0,\quad d_i=R_i d_i^0, \quad |\phi_i|<\pi . ]

Newton correction은 q←q+Δq이고, orientation은 그 total phi에서 다시 평가한다. Incremental rotation을 별도로 누적하거나 accepted step마다 chart/reference를 재설정하지 않는다. 기존 principal chart, director opposition, swing/twist, 양의 Jacobian 및 polar factor 유효성 경계를 보존한다. 실패 trial을 clamp하거나 재정규화하여 통과시키지 않는다.

[ x(\xi,\eta,\zeta;q)=\sum_i N_i(X_i+u_i+z d_i),\quad z=t\zeta/2, \quad dV_0=J_0,d\xi d\eta d\zeta . ]

N_i는 기존 bilinear MITC4 shape function이다. Interpolated director ΣN_i d_i를 단위화하지 않는다. 정확한 각속도와 각가속도는

[ \omega_i=J_l(\phi_i)v_{\phi i},\qquad \dot\omega_i=J_l(\phi_i)a_{\phi i}+DJ_l(\phi_i)[v_{\phi i}]v_{\phi i}, ] [ \dot d_i=\omega_i\times d_i,\qquad \ddot d_i=\dot\omega_i\times d_i+\omega_i\times(\omega_i\times d_i). ]

따라서 v_phi와 a_phi는 공간 각속도·각가속도가 아니다. A_phi는 이 물리 회복식에 넣지 않는다. J_l과 미분의 구현 가능한 정의는 상세 §1에 있다.

물리 내력은 기존 비선형 정식화 §§6–13를 따른다. Green–Lagrange strain의 MITC assumed shear와 fixed-reference local engineering strain ē=[E11,E22,2E12,2E13,2E23]에 대해

[ U_{phys}=\frac12\int_{\Omega_0}\bar e^TC_5\bar e,dV_0,\quad C_5=\operatorname{diag}(C_{ps},(5/6)G I_2),\quad G=E/[2(1+\nu)], ] [ C_{ps}=\frac{E}{1-\nu^2} \begin{bmatrix}1&\nu&0\\nu&1&0\0&0&(1-\nu)/2\end{bmatrix}. ]

Objective numerical drilling potential은 기존 polar-corotation / swing-twist angle χ_ei와 초기 고정 stiffness k_de를 그대로 사용한다:

[ U_{drill}=\frac12\sum_{e,i}k_{de}\chi_{ei}(q)^2,\quad f_{int}=\nabla_q(U_{phys}+U_{drill}),\quad K_{int}=\nabla_q^2(U_{phys}+U_{drill}). ]

Drilling stiffness는 물리적 회전 관성이 아니다. U_phys와 U_drill을 구별하여 에너지 evidence를 계산하고, drill stress/resultant를 만들지 않는다. 재료 내부변수·소성 update는 없다. 재사용 대상 evaluator의 PhysicalEnergy, DrillingEnergy, EvaluateInternalResidualAndTangent도 이 구분에 대조한다. 현재 공개 equilibrium/recovery 결과는 physical_strain_energy만 노출하므로 drill energy evidence의 전달은 후속 구현 작업이며 현 API가 이를 이미 제공한다고 보지 않는다.

3. 질량·관성·연속 평형

B_i=D_phi d_i=−[d_i]_× J_l(phi_i)이고 C_ij^(k)=∫rho N_iN_j z^k dV_0이다. 채택한 consistent T_C와 비교용 lumped T_L은 다음과 같다.

[ T_C=\frac12\int_{\Omega_0}\rho\left|\sum_iN_i(v_{ui}+zB_iv_{\phi i})\right|^2dV_0, ] [ T_L=\frac12\sum_i\left(m_i|v_{ui}|^2+2s_i v_{ui}\cdot B_iv_{\phi i} +j_i|B_iv_{\phi i}|^2\right),\quad (m_i,s_i,j_i)=\sum_j(C_{ij}^{(0)},C_{ij}^{(1)},C_{ij}^{(2)}). ]

비교용 lumped 식은 scalar moment의 row sum 후 tensor를 재구성하며 global matrix row sum이 아니다. 이번 HHT에는 T=T_C, M=M_C를 사용하고 artificial drilling inertia는 허용하지 않는다. T=½v^T M(q)v의 질량행렬은 회전에 의존한다. 관성 covector는

[ \boxed{I(q,v,a)=\frac{d}{dt}\partial_vT-\partial_qT=M(q)a+c(q,v)}, \qquad c_i=\sum_{j,k}\frac12(M_{ij,k}+M_{ik,j}-M_{jk,i})v_jv_k. ]

고정 구속을 제외한 모든 admissible δq에 대해

[ \delta q^T E=0,\quad E=I+f_{int}-F(t), \qquad F(t)=F^\star\ (t\ge0+). ]

이는 parent rho_0 ẍ=Div_0 P, P=F_def S의 가상일 이산화다. 강형의 F_def는 deformation gradient이며 nodal force F와 구별한다. 물리적 damping은 0이다. Tangent에는 DM·a, convective force의 configuration/velocity 미분을 모두 포함한다. 구체적인 matrix block, quadrature, derivative action과 null 조건은 상세 §§1–4에 고정한다.

4. 초기화와 HHT 시간적분 후보

음수 alpha 관례를 사용한다:

[ -\frac13\le\alpha\le0,\quad \gamma=\frac12-\alpha,\quad \beta=\frac14(1-\alpha)^2. ]

기본값은 alpha=−0.05, beta=0.275625, gamma=0.55다. 공개 기본값 근거와 Sandia의 결합계수를 따른다. alpha=0은 평균가속도 Newmark다. 이 범위는 일반 유한회전 문제의 무조건 안정성이나 큰 시간간격의 정확도를 보증하지 않는다. 별도 beta/gamma 입력을 승인한 의미도 아니다.

FESA 초기화 후보 계약: frame 0은 t=0+, q_0=v_0=0, 전값 F_0=F*다. M_ff(q_0)a_f0=F_f0−f_int,f(q_0)와 mass-null 일관성을 함께 풀어 a_0를 구한다. A_0=a_0, g_0=f_int(q_0)−F_0를 정의한다. g_0=0 또는 a_0=0을 임의 대입하지 않는다. 이는 공개 STEP 의미에 맞춘 FESA 계약 제안이며 Abaqus의 숨은 첫 history 생성 절차를 확인한 주장이 아니다.

Brüls–Arnold (2008), 식 (5)–(8)의 보조가속도 형식을 alpha_m=0, alpha_f=−alpha로 전문화한다. 아래 a는 물리적 좌표 가속도, A는 시간적분 보조변수다:

[ \boxed{A_{n+1}=(1+\alpha)a_{n+1}-\alpha a_n}, ] [ q_{n+1}=q_n+h v_n+h^2[(1/2-\beta)A_n+\beta A_{n+1}], \qquad v_{n+1}=v_n+h[(1-\gamma)A_n+\gamma A_{n+1}], ] [ \boxed{E_f(q_{n+1},v_{n+1},a_{n+1},t_{n+1})=0.} ]

같은 fixed-global chart의 a_n만 history로 사용한다. 정확도 적용 조건과 source의 singular-mass 한계는 상세 §5에 명시한다. Shell drilling DAE까지 증명되었다고 보지 않는다.

Newton 미지수는 전체 좌표의 free 성분 q_f,n+1이다. q_c=v_c=a_c=A_c=0를 항상 유지한다. Predictor와 trial 종속량은

[ q_p=q_n+h v_n+h^2(1/2-\beta)A_n,\quad v_p=v_n+h(1-\gamma)A_n, ] [ A(q)=\frac{q-q_p}{\beta h^2},\quad v(q)=v_p+\frac{\gamma}{\beta h}(q-q_p),\quad a(q)=\frac{A(q)+\alpha a_n}{1+\alpha}. ]

동일 trial에서 E를 평가하고 매 반복 full tangent를 갱신한다:

[ \boxed{K_{eff}p=D_qI[p]+\frac{\gamma}{\beta h}D_vI[p] +\frac{M(q)p}{(1+\alpha)\beta h^2}+K_{int}p},\quad (K_{eff})_{ff}\Delta q_f=-E_f,\quad\Delta q_c=0. ]

여기서 D_qI는 v,a를 고정한 편미분이다. 외력은 고정 global force이므로 D_qF=0이다. 관성 convective 항 때문에 K_eff 대칭을 일반적으로 가정하지 않는다. 모든 DOF가 구속되면 free solve는 비어 있고 full reaction을 회복한다. Singular mass만으로 모델을 거부하거나 epsilon mass를 추가하지 않으며, 유효 접선의 singularity·대수 regularity 실패는 별도로 진단한다.

5. 잔차·반력·에너지의 구별

일정 M에서 g=f_int−F라 두면, 앞선 endpoint physical equilibrium M a_n+g_n=0을 이용하여

[ R_{HHT}:=M A_{n+1}+(1+\alpha)g_{n+1}-\alpha g_n =(1+\alpha)E_{n+1}-\alpha E_n. ]

따라서 free 방정식은 표준 force-weighted HHT와 같다. 표준 HHT의 Newmark 가속도를 이 문서에서는 A라 부른다. 이를 a로 오인하면 출력·반력이 달라진다. 구성 의존 질량에는 위 M A+weighted g 식을 그대로 적용하지 않는다. 추천 후보의 Newton residual은 endpoint E_f이며, HHT filter 결함 A_{n+1}−(1+alpha)a_{n+1}+alpha a_n은 별도 운동학 evidence다.

구속 selector C_u,C_phi에 대한 endpoint 반력은

[ F_i^{react}=C_{u,i}E_{u,i},\quad m_{\phi,i}^{react}=C_{\phi,i}E_{\phi,i},\quad M_i^{react}=J_l(\phi_i)^{-T}m_{\phi,i}^{react}. ]

부호는 structure에 작용하는 support force/moment다. Free E_f는 평형 evidence이고 반력이 아니다. R_HHT,c는 이전 반력이 섞인 알고리즘량이며 endpoint 반력으로 저장하면 안 된다. Moment는 m_phi·δphi=M_spatial·δtheta에 의해 변환한다. Numerical drill의 force/work는 물리적 section stress와 분리한다.

연속계의 에너지 관계는 d(T+U_phys+U_drill)/dt=v^TF다. 일정 dead force와 고정 경계에서는

[ W_{0\to n}=F^{\star T}(q_n-q_0),\quad \Delta\mathcal E_n=T(q_n,v_n)+U_{phys}(q_n)+U_{drill}(q_n)-\mathcal E_0-W_{0\to n}. ]

ΔE_n은 이산 에너지 결함이며 물리 damping energy가 아니다. HHT의 일반 비선형 문제에서 매 frame 단조 감소나 정확한 momentum 보존을 요구하지 않는다. 무감쇠 선형 alpha=0의 보존성, 음수 alpha의 선형 고주파 소산과 비선형 work/energy 수렴은 따로 검증한다.

6. 증분·실패·회복 순서

  1. 승인 subset, source identity, 재료·geometry·quadrature 및 고정 구속을 확인한다. Consistent mass의 moment·rank/positivity 조건과 초기 대수 일관성을 검사한다.
  2. 초기 상태 및 F_0, a_0, A_0를 완성하여 frame 0을 검증한다. 초기화 실패를 zero acceleration으로 대체하지 않는다.
  3. t_n=n h의 직전 committed q,v,a,A에서 predictor를 만들고, 고정 구속을 적용한다. INC와 정확한 종료 increment를 검사한다.
  4. 각 Newton trial에서 total rotation, 내력·tangent, 관성·tangent와 E를 같은 state로 계산하고 free correction을 푼다. 병진/회전 correction 및 force/moment residual을 각각 scale-aware로 검사한다. 수치 한계·iteration 상한은 Numerical Review/I/O 소유다.
  5. 수렴 후 같은 trial에서 physical acceleration, triad/director, endpoint reaction, 물리적 strain/stress/resultant, 에너지 및 null-constraint evidence를 회복·검증한다. 고정 reference frame과 element Gauss/section-point identity는 기존 정적 회복을 보존하며 nodal averaging은 하지 않는다.
  6. time,increment,q,v,a,A,orientation,forces,recovery,evidence가 한 identity를 이루는 candidate만 원자적으로 commit한다. 모든 101 accepted frames를 보존한다. Newton iteration은 time frame이 아니다.
  7. 비수렴·invalid chart/geometry·nonfinite·선형 solve·회복·출력 실패는 trial을 폐기하고 마지막 committed 상태/history와 기존 최종 HDF5를 보존한다. 고정 DIRECT를 cutback/adaptive로 바꾸거나 실패 frame을 성공으로 내보내지 않는다.

공식 출력은 lexical sibling <stem>_results.h5다. U/UR, 물리적 v/a, omega/omega_dot, algorithmic A, 반력과 에너지의 공개 여부 및 영속 schema는 미결정이다. Algorithmic A를 physical acceleration이라는 이름으로 공개해서는 안 된다. 출력 선택과 무관하게 검증에는 서로 구별되는 내부 quantity가 필요하다.

7. 검증 인계와 미해결 결정

항목 필수 evidence / 남은 경계
HHT 계수·initialization 표준 일정 질량 HHT와 algebraic reduction; SDOF 자유진동/STEP onset 및 a_0=F/m; alpha=0 한계
관성·회전 Consistent mass의 총질량, uniform translation, 유한회전 angular motion, d·ḋ=0, `d·d̈+
mass-null 구성 의존 null, ψ=0과 초기 1·2차 일관성; coupled massless DAE에서 q/v/a 시간 수렴, finite-rotation drill coupling. 이 source theorem의 직접 적용은 미확인
정확도 Smooth M(q) manufactured case의 q/v/a를 h,h/2,h/4로 검산; constant-M와 variable-M를 분리. STEP의 시작 불연속과 이후 smooth 구간을 분리
에너지·반력 관성 구속 평형, generalized/spatial moment work, algorithmic/endpoint reaction 구별, T/U_phys/U_drill/W와 momentum evidence
구현 lifecycle 동일 state tangent, 실패 rollback, time/frame 순서, deterministic assembly와 반복 결과
외부 reference 지정 원본 101 frames × 9 nodes, frame/instance/node/component의 exact mapping; U blocking / UR warning-only 및 공통 TOLERANCE. V/A/energy correctness는 해당 CSV에 coverage 없음

수식 후보의 독립 검토는 완료되었으며 표시된 수식에서 새 오류는 발견되지 않았다. 질량 선택은 consistent로 확정되었다. 구현 전 남은 계약은 (1) HHT 확장·초기화의 FESA 채택과 나머지 physics 범위, (2) initial algebraic regularity와 trial 실패 판정의 수치 기준, (3) 공식 output 의미의 공개 범위/I/O다. 이 선택 후 구현할 수 있는 수식 후보이며, 실제 유한회전 shell의 drilling DAE·관성 tangent·q/v/a 시간 수렴은 구현 후 기능 acceptance에 필요한 검증이다. 부모 문헌의 nonsingular-mass 조건을 shell mass-null에 자동 적용하지 않되, 모든 shell DAE의 보편적 수렴 증명을 코드 작성의 선행 조건으로 요구하지 않는다. 현재 문서 완료나 제한된 probe를 전체 solver 승인·검증 통과로 표시하지 않는다.

이번 owner 작업은 기존 정식화·source 대조, 일차 문헌 대조와 수식 유도 및 문서 정적 검토다. 독립 시간정확도 검산은 variable-mass SDOF와 moving-null 예제에서 추천 후보의 약 2차 q/v/physical-a 수렴을 관찰했고 naive 확대의 약 1차 오차를 구별했다. 이는 shell solver 검증이나 일반 drilling DAE 증명을 대신하지 않는다. Build/test·solver/reference 실행은 문서 작업 범위가 아니므로 미실행이다. 독립 검산과 Numerical Review 결과는 해당 owner의 보고서가 소유하며 본 문서가 스스로 승인하지 않는다.