Aerodynamika, struktura i obciążenia pocisku

Pocisk jest jednym układem dynamicznym, zanim stanie się zbiorem dokumentów branżowych. Rozkład ciśnienia wytwarza siły i momenty; masa oraz bezwładność zamieniają je w ruch; ten ruch tworzy obciążenia bezwładnościowe; powłoka, ramy, zawiasy i połączenia prowadzą je przez konstrukcję; ugięcia wracają do aerodynamiki i pomiaru żyroskopowego; temperatura zmienia sztywność, docisk i oś sensora. Tezą tej monografii jest więc zdanie mocniejsze niż zwykłe „dyscypliny trzeba integrować”: zdolność płatowca można ocenić wiarygodnie tylko wtedy, gdy jeden kontrolowany opis konfiguracji przenosi te same stany masowe, osie, obciążenia, niepewności i kryteria między modelami oraz próbami aż do decyzji eksploatacyjnej.

Tekst pokazuje ten łańcuch na fikcyjnym obiekcie dydaktycznym OR-D1. Wszystkie jego dane są syntetyczne, nieoperacyjne i dobrane wyłącznie po to, aby rachunki dało się prześledzić. OR-D1 nie odwzorowuje istniejącego uzbrojenia. Podane wartości nie tworzą kompletnej geometrii wykonawczej, bazy osiągów ani nastaw sterowania i nie pozwalają odtworzyć współczesnego systemu. Tam, gdzie potrzebny jest przykład próby, chodzi o obojętny artykuł badawczy pozbawiony materiałów energetycznych i funkcji bojowych. Trzy studia wykonują ilościowo wybrane przekazania danych i decyzje, ale nie są pełną symulacją misji od atmosfery po trwałość. Monografia dowodzi spójności interfejsów i reguły decyzji; nie przypisuje OR-D1 kompletnego, numerycznie zamkniętego dowodu kwalifikacyjnego.

Płatowiec jako granica dowodu

Płatowiec oznacza tu wszystko, co ustala geometrię zewnętrzną, przenosi obciążenia albo wyznacza ruch mechaniczny: korpus i jego powłokę, ramy, owiewkę sensora, powierzchnie nośne i sterowe, zawiasy, mocowania urządzeń, interfejs napędu oraz punkty podparcia podczas transportu i startu. Jeden element może należeć do dwóch domen. Obudowa napędu może być główną częścią nośną; osłona sensora jest zarazem powierzchnią aerodynamiczną, przegrodą cieplną i odniesieniem osiowania; wspornik żyroskopu należy do awioniki funkcjonalnie, lecz do struktury dynamicznie.

Konfigurację OR-D1 identyfikuje się jako OR-D1/A-07. W syntetycznym rejestrze istnieją trzy stany masowe:

Stan	Masa $m$	Stacja środka masy $s_{CG}$ od płaszczyzny odniesienia	Zastosowanie
początkowy	$72{,}0\pm0{,}3 \mathrm{kg}$	$0{,}820\pm0{,}010 \mathrm{m}$	obciążenia początkowe i podparcie
punkt wysokiego $q_\infty$	$61{,}0\pm0{,}3 \mathrm{kg}$	$0{,}860\pm0{,}010 \mathrm{m}$	wspólny punkt trzech studiów
końcowy	$54{,}0\pm0{,}3 \mathrm{kg}$	$0{,}905\pm0{,}012 \mathrm{m}$	dolna masa i modalna regresja

Każdy zapis $\pm$ w tej tabeli oznacza syntetyczną niepewność standardową $u$ wyniku pomiaru albo estymaty; nie jest tolerancją geometryczną ani granicą specyfikacji.

W stanie wysokiego ciśnienia dynamicznego syntetyczny tensor bezwładności względem chwilowego środka masy, zapisany w osiach korpusu, wynosi

$$I_b= \begin{bmatrix} 0{,}32 & 0{,}00 & -0{,}02\\ 0{,}00 & 12{,}60 & 0{,}00\\ -0{,}02 & 0{,}00 & 12{,}50 \end{bmatrix} \mathrm{kg\,m^2}.$$

Przyjmujemy klasyczną definicję iloczynu bezwładności $J_{xz}=\int x_bz_b\,dm$; element tensora ma wtedy znak przeciwny, $(I_b)_{xz}=(I_b)_{zx}=-J_{xz}$. Macierz wyżej odpowiada więc $J_{xz}=+0{,}02 \mathrm{kg\,m^2}$, a nie „iloczynowi $I_{xz}=-0{,}02$”. Pozostawiono go celowo, aby pozorna symetria nie uzasadniała automatycznego wyzerowania sprzężeń. Tensor jest symetryczny i dodatnio określony: jego kolejne główne minory wynoszą $0{,}32$, $4{,}032$ i $50{,}39496$ w odpowiednich jednostkach, a wartości własne to $0{,}31997$, $12{,}50003$ i $12{,}60000 \mathrm{kg\,m^2}$. W modelu występują kontrolowane tablice $m(t)$, $s_{CG}(t)$ i $I_b(t)$, a nie tylko trzy wiersze tabeli. Liczby w tabeli są punktami kontrolnymi tych tablic.

Wspólne wymaganie demonstracyjne brzmi: konfiguracja A-07 ma zachować funkcję i dodatnie marginesy strukturalne w zadanej obwiedni, spełnić syntetyczne kryterium zapasu fazy co najmniej $35^\circ$ i zapasu wzmocnienia co najmniej $6 \mathrm{dB}$ w badanym kanale oraz utrzymać błąd osiowania sensora poniżej $0{,}25 \mathrm{mrad}$ podczas przejściowego nagrzewania. Te liczby są kryteriami ćwiczenia, nie zaleceniem dla realnego pocisku.

Cykl dowodu obejmuje także transport, magazynowanie, obsługę, noszenie, start, manewr, zdarzenia przejściowe i starzenie. Maksimum jednej zmiennej nie musi być przypadkiem krytycznym wszystkiego: najwyższe $q_\infty$ może sterować obciążeniem powierzchni, inny punkt kombinacją ściskania i zginania powłoki, a jeszcze inny temperaturą wewnętrznego mocowania. NASA-STD-5002B wymaga rozpoznania wszystkich krytycznych reżimów, odrębnych konfiguracji i źródeł obciążenia; jest to standard statków i ładunków kosmicznych, lecz jego logika kontroli przypadków oraz modeli jest użyteczna także tutaj, po świadomym dostosowaniu domeny.¹

Konwencja osi, znaków i wielkości odniesienia

Na potrzeby przykładu przyjmujemy jedną, jawną konwencję. Inne dokumentacje lotnicze i programy symulacyjne mogą stosować inne zwroty osi, definicje kąta ślizgu, kolejność obrotów albo punkt redukcji momentu. Wartość bez metadanych konwencji nie jest przenośna.

Płaszczyzna odniesienia $s=0$ jest umowną poprzeczną płaszczyzną w przedniej części obiektu. Nie jest czubkiem rzeczywistego pocisku ani informacją wykonawczą. Niezależna współrzędna stacyjna $s$ rośnie od nosa ku tyłowi. Prawoskrętny układ korpusu $b=(x_b,y_b,z_b)$ ma natomiast oś $+x_b$ do przodu, $+y_b$ w prawo, a $+z_b$ w dół; zachodzi $\boldsymbol e_{x_b}\times\boldsymbol e_{y_b}=\boldsymbol e_{z_b}$. Wzdłuż osi korpusu oba zwroty są przeciwne. Dla punktu osiowego $P$ o stacji $s_P$, przy początku układu $b$ w płaszczyźnie $s=0$, zachodzi

$$\boldsymbol r_P^b=-s_P\boldsymbol e_{x_b}, \qquad x_b(P)=-s_P.$$

Dodatnie momenty wynikają z reguły prawej dłoni: $L$ i prędkość $p_b$ oznaczają przechylenie wokół $+x_b$, $M$ i $q_b$ — obrót wokół $+y_b$, a $N$ i $r_b$ — obrót wokół $+z_b$. Określenia słowne, takie jak „uniesienie nosa”, zależą od konwencji orientacji i nie zastępują tego zapisu wektorowego.

Wektor prędkości względem powietrza w osiach korpusu zapisujemy jako $\boldsymbol V_a^b=[u,v,w]^T$, gdzie wszystkie składowe mają jednostkę $\mathrm{m/s}$. Dla $u>0$ definiujemy

$$V_a=\sqrt{u^2+v^2+w^2},\qquad \alpha=\operatorname{atan2}(w,u),\qquad \beta=\arcsin\left(\frac{v}{V_a}\right).$$

Dodatni kąt natarcia $\alpha$ odpowiada dodatniej składowej $w$; dodatni kąt ślizgu $\beta$ — dodatniej składowej $v$. Są to kąty wektora prędkości obiektu względem powietrza, nie kierunku, z którego napływa struga. Tak jawna definicja usuwa częsty błąd znaku wynikający z zamiany „wektora prędkości” na „wiatr względny”.

Układ inercjalny $i=(x_i,y_i,z_i)$ służy do opisu położenia i prędkości środka masy. Na krótkim odcinku demonstracyjnym uznajemy go za lokalnie nierotujący i płaski; model globalnej Ziemi wymagałby osobnego układu związanego z obracającą się Ziemią, modelu geodezji oraz konsekwentnych transformacji. NASA wykazała na porównaniach niezależnych symulatorów 6-DOF, że niejednoznaczność między prędkością kątową inercjalną i względem obracającej się planety oraz różne punkty redukcji momentów dają mierzalnie różne wyniki mimo pozornie tych samych danych wejściowych.²

Układ aerodynamiczny $w=(x_w,y_w,z_w)$ ma $x_w$ wzdłuż $\boldsymbol V_a$, $y_w$ w prawo w płaszczyźnie poprzecznej, a $z_w$ domyka układ prawoskrętny. Układ korpusu obraca się wraz z płatowcem; aerodynamiczny wraz z kierunkiem prędkości względem powietrza; inercjalny pozostaje odniesieniem ruchu postępowego. Macierz cosinusów kierunkowych $C_{ib}$ przekształca składowe z układu korpusu do inercjalnego, a $C_{bi}=C_{ib}^T$ działa odwrotnie. Dla dowolnego wektora $\boldsymbol a$ obowiązuje zatem $\boldsymbol a^i=C_{ib}\boldsymbol a^b$.

Siły aerodynamiczne i napędowe redukujemy najpierw do chwilowego środka masy. Jeśli baza aerodynamiczna przechowuje moment względem stałego punktu $R$ o położeniu $\boldsymbol r_R$, to

$$\boldsymbol M_{CG}=\boldsymbol M_R+(\boldsymbol r_R-\boldsymbol r_{CG})\times\boldsymbol F.$$

W OR-D1 pole odniesienia wynosi $S=0{,}0250 \mathrm{m^2}$, długość odniesienia $\ell=1{,}60 \mathrm m$, a stały punkt redukcji bazy leży na osi korpusu w stacji $s_R=0{,}800 \mathrm m$, czyli $\boldsymbol r_R^b=-0{,}800\boldsymbol e_{x_b} \mathrm m$. Dla osiowych $R$ i $CG$ transformacja przyjmuje zatem jednoznaczną postać

$$\boldsymbol M_{CG}=\boldsymbol M_R+(s_{CG}-s_R)\boldsymbol e_{x_b}\times\boldsymbol F.$$

Przykładowo w stanie wysokiego $q_\infty$ jest $s_{CG}-s_R=+0{,}060 \mathrm m$, mimo że składowa wektora $\boldsymbol r_R-\boldsymbol r_{CG}$ wzdłuż $x_b$ jest dodatnia. Są to umowne wielkości do bezwymiarowania; $S$ nie musi być fizycznym polem konkretnej powierzchni, a $\ell$ nie ujawnia kompletnej długości obiektu. Moment względem $R$ i moment względem $CG$ nigdy nie mogą nosić tej samej etykiety bez transformacji.

Ciśnienie dynamiczne oznaczamy wyłącznie $q_\infty$:

$$q_\infty=\frac{1}{2}\rho V_a^2,$$

gdzie $\rho$ jest gęstością powietrza w $\mathrm{kg/m^3}$. Symbol $q_b$ jest zarezerwowany dla prędkości kątowej pochylenia w $\mathrm{rad/s}$. Indeks b odróżnia ją od $q_\infty$ także w tabelach i kodzie. Brak tego rozróżnienia jest błędem interfejsu danych, nie kosmetyką notacji.

Od atmosfery do kontrolowanej bazy aerodynamicznej

Atmosfera jest wejściem do aerodynamiki, termiki i napędu. Dla gazu doskonałego lokalna prędkość dźwięku wynosi $a=\sqrt{\gamma R T}$, liczba Macha $M=V_a/a$, a liczba Reynoldsa

$$Re=\frac{\rho V_a\ell}{\mu},$$

gdzie $\gamma$ jest stosunkiem ciepeł właściwych, $R$ — indywidualną stałą gazową w $\mathrm{J/(kg\,K)}$, $T$ — temperaturą statyczną w kelwinach, a $\mu$ — lepkością dynamiczną w $\mathrm{Pa\,s}$. Standardowa atmosfera jest modelem odniesienia, nie zapisem pogody podczas próby; U.S. Standard Atmosphere 1976 definiuje spójne profile wielkości termodynamicznych, ale rzeczywisty punkt wymaga pomiaru lub zaakceptowanego modelu odchyłek.³

W syntetycznym punkcie wysokiego $q_\infty$ przyjmujemy $q_\infty=45{,}0 \mathrm{kPa}$, $\rho=0{,}750 \mathrm{kg/m^3}$, $T=250 \mathrm K$, $\gamma=1{,}40$, $R=287{,}05 \mathrm{J/(kg\,K)}$ i $\mu=1{,}60\cdot10^{-5} \mathrm{Pa\,s}$. Prędkość wynika z ciśnienia dynamicznego:

$$V_a=\sqrt{\frac{2q_\infty}{\rho}} =\sqrt{\frac{2(45\,000 \mathrm{Pa})}{0{,}750 \mathrm{kg/m^3}}} =346 \mathrm{m/s}.$$

Sprawdzenie wymiaru daje $\sqrt{\mathrm{Pa}/(\mathrm{kg/m^3})}=\sqrt{\mathrm{m^2/s^2}}=\mathrm{m/s}$. Dalej

$$a=\sqrt{1{,}40(287{,}05 \mathrm{J/(kg\,K)})(250 \mathrm K)}=317 \mathrm{m/s},$$

$$M=\frac{346}{317}=1{,}09, \qquad Re=\frac{(0{,}750)(346)(1{,}60)}{1{,}60\cdot10^{-5}}=2{,}60\cdot10^7.$$

Jest to celowo punkt lekko naddźwiękowy, w którym podobne $q_\infty$ nie gwarantuje podobnej aerodynamiki do punktu poddźwiękowego. Fale uderzeniowe, przesunięcie separacji, lokalna grubość warstwy przyściennej i interferencja powierzchnia–korpus mogą zmieniać współczynniki mimo tej samej skali siły $q_\infty S$. Dopasowanie modelu tunelowego tylko w $M$ bez kontroli $Re$, chropowatości, szczelin i podpory również nie zapewnia podobieństwa.

Rozkład ciśnienia $p_s$ i naprężenia stycznego $\boldsymbol\tau_s$ na powierzchni $A$ daje siłę i moment przez całkowanie:

$$\boldsymbol F_a=\int_A\left[-p_s\boldsymbol n+\boldsymbol\tau_s\right]dA, \qquad \boldsymbol M_R=\int_A(\boldsymbol r-\boldsymbol r_R)\times \left[-p_s\boldsymbol n+\boldsymbol\tau_s\right]dA.$$

W praktyce wynik zapisuje się bezwymiarowo, na przykład $X=q_\infty S C_X$, $Y=q_\infty S C_Y$, $Z=q_\infty S C_Z$ oraz $L=q_\infty S\ell C_l$, $M=q_\infty S\ell C_m$, $N=q_\infty S\ell C_n$. Każdy współczynnik zależy od argumentów bazy:

$$C_k=C_k(M,Re,\alpha,\beta,\hat p,\hat q,\hat r,\boldsymbol\delta,\mathcal C),$$

gdzie $k$ wskazuje siłę albo moment, $\boldsymbol\delta$ jest wektorem wychyleń, $\mathcal C$ identyfikuje dyskretną konfigurację, a znormalizowane prędkości kątowe wynoszą

$$\hat p=\frac{p_b\ell}{2V_a},\qquad \hat q=\frac{q_b\ell}{2V_a},\qquad \hat r=\frac{r_b\ell}{2V_a}.$$

Pochodna $C_{m_q}$ oznacza zatem $\partial C_m/\partial\hat q$, a nie pochodną względem ciśnienia dynamicznego. W małym otoczeniu punktu pracy można napisać

$$C_m\approx C_{m0}+C_{m_\alpha}\alpha+C_{m_q}\hat q+C_{m_\delta}\delta,$$

lecz zakres linearyzacji musi być opisany. Przejście przez falę, duże wychylenie, otwarcie elementu albo nieliniowy luz zawiasu są zmianą modelu, nie kolejnym punktem tej samej prostej.

Kontrolowana baza aerodynamiczna jest wyrobem konfiguracyjnym: zawiera wartości, metody interpolacji, korekty podpory, identyfikatory geometrii i siatek, zakres ważności oraz niepewność. Missile DATCOM jest oficjalnie opisanym narzędziem do szybkich, ekonomicznych predykcji szerokiej klasy konfiguracji na etapie wstępnym; instrukcja wylicza osobno warunki lotu, wielkości odniesienia, geometrię, wychylenia, dane doświadczalne i pochodne dynamiczne. Ta szerokość jest zaletą, ale raport nie zmienia metody w dowód kwalifikacyjny ani nie usuwa ograniczeń konfiguracji nietypowej.⁴

W OR-D1 baza ADB-A07-R3 powstaje z trzech warstw. Model półempiryczny wypełnia obwiednię, obliczenia przepływowe rozstrzygają interferencję w wybranych obszarach, a tunel dostarcza zintegrowanych sił, momentów i ciśnień do korelacji. Dane z próby w locie służą dopiero do identyfikacji odchyłek kompletnego obiektu. Zgodność oporu przy $\alpha=0$ nie waliduje momentu przy manewrze ani lokalnego momentu zawiasowego. Dla każdego pola bazy zapisuje się: najlepszą estymatę, standardową niepewność, korelacje z sąsiednimi współczynnikami oraz flagę interpolacji lub ekstrapolacji.

Niepewności nie wolno dodawać bez rozpoznania zależności. Jeśli $F=q_\infty S C$, a $q_\infty$, $S$ i $C$ są niezależne, względna standardowa niepewność w przybliżeniu pierwszego rzędu wynosi

$$\left(\frac{u_F}{F}\right)^2= \left(\frac{u_q}{q_\infty}\right)^2+ \left(\frac{u_S}{S}\right)^2+ \left(\frac{u_C}{C}\right)^2.$$

Jeżeli ten sam błąd ustawienia modelu wpływa na $\alpha$ i moment, potrzebny jest człon kowariancji. Naturalny rozrzut egzemplarzy i atmosfery jest niepewnością losową; brak wiedzy o modelu separacji lub przejścia jest niepewnością poznawczą. Tej drugiej nie należy zamieniać w wygodny, wąski rozkład normalny bez danych.

Ruch bryły sztywnej, lokalny model zamrożony i granice sterowania

Stan 6-DOF obejmuje położenie środka masy $\boldsymbol r^i$, prędkość $\boldsymbol v^b=[u,v,w]^T$, orientację oraz prędkość kątową $\boldsymbol\omega^b=[p_b,q_b,r_b]^T$. W tej monografii nie pozorujemy pełnego wyprowadzenia układu otwartego. Rozpatrujemy tylko krótki odcinek wokół punktu A-07, na którym masa, stacja środka masy i tensor są zamrożone: $m_A=61{,}0 \mathrm{kg}$, $s_{CG,A}=0{,}860 \mathrm m$ i $I_{b,A}$ jak wyżej. Równania translacji, zapisane w osiach korpusu, mają wtedy postać

$$\dot{\boldsymbol r}^{i}=C_{ib}\boldsymbol v^b,$$

$$m_A\left(\dot{\boldsymbol v}^{b}+\boldsymbol\omega^b\times\boldsymbol v^b\right) =\boldsymbol F_a^b+\boldsymbol F_T^b+m_AC_{bi}\boldsymbol g^i+\boldsymbol F_{inne}^b.$$

$\boldsymbol F_a^b$ jest siłą aerodynamiczną, $\boldsymbol g^i$ przyspieszeniem grawitacyjnym, a $\boldsymbol F_{inne}^b$ sumą pozostałych sił zewnętrznych. $\boldsymbol F_T^b$ nie jest wyprowadzana z powyższego równania: dostarcza ją odrębny, kontrolowany model objętości kontrolnej strumienia napędowego jako wypadkową już uwzględniającą transport pędu przez granicę. W lokalnym modelu zamrożonym nie występuje $\dot m$.

Równanie obrotu względem zamrożonego środka masy A-07 brzmi

$$I_{b,A}\dot{\boldsymbol\omega}^{b} +\boldsymbol\omega^b\times\left[I_{b,A}\boldsymbol\omega^b\right] =\boldsymbol M_a^b+\boldsymbol M_T^b+\boldsymbol M_{inne}^b.$$

Człon $\boldsymbol\omega\times(I\boldsymbol\omega)$ wyraża sprzężenie żyroskopowe. $\boldsymbol M_a^b$ jest momentem aerodynamicznym, a $\boldsymbol M_{inne}^b$ sumą pozostałych momentów zewnętrznych. Wypadkowy $\boldsymbol M_T^b$, podobnie jak $\boldsymbol F_T^b$, jest wejściem z odrębnego modelu strumienia i jest zredukowany do tego samego zamrożonego $CG$; obejmuje moment sił na granicy oraz transport momentu pędu przez granicę przy założeniach tego modelu. Nie jest domyślnym „momentem ciągu” bez metadanych punktu redukcji.

Pełny przypadek zmiennej masy wymaga bilansu na jawnie wybranej objętości kontrolnej. Oprócz zmiany pędu i momentu pędu materii pozostającej wewnątrz trzeba uwzględnić strumienie obu wielkości przez każdy wlot i wylot, ruch przyjętego punktu redukcji, rozkład prędkości wewnętrznej oraz ewentualny względny ruch i deformację masy. Zależnie od granicy dochodzą siły i pary na tej granicy. Samo dopisanie $\dot I_b\boldsymbol\omega$ albo $\dot m\boldsymbol v$ nie zamyka takiego bilansu. Tablice $m(t)$, $s_{CG}(t)$ i $I_b(t)$ służą w pełnym modelu misji do wyboru kolejnych punktów lokalnych, lecz nie są tu podstawą ogólnego równania układu otwartego. Ani przypadki kontrolne 6-DOF z ², ani model elastycznego samolotu z ⁶ nie są źródłem takiego wyprowadzenia.

Orientację całkujemy kwaternionem jednostkowym w kolejności skalar–wektor, $\boldsymbol e=[e_0,e_1,e_2,e_3]^T=(e_0,\boldsymbol e_v)$. Jest to kwaternion aktywnej rotacji b→i: wektor o składowych $\boldsymbol a^b$ reprezentujemy kwaternionem czystym $(0,\boldsymbol a^b)$ i przekształcamy jako $(0,\boldsymbol a^i)=\boldsymbol e\otimes(0,\boldsymbol a^b)\otimes\boldsymbol e^*$, gdzie $\boldsymbol e^*=(e_0,-\boldsymbol e_v)$. Stosujemy iloczyn Hamiltona

$$(a_0,\boldsymbol a)\otimes(b_0,\boldsymbol b)= (a_0b_0-\boldsymbol a\mathbin{\cdot}\boldsymbol b, a_0\boldsymbol b+b_0\boldsymbol a+\boldsymbol a\mathbin{\times}\boldsymbol b).$$

Kolejność czynników jest częścią definicji: poniżej występuje prawostronne mnożenie $\boldsymbol e\otimes\boldsymbol\omega_q$ przez czysty kwaternion prędkości $\boldsymbol\omega_q=(0,p_b,q_b,r_b)$, a nie lewostronne $\boldsymbol\omega_q\otimes\boldsymbol e$. Równanie kinematyczne ma zatem postać

$$\dot{\boldsymbol e}=\frac{1}{2} \begin{bmatrix} 0&-p_b&-q_b&-r_b\\ p_b&0&r_b&-q_b\\ q_b&-r_b&0&p_b\\ r_b&q_b&-p_b&0 \end{bmatrix}\boldsymbol e, \qquad \boldsymbol e^T\boldsymbol e=1.$$

Odpowiada jej dokładnie macierz działająca na kolumnę składowych, $\boldsymbol a^i=C_{ib}(\boldsymbol e)\boldsymbol a^b$:

$$C_{ib}(\boldsymbol e)= \begin{bmatrix} e_0^2+e_1^2-e_2^2-e_3^2 & 2(e_1e_2-e_0e_3) & 2(e_1e_3+e_0e_2)\\ 2(e_1e_2+e_0e_3) & e_0^2-e_1^2+e_2^2-e_3^2 & 2(e_2e_3-e_0e_1)\\ 2(e_1e_3-e_0e_2) & 2(e_2e_3+e_0e_1) & e_0^2-e_1^2-e_2^2+e_3^2 \end{bmatrix}, \qquad C_{bi}=C_{ib}^T.$$

Zgodność znaków można sprawdzić bez odwołania do pamięci. Z definicji iloczynu Hamiltona $\dot{\boldsymbol e}=\tfrac12\boldsymbol e\otimes\boldsymbol\omega_q$ daje kolejno $\dot e_0=-\tfrac12\boldsymbol e_v\cdot\boldsymbol\omega^b$ oraz $\dot{\boldsymbol e}_v=\tfrac12(e_0\boldsymbol\omega^b+\boldsymbol e_v\times\boldsymbol\omega^b)$, czyli dokładnie cztery wiersze wydrukowanej macierzy. Przy $\boldsymbol e(0)=(1,0,0,0)$ otrzymujemy ponadto $\dot C_{ib}(0)=[\boldsymbol\omega^b]_\times$, a ogólnie $\dot C_{ib}=C_{ib}[\boldsymbol\omega^b]_\times$; jest to znak właściwy dla aktywnego mapowania b→i i prędkości wyrażonej w osiach $b$.

Test jednostkowy dodatniego obrotu aktywnego $+90^\circ$ wokół osi $+z$ przyjmuje $\boldsymbol e=(\sqrt2/2,0,0,\sqrt2/2)$. Bezpośrednie podstawienie daje

$$C_{ib}= \begin{bmatrix}0&-1&0\\1&0&0\\0&0&1\end{bmatrix},\qquad C_{ib}\boldsymbol e_{x_b}=\boldsymbol e_{y_i},\qquad C_{ib}\boldsymbol e_{y_b}=-\boldsymbol e_{x_i},\qquad C_{ib}\boldsymbol e_{z_b}=\boldsymbol e_{z_i}.$$

Transpozycja dałaby obrót $-90^\circ$ i oblałaby ten test. Testy zachowania normy, zachowania energii w ruchu bez sił i porównania niezależnych implementacji pozostają dalszą częścią weryfikacji. NASA-owski zestaw przypadków kontrolnych dla symulacji lotu pokazuje właśnie potrzebę wspólnych definicji układów, początkowych prędkości kątowych, grawitacji i punktów momentowych.²

Model aktuatora jest częścią płatowca dynamicznego, nawet jeśli jego elektronikę opisuje inny artykuł. Dla celów poglądowych stosujemy opóźniony człon pierwszego rzędu

$$\tau_a\dot\delta+\delta=K_a\delta_c(t-T_d), \qquad |\delta|\le\delta_{max},\qquad |\dot\delta|\le\dot\delta_{max},$$

gdzie $\delta_c$ jest komendą, $\delta$ rzeczywistym wychyleniem, $\tau_a$ stałą czasową, $T_d$ opóźnieniem, $K_a$ wzmocnieniem, a ograniczenia opisują nasycenie położenia i prędkości. Liczb $K_a$, $\delta_{max}$ i $\dot\delta_{max}$ dla żadnego realnego systemu tu nie podajemy. Luz, tarcie, podatność zasilania i zależność momentu od temperatury mogą wymagać modelu nieliniowego.

Pomiar sensora zapisujemy jako

$$\boldsymbol y_s(t)=H_s\boldsymbol x(t-T_s)+\boldsymbol b_s(t)+\boldsymbol n_s(t)+\boldsymbol y_{flex}(t)+\boldsymbol y_{therm}(t),$$

gdzie $H_s$ wybiera mierzone składowe stanu $\boldsymbol x$, $T_s$ jest opóźnieniem, $\boldsymbol b_s$ błędem systematycznym, $\boldsymbol n_s$ szumem, a dwa ostatnie człony są odpowiednio ruchem lokalnym od elastyczności i dryfem cieplnym. Filtr przeciwaliasingowy i częstotliwość próbkowania muszą wejść do modelu czasu dyskretnego. Sensor nie mierzy abstrakcyjnej osi nominalnej, lecz ruch i temperaturę własnego miejsca montażu.

Prawo sterowania otrzymuje więc nie „prawdziwy kąt”, lecz opóźnioną estymatę z błędami. Sterowność wymaga jednocześnie wystarczającego momentu, pasma aktuatora, zapasu przed nasyceniem, obserwowalności stanu i odporności na niepewność. Ogranicznik komendy może chronić strukturę, ale zmienia osiągalny ruch; filtr wycinający postać drgań może chronić zapas stabilności, ale wnosi fazę; przesunięcie środka masy zmienia zarówno moment aerodynamiczny, jak i potrzebną pracę aktuatora. Z tego powodu model 6-DOF jest interfejsem między aerodynamiką, masą, strukturą, sensorami oraz naprowadzaniem, nawigacją i sterowaniem (GNC), a nie tylko „symulatorem trajektorii”.

Od ruchu do obciążeń przekrojowych

Wynik symulacji ruchu staje się obciążeniem konstrukcji dopiero po rozłożeniu wypadkowych na drogę przenoszenia sił. Rozważmy wydzielony segment $s\in[s_1,s_2]$, gdzie $s_1<s_2$. W przekroju $s$ definiujemy $\boldsymbol R(s)$ i $\boldsymbol M(s)$ jako siłę i parę, które materiał po stronie mniejszych stacji wywiera na materiał po stronie większych stacji; $\boldsymbol M(s)$ jest zredukowany do punktu osi przekroju $\boldsymbol r(s)$. Na lewym brzegu wybranego segmentu działają więc $+\boldsymbol R(s_1),+\boldsymbol M(s_1)$, a na prawym ich odpowiedniki ze znakiem przeciwnym. Równowaga dynamiczna ma postać

$$\boldsymbol R(s_1)-\boldsymbol R(s_2) +\int_{s_1}^{s_2}\boldsymbol w_a(s,t)\,ds +\sum_k\boldsymbol F_k(t) =\int_{s_1}^{s_2}\lambda(s,t)\boldsymbol a(s,t)\,ds,$$

$$\boldsymbol M(s_1)-\boldsymbol M(s_2) +[\boldsymbol r(s_1)-\boldsymbol r_R]\times\boldsymbol R(s_1) -[\boldsymbol r(s_2)-\boldsymbol r_R]\times\boldsymbol R(s_2) +\int_{s_1}^{s_2}[\boldsymbol r(s)-\boldsymbol r_R]\times\boldsymbol w_a(s,t)\,ds +\sum_k[\boldsymbol r_k-\boldsymbol r_R]\times\boldsymbol F_k(t) +\sum_k\boldsymbol Q_k(t) =\int_{s_1}^{s_2}[\boldsymbol r(s)-\boldsymbol r_R]\times\lambda(s,t)\boldsymbol a(s,t)\,ds +\int_{s_1}^{s_2}\dot{\boldsymbol h}_{loc}(s,t)\,ds.$$

$\boldsymbol r_R$ jest jednym, dowolnie wybranym punktem redukcji całego bilansu. $\boldsymbol w_a$ oznacza aerodynamiczne obciążenie rozłożone w $\mathrm{N/m}$, $\lambda$ masę liniową w $\mathrm{kg/m}$, a $\boldsymbol a$ przyspieszenie osi rozkładu w $\mathrm{m/s^2}$. Siły skupione $\boldsymbol F_k$ obejmują reakcje podpór, zespołów wewnętrznych i aktuatorów; $\boldsymbol Q_k$ są skupionymi parami momentów. Ostatni człon reprezentuje zmianę lokalnego momentu pędu przekroju, jeżeli model liniowy zachowuje bezwładność obrotową; w modelu czysto belkowym może być zerowy. Ramiona sił przekrojowych znikają tylko przy redukowaniu każdej pary $(\boldsymbol R,\boldsymbol M)$ do wspólnego $\boldsymbol r_R$, czego tutaj nie robimy.

Znak sprawdza prosty wspornik. Niech segment od przekroju $s$ do wolnego końca $L$ obciąża siła $+P\boldsymbol e_{z_b}$ w $L$, bez bezwładności i innych obciążeń. Ponieważ $\boldsymbol r(L)-\boldsymbol r(s)=-(L-s)\boldsymbol e_{x_b}$, równania dają

$$\boldsymbol R(s)=-P\boldsymbol e_{z_b},\qquad \boldsymbol M(s)=-P(L-s)\boldsymbol e_{y_b}.$$

Skalarne składowe są rzutami wcześniej zdefiniowanych wektorów: $M_y(s)=\boldsymbol e_{y_b}\mathbin{\cdot}\boldsymbol M(s)$, analogicznie dla $M_z,T,V_y,V_z$ i $N$. Nie są modułami. W teście powyżej $M_y=-P(L-s)$, a jego moduł rośnie od zera na końcu wolnym do $PL$ przy utwierdzeniu, zgodnie z elementarnym testem wspornika.⁵ Znaki wszystkich składowych dziedziczą zdefiniowane osie i stronę cięcia.

W praktyce obciążenia rozłożone pozyskuje się z bazy aerodynamicznej lub z ciśnień powierzchniowych, a człon bezwładnościowy z historii 6-DOF i modelu masy. Oba muszą dotyczyć tej samej chwili i konfiguracji. Użycie maksymalnego ciśnienia dynamicznego z jednej chwili, maksymalnego kąta natarcia z drugiej i najniekorzystniejszego środka masy z trzeciej tworzy kombinację konserwatywną tylko wtedy, gdy jawnie wykaże się jej fizyczną osiągalność albo świadomie ustanowi ją jako obwiednię projektową. W przeciwnym razie jest to niezidentyfikowany przypadek, którego ani nie można odtworzyć, ani sensownie skorelować z próbą.

Droga obciążenia jest grafem, nie pojedynczą strzałką. Na przykład siła na ruchomej powierzchni przechodzi przez jej dźwigar i łożyska do osi zawiasu; część momentu równoważy aktuator, reakcje zawiasu trafiają do ramy, rama rozprowadza je do powłoki, a powłoka przenosi ścinanie i zginanie do sąsiednich przedziałów. Globalny moment powierzchni względem środka masy nie jest momentem zawiasowym. Pierwszy służy do równania ruchu całego obiektu, drugi wynika z rozkładu ciśnienia względem osi obrotu i obciąża napęd oraz mocowanie. Pomylenie tych dwóch wielkości może zawyżyć jedną reakcję o rząd wielkości, a inną pominąć.

Reakcje wyznacza się dla wszystkich stanów masowych i istotnych faz przejściowych. Stan A-07 przy wysokim $q_\infty$ jest krytyczny dla badanego mocowania, lecz stan cięższy może być krytyczny dla podparcia wyposażenia wskutek większej siły bezwładności, a stan końcowy — dla sterowania wskutek przesunięcia stacji $s_{CG}$. Lista przypadków obciążenia powinna zatem wiązać identyfikator konfiguracji, czas, atmosferę, stan ruchu, wariant bazy aerodynamicznej, rozkład masy, model połączeń i wersję oprogramowania. NASA-STD-5002B wymaga analogicznej identyfikowalności źródeł, konfiguracji, wszystkich krytycznych reżimów oraz niepewności; zastosowanie tego standardu tutaj jest zasadą organizacji dowodu, nie stwierdzeniem zgodności rakietowego programu z wymaganiami statków kosmicznych.¹

Odpowiedź konstrukcji i kryterium marginesu

Dla smukłego elementu w zakresie liniowym pierwszą kontrolę przekroju daje naprężenie normalne

$$\sigma_x=\frac{N}{A}+\frac{M_y z}{I_y}-\frac{M_z y}{I_z},$$

gdzie $A$ jest polem przekroju, $I_y,I_z$ jego momentami bezwładności pola, a $y,z$ współrzędnymi badanego punktu. W skrajnym włóknie jednokierunkowego zginania zapis wartości maksymalnej upraszcza się do $\sigma_b=|M_y|/W$ (albo odpowiednio $|M_z|/W$), gdzie $W=I/c$ jest wskaźnikiem wytrzymałości, a $c$ odległością włókna od osi obojętnej. Skręcanie zamkniętej cienkościennej komórki wygodniej opisać przepływem ścinającym $q_s=T/(2A_m)$ niż wzorem dla pełnego wału; $A_m$ jest polem zamkniętym przez linię środkową ścianki, a naprężenie wynosi $\tau=q_s/t$ przy grubości $t$. Te wzory są modelami belkowymi. Otwory, zakończenia żeber, koncentracje przy łącznikach i lokalny docisk wymagają modelu lokalnego albo zweryfikowanego współczynnika koncentracji.

W stanie wieloosiowym dla materiału ciągliwego stosujemy jako jedną z kontroli naprężenie zastępcze Hubera–Misesa

$$\sigma_{VM}=\sqrt{\frac{1}{2}\left[(\sigma_1-\sigma_2)^2+(\sigma_2-\sigma_3)^2+(\sigma_3-\sigma_1)^2\right]},$$

gdzie $\sigma_1,\sigma_2,\sigma_3$ są naprężeniami głównymi. Nie jest to uniwersalne kryterium wszystkich uszkodzeń: dla kompozytu, kruchej powłoki, kleju, łożyska czy połączenia śrubowego trzeba sprawdzić właściwy mechanizm. Misiak wyprowadza związek obciążenia z siłami przekrojowymi, naprężeniami, ugięciem, hipotezami wytężenia oraz statecznością; użyte tutaj równania są klasycznymi redukcjami tych relacji, a nie substytutem analizy szczegółu.⁵

W płaszczyźnie $x_b$–$z_b$ definiujemy $v_z(s)$ jako podpisane ugięcie osi belki: $v_z>0$ oznacza przemieszczenie w $+z_b$. Ponieważ stacja $s$ rośnie przeciwnie do $+x_b$, a $M_y$ jest podpisanym rzutem wektora momentu przekrojowego zdefiniowanego wyżej, relacja Eulera–Bernoulliego ma w tej konwencji postać

$$E I_y(s)\frac{d^2v_z}{ds^2}=-M_y(s).$$

W tym zapisie $E$ jest modułem Younga, a $I_y$ momentem bezwładności pola względem osi $y_b$. Ten sam test wspornika domyka znak: dla $M_y=-P(L-s)$ otrzymujemy $EI_yv_z''=P(L-s)$; przy $v_z(0)=v_z'(0)=0$ całkowanie daje $v_z(s)=P(Ls^2/2-s^3/6)/(EI_y)$ oraz $v_z(L)=PL^3/(3EI_y)>0$, zgodnie z kierunkiem przyłożonej siły $+P\boldsymbol e_{z_b}$. Kryterium sztywności bywa ważniejsze od granicy plastyczności: małe ugięcie mocowania może zmienić oś sensora, szczelinę uszczelnienia albo luz mechanizmu. Z kolei stateczność nie może być oceniona tylko przez porównanie naprężenia z granicą plastyczności. Idealna siła krytyczna pręta $P_{cr}=\pi^2EI/(K L)^2$, w której $L$ jest długością niepodpartą, a $K$ współczynnikiem długości wyboczeniowej zależnym od zamocowania, wskazuje zależności, lecz cienka powłoka jest wrażliwa na owalizację, wycięcia, mimośród i imperfekcje. Nośność z analizy własnej wyboczenia należy obniżyć współczynnikiem niedoskonałości uzasadnionym pomiarem, literaturą lub próbą reprezentatywną.

Dla każdego mechanizmu zapisujemy margines

$$MS=\frac{R_{dop}}{L_{obl}}-1,$$

gdzie $R_{dop}$ jest dopuszczalną nośnością po wszystkich redukcjach materiałowych, temperaturowych, środowiskowych i produkcyjnych, a $L_{obl}$ odpowiada obciążeniu oraz temu samemu sposobowi jego wyrażenia. Dodatni $MS$ oznacza spełnienie tylko wskazanego kryterium. Nie wolno dzielić naprężenia dopuszczalnego przez siłę ani porównywać nośności niszczącej z obciążeniem granicznym bez właściwego współczynnika. W tej monografii „obciążenie graniczne” (limit load) oznacza największe obciążenie przewidywane w zdefiniowanej eksploatacji, a „obciążenie niszczące” (ultimate load) — obciążenie graniczne pomnożone przez przyjęty współczynnik bezpieczeństwa. Są to definicje robocze; wartość współczynnika musi wynikać z programu i nie jest automatycznie przenoszona z norm kosmicznych.¹⁴

Studium I: manewr przy wysokim ciśnieniu dynamicznym

Wymaganie i konfiguracja. Mocowanie powierzchni w konfiguracji A-07 ma przenieść obciążenie graniczne bez uplastycznienia, utraty stateczności i poślizgu połączenia oraz obciążenie niszczące bez zerwania. Wymagamy $MS\geq0$ dla każdego mechanizmu. Przypadek LC-A07-184 jest jawnym przekazaniem danych z lokalnego modelu 6-DOF: w chwili umownej $t_*=18{,}4 \mathrm s$ pakiet stanu wskazuje $q_\infty=45{,}0 \mathrm{kPa}$, masę $61{,}0 \mathrm{kg}$, $s_{CG}=0{,}860 \mathrm{m}$, małe prędkości kątowe oraz identyfikatory ADB-A07-R3 i tablicy rozkładu masy. Zapisane w tej samej chwili rozkłady $\boldsymbol w_a(s,t_*)$ i $\lambda(s,t_*)\boldsymbol a(s,t_*)$ są wejściem do bilansu przekrojowego. Ich redukcja dla badanej powierzchni daje $C_N=0{,}70$ i odległość środka parcia od $CG$ równą $l_N=0{,}550 \mathrm m$. Nie są to wymiary ani parametry żadnego uzbrojenia, lecz zredukowany przypadek do sprawdzenia przepływu danych. Monografia nie oblicza poprzedzającej trajektorii; nie twierdzi więc, że wykonała pełny łańcuch 6-DOF.

Konkurują cztery hipotezy anomalii obserwowanej wcześniej w odkształceniu ramy: H1 — baza aerodynamiczna zaniża siłę normalną; H2 — błędne położenie środka masy zawyża globalne ramię; H3 — baza momentu zawiasowego zaniża moment lokalny; H4 — połączenie ramy z powłoką traci napięcie wstępne i zmienia drogę obciążenia. Model uproszczony zakłada quasi-statyczne maksimum, liniową sprężystość przed poślizgiem, dwie reakcje zawiasu, zamknięty przekrój powłoki oraz brak równoczesnego rezonansu. Założenia te są ważne tylko dla przedstawionego zdarzenia; badanie dynamiczne pozostaje osobnym przypadkiem.

Siła i moment globalny wynoszą

$$N_a=q_\infty S C_N =(45\,000 \mathrm{N/m^2})(0{,}0250 \mathrm{m^2})(0{,}70) =787{,}5 \mathrm{N},$$

$$M_{CG}=N_a l_N=(787{,}5 \mathrm{N})(0{,}55 \mathrm{m})=433{,}1 \mathrm{N\,m}.$$

Standardowe niepewności względne przyjmujemy jawnie: $u_q/q=2{,}0\%$, $u_S/S=0{,}8\%$, $u_C/C=7{,}0\%$ oraz $u_l/l=1{,}8\%$. Są traktowane jako niezależne wyłącznie na potrzeby tego rachunku. Zatem

$$u_N=787{,}5\sqrt{0{,}020^2+0{,}008^2+0{,}070^2}=57{,}7 \mathrm{N},$$

$$u_M=433{,}1\sqrt{(57{,}7/787{,}5)^2+0{,}018^2}=32{,}7 \mathrm{N\,m}.$$

Programowy mnożnik przejściowej nadwyżki dynamicznej ma syntetyczną wartość $k_d=1{,}10$. Obciążenie graniczne definiujemy w tym studium jako wartość nominalną powiększoną o dwa odchylenia standardowe i pomnożoną przez $k_d$:

$$N_{lim}=1{,}10(787{,}5+2\cdot57{,}7)=993 \mathrm{N},$$

$$M_{lim}=1{,}10(433{,}1+2\cdot32{,}7)=548 \mathrm{N\,m}.$$

Przy przyjętym tylko dla OR-D1 współczynniku $1{,}50$ wartości niszczące wynoszą odpowiednio $N_{ult}=1{,}49 \mathrm{kN}$ i $M_{ult}=822 \mathrm{N\,m}$. Zaokrąglenie do trzech cyfr jest współmierne do danych wejściowych; liczby nie są dopuszczeniami dla realnego sprzętu.

Kontrolowana baza lokalna daje nominalny moment zawiasowy $H=34\pm5 \mathrm{N\,m}$, gdzie $5 \mathrm{N\,m}$ jest niepewnością standardową. Jego wartość graniczna $H_{lim}=1{,}10(34+2\cdot5)=48{,}4 \mathrm{N\,m}$ jest ponad jedenastokrotnie mniejsza od $M_{lim}$, co ilustruje, dlaczego wielkości nie wolno utożsamiać. Reakcje dwóch łożysk, z rozdziałem $55/45$ z modelu sztywności, wynoszą $546$ i $447 \mathrm{N}$. Rama przekazuje do powłoki ich sumę $993 \mathrm{N}$ oraz lokalne skręcanie $48{,}4 \mathrm{N\,m}$, natomiast przekrój korpusu przy środku masy równoważy zginanie $548 \mathrm{N\,m}$. Każda liczba ma więc przypisaną krawędź drogi obciążenia.

Model ramy i powłoki daje zredukowany wskaźnik przekroju $W_{eff}=3{,}60\cdot10^{-6} \mathrm{m^3}$ po uwzględnieniu otworów. Naprężenie zginające przy obciążeniu granicznym wynosi

$$\sigma_b=\frac{548 \mathrm{N\,m}}{3{,}60\cdot10^{-6} \mathrm{m^3}} =152 \mathrm{MPa}.$$

Po dodaniu lokalnego ścinania i docisku model szczegółowy daje $\sigma_{VM}=157 \mathrm{MPa}$. Dopuszczalna granica plastyczności po redukcji temperaturowej i produkcyjnej wynosi syntetycznie $240 \mathrm{MPa}$, stąd $MS_y=240/157-1=0{,}53$. Zredukowana nośność wyboczeniowa odpowiada $190 \mathrm{MPa}$ naprężenia ściskającego, więc $MS_b=190/152-1=0{,}25$. Przy obciążeniu niszczącym model liniowy daje $1{,}50\cdot157=236 \mathrm{MPa}$ wobec nośności zerwania $330 \mathrm{MPa}$, czyli $MS_u=0{,}40$. Krytyczny jest zatem margines stateczności, a nie uplastycznienie.

Połączenie cierne ma nominalną sumę sił zacisku $F_c=12{,}0 \mathrm{kN}$, zmierzony współczynnik tarcia $\mu=0{,}22\pm0{,}01$ i efektywny promień $r_e=0{,}070\pm0{,}002 \mathrm{m}$; oba zapisy $\pm$ oznaczają syntetyczne niepewności standardowe. Nośność centralna przed poślizgiem wynosi $\mu F_c=2{,}64 \mathrm{kN}$. Przy centralnym $r_e$ moment zawiasowy odpowiada składowej $H_{lim}/r_e=0{,}691 \mathrm{kN}$, a zredukowane żądanie

$$L_{slip}=\sqrt{(0{,}993 \mathrm{kN})^2+(0{,}691 \mathrm{kN})^2}=1{,}21 \mathrm{kN}.$$

Centralny margines poślizgu jest wysoki: $MS_{slip,cent}=2{,}64/1{,}21-1=1{,}18$. Regułę $2u$ stosujemy jednak również do tarcia i ramienia: $\mu^-=0{,}22-2(0{,}01)=0{,}20$, $r_e^-=0{,}070-2(0{,}002)=0{,}066 \mathrm m$. Prostokątna, jednoczesna obwiednia bierze $\mu^-$ po stronie nośności i $r_e^-$ po stronie popytu, bez przypisywania jej dokładnego prawdopodobieństwa. Daje ona

$$L_{slip}^{+}=\sqrt{0{,}993^2+(48{,}4/66)^2}=1{,}234 \mathrm{kN},\qquad R_{slip}^{-}=0{,}20(12{,}0)=2{,}400 \mathrm{kN},$$

czyli $MS_{slip}^{cov}=2{,}400/1{,}234-1=0{,}944$. Po utracie połowy centralnego napięcia wstępnego pokryciowa nośność wynosi już $0{,}20(6{,}0)=1{,}200 \mathrm{kN}$, a margines $1{,}200/1{,}234-1=-0{,}028$. Uwzględnienie niepewności nie osłabia zatem konserwatyzmu i ujawnia mechanizm, który analiza samej powłoki pominęłaby.

Sygnatury i pomiar. H1 powinna podnieść równocześnie siłę z wagi aerodynamicznej, ciśnienia na powierzchni oraz globalne przyspieszenie. H2 zmieni moment względem CG bez porównywalnej zmiany momentu zawiasowego. H3 podniesie prąd lub siłę aktuatora i odkształcenie lokalne przy zgodnej sile normalnej. H4 da histerezę przemieszczenia ramy, redystrybucję odkształceń wokół łączników i spadek zmierzonej siły zacisku, przy możliwie poprawnych wielkościach aerodynamicznych.

Plan rozstrzygający łączy próbę aerodynamiczną reprezentatywnego kształtu z wagą sześcioskładową, pomiar CG metodą dwóch podparć, kalibrowany tor momentu aktuatora, tensometry na ramie i powłoce, znaczniki przemieszczenia złącza oraz bezpośredni pomiar napięcia łączników przed i po cyklu cieplnym. Kanały synchronizuje wspólna podstawa czasu. Kalibracje obejmują krzyżową czułość, temperaturę i powtarzalność montażu. Budżet wyniku siły ma $57{,}7 \mathrm{N}$ niepewności modelowej; pomiar siły ma niezależnie $25 \mathrm{N}$, CG $3 \mathrm{mm}$, moment zawiasowy $2 \mathrm{N\,m}$, a siła zacisku $0{,}6 \mathrm{kN}$, wszystkie jako standardowe niepewności syntetyczne.

Wynik i ponowny rachunek. W fikcyjnej serii uzyskano $N=806\pm25 \mathrm{N}$ wobec predykcji $788\pm58 \mathrm{N}$, $s_{CG}=0{,}868\pm0{,}003 \mathrm{m}$ wobec $0{,}860\pm0{,}010 \mathrm{m}$, moment zawiasowy $46\pm2 \mathrm{N\,m}$ wobec $34\pm5 \mathrm{N\,m}$ oraz siłę zacisku po cyklu $6{,}4\pm0{,}6 \mathrm{kN}$ wobec nominalnych $12{,}0 \mathrm{kN}$. Każde $\pm$ w tym zdaniu oznacza niepewność standardową. Nie wykonano pomiaru siły zacisku przed cyklem, dlatego porównanie z wartością nominalną nie jest pomiarem spadku przed–po. Lokalna deformacja była o $28\%$ większa od pierwotnej predykcji. Po wprowadzeniu zmierzonego momentu i podatności złącza model odtworzył ją w granicy tolerancji porównania $\pm6\%$. Nie wyznaczono jednak ilościowej miary pętli histerezy ani jej niepewności. H1 i H2 nie mają sygnatury większej od połączonej niepewności; H3 jest potwierdzona, natomiast niższa od nominalnej siła zacisku po cyklu i większa deformacja jedynie wspierają H4. Do jej rozstrzygnięcia potrzeba pomiaru przed–po z kowariancją wyników powtórzonych oraz ilościowej sygnatury histerezy.

Reguła pokrycia pozostaje taka jak przed próbą: dla popytu bierzemy jednostronną górną granicę $\hat y+2u_y$ i mnożnik $k_d=1{,}10$, a dla nośności jednostronną dolną granicę $\hat y-2u_y$. Wartość $2u$ jest tu programową regułą pokrycia, nie twierdzeniem o dokładnym prawdopodobieństwie bez znajomości rozkładu. Dla niezależnych wielkości obliczanych stosujemy propagację RSS do wyniku, natomiast dla granicznego sprawdzenia poślizgu — konserwatywną prostokątną obwiednię równoczesnych granic wejść. Z rozkładu ciśnienia i zmierzonego $s_{CG}$ otrzymano ramię $l_N=0{,}550\pm0{,}010 \mathrm m$ z niepewnością standardową. Przy niezależności $N$ i $l_N$:

$$\hat M=806(0{,}550)=443{,}3 \mathrm{N\,m},\qquad u_M=443{,}3\sqrt{(25/806)^2+(0{,}010/0{,}550)^2}=15{,}94 \mathrm{N\,m}.$$

Zaktualizowane obciążenia graniczne wynoszą

$$N_{lim}^{post}=1{,}10(806+2\cdot25)=941{,}6 \mathrm N,$$

$$M_{lim}^{post}=1{,}10(443{,}3+2\cdot15{,}94)=522{,}7 \mathrm{N\,m},\qquad H_{lim}^{post}=1{,}10(46+2\cdot2)=55{,}0 \mathrm{N\,m}.$$

Reakcje w podziale $55/45$ to $518$ i $424 \mathrm N$. Naprężenie belkowe jest równe $522{,}7/(3{,}60\cdot10^{-6})=145{,}2 \mathrm{MPa}$. Skorelowany model lokalny reprezentujemy jawnym mnożnikiem $k_{loc}=1{,}06$, zatem $\sigma_{VM}^{post}=1{,}06(145{,}2)=153{,}9 \mathrm{MPa}$. Marginesy po pomiarze pozostają dodatnie:

$$MS_y=\frac{240}{153{,}9}-1=0{,}56,\qquad MS_b=\frac{190}{145{,}2}-1=0{,}31,\qquad MS_u=\frac{330}{1{,}50(153{,}9)}-1=0{,}43.$$

Kryterium odrzucające ujawnia dopiero poślizg. Dla wartości centralnych oraz pokryciowej obwiedni, która używa już wyznaczonych $\mu^-=0{,}20$ i $r_e^-=0{,}066 \mathrm m$, popyt i nośności cierne wynoszą

$$L_{slip,cent}^{post}=\sqrt{(0{,}9416)^2+(55{,}0/70)^2}=1{,}226 \mathrm{kN},\qquad L_{slip,+}^{post}=\sqrt{(0{,}9416)^2+(55{,}0/66)^2}=1{,}257 \mathrm{kN},$$

$$R_{slip,cent}=0{,}22(6{,}4)=1{,}408 \mathrm{kN},\qquad R_{slip,-2u}=[0{,}22-2(0{,}01)][6{,}4-2(0{,}6)]=1{,}040 \mathrm{kN}.$$

Stąd $MS_{slip,cent}=1{,}408/1{,}226-1=+0{,}148$, lecz margines według przyjętej reguły pokrycia wynosi $MS_{slip,-2u}=1{,}040/1{,}257-1=-0{,}173$. Konfiguracja procesu połączenia A-07 zostaje więc odrzucona przez ilościowe kryterium poślizgu, mimo dodatnich marginesów powłoki i dodatniego marginesu przy wartościach centralnych. Należy wykazać retencję napięcia wstępnego; równocześnie baza momentu zawiasowego otrzymuje skorygowaną obwiednię. Decyzja dotyczy tego materiału, powierzchni styku, procedury montażu i badanego zakresu temperatur, nie całej klasy płatowców. Regresja po zmianie obejmuje pomiar napięcia po cyklach środowiskowych, powtórzenie obciążenia do wartości granicznej i niszczącej na egzemplarzu kwalifikacyjnym, korelację tensometryczną drogi obciążenia oraz ponowne przeliczenie lokalnych punktów 6-DOF dla trzech stanów masowych. Nie wolno zamknąć niezgodności samym podniesieniem współczynnika bezpieczeństwa: zmieniła się fizyka interfejsu.

Gdy płatowiec przestaje być bryłą sztywną

Obliczenia statyczne odpowiadają na pytanie o nośność i ugięcie pod powoli zmieniającym się obciążeniem. Nie odpowiadają, czy obciążenie, sterowanie albo drganie wymuszone pobudzi postać drgań konstrukcji. W liniowym modelu elementów skończonych (MES) ruch opisuje równanie

$$M_s\ddot{\boldsymbol u}+C_s\dot{\boldsymbol u}+K_s\boldsymbol u =B_a\boldsymbol f_a+B_c\boldsymbol f_c,$$

gdzie $M_s,C_s,K_s$ są macierzami masy, tłumienia i sztywności, $\boldsymbol u$ wektorem przemieszczeń węzłowych, $B_a\boldsymbol f_a$ rozłożonym wymuszeniem aerodynamicznym, a $B_c\boldsymbol f_c$ siłami wprowadzenia przez aktuator i wyposażenie. Przy redukcji modalnej $\boldsymbol u=\Phi\boldsymbol\eta$, gdzie kolumnami macierzy modalnej $\Phi$ są postacie drgań $\boldsymbol\phi_j$, normalizacji postaci do masy jednostkowej i pominięciu sprzężenia tłumienia otrzymujemy dla postaci $j$

$$\ddot\eta_j+2\zeta_j\omega_j\dot\eta_j+\omega_j^2\eta_j =\boldsymbol\phi_j^T B_a\boldsymbol f_a+\boldsymbol\phi_j^T B_c\boldsymbol f_c.$$

$\eta_j$ jest współrzędną modalną, $\omega_j=2\pi f_j$ częstością kołową, $f_j$ częstotliwością w hercach, $\zeta_j$ bezwymiarowym współczynnikiem tłumienia, a $\boldsymbol\phi_j$ postacią drgań. Wartość $\boldsymbol\phi_j$ w miejscu sensora określa obserwowalność postaci; wartość i kierunek w miejscu aktuatora — jej sterowalność. Sensor w węźle drgań może nie zobaczyć groźnej postaci drgań, a sensor przy brzuchu drgań może raportować głównie lokalne ugięcie zamiast ruchu bryły. Podobnie dwa aktuatory o tej samej sile mogą bardzo różnie pobudzać konstrukcję.

Transmitancja częstotliwościowa musi odpowiadać rzeczywistej wielkości mierzonej. W rozpatrywanym kanale jednoskładowy żyroskop mierzy lokalną prędkość kątową $y_g$ w $\mathrm{rad/s}$. Dla zginania w płaszczyźnie pomiaru

$$y_g(t)=q_b(t)+\sum_{j=1}^{n}\phi'_j(s_g)\dot\eta_j(t),$$

gdzie $s_g$ jest stacją żyroskopu, kreska oznacza pochodną postaci drgań po stacji, a $q_b$ prędkość kątową bryły. Jeżeli rzeczywiste wychylenie aktuatora $\delta$ pobudza postać ze współczynnikiem uczestnictwa $\Gamma_j$, transmitancja od $\delta$ w radianach do $y_g$ ma postać

$$P_g(s)=P_{rig}(s)+\sum_{j=1}^{n} \frac{sR_j}{s^2+2\zeta_j\omega_j s+\omega_j^2}, \qquad R_j=\phi'_j(s_g)\Gamma_j.$$

Zmienna Laplace'a $s$ ma jednostkę $\mathrm{s^{-1}}$; czynnik $s$ wynika z pomiaru prędkości, a znak i wartość $R_j$ wynikają z położenia oraz kierunku aktuatora i sensora. Aktuator, filtr sensora i opóźnienie nie są ukryte w stałym residuum. Jest to odpowiedź częstotliwościowa FRF (frequency-response function), nie lista samych częstotliwości własnych. Dwa modele o zgodnych $f_j$ mogą dać przeciwne skutki w pętli, jeśli różnią się residuum, tłumieniem albo położeniem zera antyrezonansu.

Aerodynamiczne obciążenie zależne od ruchu dodaje efektywną sztywność i tłumienie. W formie liniowej

$$M_s\ddot{\boldsymbol u}+\left(C_s-C_{aero}(q_\infty,M)\right)\dot{\boldsymbol u} +\left(K_s-K_{aero}(q_\infty,M)\right)\boldsymbol u =B_c\boldsymbol f_c,$$

gdzie $C_{aero}$ i $K_{aero}$ oznaczają macierze uogólnionego tłumienia i sztywności aerodynamicznej. Znak zależy od konwencji i nie może być odgadywany. Utrata dodatniego tłumienia efektywnego prowadzi do niestateczności aeroelastycznej; sprzężenie z sensorem, aktuacją, filtracją i prawem sterowania tworzy problem aeroservoelastyczny. Model samolotu elastycznego Meirovitcha i Tuzcu pokazuje wspólną formulację ruchu bryły, postaci drgań, aerodynamiki i sterowania. Jest źródłem metodycznym, lecz dotyczy manewrującego samolotu; dla pocisku stanowi analogię pośrednią, wymagającą własnych warunków brzegowych, geometrii i dowodu.⁶

Próba drgań naziemnych GVT (ground vibration test) dostarcza częstotliwości, tłumienia, postaci drgań i FRF w kontrolowanej konfiguracji. Nie jest „pomiarowym modelem lotu”: zawieszenie, przewody, imitatory masy, temperatura i brak obciążenia aerodynamicznego zmieniają układ. Celem korelacji nie jest dopasowanie jednej częstotliwości dowolną zmianą modułu materiału, lecz znalezienie fizycznie uzasadnionych parametrów, które jednocześnie poprawiają kilka postaci drgań, kilka FRF i przewidywanie innych konfiguracji. NASA-STD-5002B wskazuje, że tłumienie ma opierać się na badaniach podobnych konstrukcji, a krytyczne modele dynamiczne podlegają walidacji i korelacji; znów wykorzystujemy logikę procesu, nie przenosimy automatycznie kryteriów statku kosmicznego.¹

Próba sprzętowa w pętli HIL (hardware-in-the-loop) zamyka rzeczywisty sensor, komputer i aktuator wokół modelu czasu rzeczywistego. Musi odtwarzać opóźnienia, częstotliwość próbkowania, nasycenia, szumy i elastyczne FRF w zakresie ważnym dla decyzji. HIL nie waliduje modelu modalnego, jeśli ten sam niezweryfikowany model generuje bodziec i stanowi wzorzec odpowiedzi. Łańcuch dowodu biegnie od analizy modalnej przez GVT i korelację do HIL, a nie na skróty.

Studium II: postać drgań zginających w kanale żyroskopowym

Wymaganie i stan początkowy. Dla wszystkich trzech stanów masowych OR-D1 i zdefiniowanego zakresu temperatur kanał sterowania ma zachować zapas fazy co najmniej $35^\circ$ i zapas wzmocnienia co najmniej $6 \mathrm{dB}$ po uwzględnieniu elastyczności, opóźnienia i niepewności. W konfiguracji A-07 model MES przewidywał pierwszą istotną postać drgań zginających $f_{MES}=49{,}0 \mathrm{Hz}$ z $\zeta_{MES}=0{,}018$. GVT ujawniła $f_{GVT}=42{,}3\pm0{,}4 \mathrm{Hz}$ i $\zeta_{GVT}=0{,}008\pm0{,}001$; oba zapisy $\pm$ oznaczają niepewności standardowe. Są to dane szkoleniowe, niezwiązane z pasmem żadnego systemu.

Rozważamy H1 — niezamodelowaną podatność połączenia ramy z powłoką, spójną z wynikiem pierwszego studium; H2 — lokalną postać drgań uchwytu żyroskopu omyłkowo rozpoznaną jako zginanie korpusu; H3 — zbyt sztywne albo zbyt miękkie podparcie GVT; H4 — aliasing, czyli odwzorowanie wyższej częstotliwości poniżej połowy częstotliwości próbkowania. Każda hipoteza ma inną sygnaturę. H1 zmienia globalną krzywiznę i FRF między odległymi punktami oraz reaguje na zmianę sztywności złącza. H2 ma duże lokalne przyspieszenie uchwytu, małą energię globalną i niewielki udział w odległych punktach. H3 przesuwa częstotliwość po zmianie zawieszenia. H4 znika lub przemieszcza się po zmianie próbkowania i filtru przeciwaliasingowego.

Do porównania postaci drgań używamy kryterium MAC (modal assurance criterion)

$$MAC(\boldsymbol\phi_t,\boldsymbol\phi_a)= \frac{|\boldsymbol\phi_t^T\boldsymbol\phi_a|^2} {(\boldsymbol\phi_t^T\boldsymbol\phi_t)(\boldsymbol\phi_a^T\boldsymbol\phi_a)},$$

gdzie $\boldsymbol\phi_t$ jest zmierzoną, a $\boldsymbol\phi_a$ obliczoną postacią drgań w tych samych punktach i kierunkach. MAC bliski jedności wskazuje podobieństwo kształtu, ale nie dowodzi poprawnego tłumienia ani residuum. Plan GVT obejmuje wielopunktowe wymuszenie o małej amplitudzie, akcelerometry na korpusie i uchwycie, bezkontaktowy pomiar ruchu wybranych punktów, dwa warianty zawieszenia oraz rejestrację siły wymuszającej. Tor podstawowy próbkuje z $1{,}0 \mathrm{kHz}$ i ma analogowy filtr przed próbkowaniem; równoległy tor diagnostyczny rejestruje z $5{,}0 \mathrm{kHz}$. Niepewność częstotliwości $0{,}4 \mathrm{Hz}$ obejmuje rozdzielczość estymatora, powtarzalność montażu i zmiany temperatury; niepewność tłumienia $0{,}001$ obejmuje rozrzut między metodami połowy mocy i dopasowania zaniku.

Wynik ma następującą strukturę dowodową. Zmierzona globalna postać ma $MAC=0{,}91$ z pierwszą postacią drgań zginających MES, podczas gdy lokalny uchwyt wykazuje osobny szczyt $118\pm2 \mathrm{Hz}$, gdzie $\pm2 \mathrm{Hz}$ jest niepewnością standardową. Zmiana zawieszenia przesuwa badaną postać z $42{,}3\pm0{,}4 \mathrm{Hz}$ do $42{,}0\pm0{,}3 \mathrm{Hz}$; oba zapisy $\pm$ oznaczają niepewności standardowe. Dla niezależnych wyników $u_\Delta=\sqrt{0{,}4^2+0{,}3^2}=0{,}5 \mathrm{Hz}$, a rozszerzona niepewność różnicy wynosi $U=ku_\Delta=2(0{,}5)=1{,}0 \mathrm{Hz}$. Zaobserwowana różnica $0{,}3 \mathrm{Hz}$ jest od niej mniejsza, co oznacza brak wykrytej sygnatury H3 w zakresie dwóch zbadanych zawieszeń, ale nie dowodzi równoważności ani nie wyklucza wpływu podparcia poza tym zakresem. Do odrzucenia H3 potrzebne byłyby jawne właściwości obu zawieszeń, przewidywana czułość częstotliwości oraz uprzednio ustalona granica równoważności. W torze $5 \mathrm{kHz}$ nie ma szczytu przy $957{,}7 \mathrm{Hz}$, który mógłby aliasować do $42{,}3 \mathrm{Hz}$ w torze $1 \mathrm{kHz}$; dodatkowo faza i koherencja obu torów są zgodne. H2 i H4 zostają odrzucone, natomiast H3 pozostaje nierozstrzygnięta. Po zastąpieniu badanego interfejsu sztywnym imitatorem częstotliwość rośnie do $47{,}6\pm0{,}5 \mathrm{Hz}$, gdzie $0{,}5 \mathrm{Hz}$ jest niepewnością standardową. Aktualizacja sztywności wyłącznie fizycznych elementów złącza daje $42{,}8 \mathrm{Hz}$ i poprawia MAC pozostałych dwóch postaci drgań, dlatego H1 zostaje przyjęta jako najlepiej wsparta hipoteza, bez twierdzenia, że jest przyczyną wyłączną.

Wpływ na sterowanie oceniamy na jawnej, lokalnej pętli szkoleniowej w paśmie $0{,}01\leq f\leq100 \mathrm{Hz}$. Nie jest to regulator ani pasmo realnego systemu. Dla $s$ w $\mathrm{s^{-1}}$ przyjmujemy

$$\begin{aligned} C(s)&=0{,}76\frac{1+s/(2\pi\,6{,}7)}{s[1+s/(2\pi\,267)]}\quad[\mathrm s],\\ A(s)&=\frac{1}{1+s/(2\pi\,194)},\qquad H_g(s)=\frac{1}{1+s/(2\pi\,281)},\\ D(s)&=e^{-0{,}00162s},\\ P_g(s;f_n,\zeta)&=\frac{80 \mathrm{s^{-1}}}{1+s/(2\pi\,8{,}7)} -\frac{(220 \mathrm{s^{-2}})s}{s^2+2\zeta(2\pi f_n)s+(2\pi f_n)^2},\\ L(s)&=C(s)A(s)P_g(s;f_n,\zeta)H_g(s)D(s). \end{aligned}$$

$C$ przekształca błąd prędkości kątowej w komendę wychylenia, $A$ jest dynamiką aktuatora, $P_g$ prowadzi od rzeczywistego wychylenia do prędkości kątowej żyroskopu, $H_g$ jest filtrem sensora, a $D$ opóźnieniem. Iloczyn $L$ jest bezwymiarowy. Ujemne residuum $-220 \mathrm{s^{-2}}$ koduje przeciwny znak udziału postaci drgań w miejscu sensora i aktuatora. Przed korelacją podstawiamy $(f_n,\zeta)=(49{,}0 \mathrm{Hz},0{,}018)$, po korelacji $(42{,}3 \mathrm{Hz},0{,}008)$; jest to jawne przekazanie wyniku GVT do pętli.

Przecięcia wyznaczamy z $|L(j2\pi f_c)|=1$ oraz $\arg L(j2\pi f_\pi)=-180^\circ$, z fazą rozwiniętą w sposób ciągły. Dla każdego przecięcia fazy liczymy $GM=-20\log_{10}|L|$ i przyjmujemy najmniejszy zapas w badanym paśmie. Wyniki, otrzymane bez dodatkowego strojenia, są następujące:

Model $P_g$	Przecięcia $|L|=1$: $f_c$ i $PM$	Przecięcia fazy $-180^\circ$: $f_\pi$ i $GM$	Wynik
MES, $49{,}0 \mathrm{Hz}$	$11{,}601 \mathrm{Hz}$; $81{,}69^\circ$	$48{,}231 \mathrm{Hz}$; $7{,}23 \mathrm{dB}$; $50{,}803 \mathrm{Hz}$; $20{,}35 \mathrm{dB}$; $72{,}160 \mathrm{Hz}$; $17{,}21 \mathrm{dB}$	oba kryteria spełnione nominalnie
GVT, $42{,}3 \mathrm{Hz}$	$11{,}616 \mathrm{Hz}$; $81{,}62^\circ$	$41{,}944 \mathrm{Hz}$; $1{,}76 \mathrm{dB}$; $43{,}503 \mathrm{Hz}$; $29{,}54 \mathrm{dB}$; $72{,}210 \mathrm{Hz}$; $17{,}05 \mathrm{dB}$	kryterium $GM\geq6 \mathrm{dB}$ niespełnione

Kontrola faz i jednostek w pierwszym przecięciu wzmocnienia daje dla modelu MES kolejno: $|C|=0{,}020828 \mathrm s$, $|A|=0{,}99822$, $|P_g|=48{,}1388 \mathrm{s^{-1}}$, $|H_g|=0{,}99915$, $|D|=1$, a fazy $-32{,}50^\circ$, $-3{,}42^\circ$, $-53{,}26^\circ$, $-2{,}36^\circ$ i $-6{,}77^\circ$. Iloczyn modułów wynosi $1{,}0000$, suma faz $-98{,}31^\circ$, zatem $PM=180^\circ-98{,}31^\circ=81{,}69^\circ$. Przy $48{,}231 \mathrm{Hz}$ faza wynosi $-180^\circ$, a $20\log_{10}|L|=-7{,}23 \mathrm{dB}$, więc najmniejszy $GM=7{,}23 \mathrm{dB}$. Po korelacji analogiczne pierwiastki dają $PM=81{,}62^\circ$ i najmniejszy $GM=1{,}76 \mathrm{dB}$. Zmiana postaci drgań prawie nie przesuwa niskiego $f_c$, lecz podnosi rezonansowy moduł przy pierwszym przecięciu fazy; dlatego rozstrzyga zapas wzmocnienia, a nie dopisana ręcznie faza przy z góry ustalonej częstotliwości.

Dla danych GVT stosujemy dodatkowo prostokątną obwiednię $\pm2u$: $41{,}5\leq f_n\leq43{,}1 \mathrm{Hz}$ i $0{,}006\leq\zeta\leq0{,}010$. Symbol $\pm$ oznacza tu rozszerzenie standardowych niepewności współczynnikiem $k=2$, a nie tolerancję. Obliczenie wszystkich dziewięciu kombinacji środka i krańców daje dolne obwiednie $PM=-71{,}82^\circ$ oraz $GM=-0{,}05 \mathrm{dB}$; ujemne wartości pojawiają się wskutek dodatkowych przecięć w najniżej tłumionym wariancie. Nawet wynik centralny odrzuca konfigurację przez $GM=1{,}76 \mathrm{dB}$, więc niepewność nie służy do „zaliczenia przez zachodzenie przedziałów”. HIL odtwarza ponadto wzrost składowej około $42 \mathrm{Hz}$ w żyroskopie bez analogicznego ruchu modelu bryły, co jest rozróżniającą sygnaturą sprzężenia elastycznego.

Decyzja i regresja. Konfiguracja A-07 nie jest dopuszczona w badanym kanale: odrzuca ją centralny i obwiedniowy zapas wzmocnienia. Najsilniej wspartą przyczyną jest podatność połączenia, wskazana niezależnie przez GVT i pierwsze studium; lokalny uchwyt i aliasing nie mają oczekiwanych sygnatur, natomiast wpływu podparcia nie rozstrzygnięto formalnym testem równoważności. Decyzją projektową jest korekta fizycznego interfejsu i modelu oraz domknięcie testu H3, a nie samo dodanie filtru do nieskorelowanego obiektu. Ewentualna zmiana filtru stanowi odrębną zmianę GNC i również wymaga dowodu zapasów.

Po zmianie interfejsu regresja obejmuje GVT wszystkich stanów masowych i skrajnych temperatur, identyfikację częstotliwości, tłumienia, postaci drgań oraz pełnych FRF sensor–aktuator, ponowną korelację MES bez pogorszenia innych postaci, analizę aeroservoelastyczną przez całą obwiednię $q_\infty$–$M$, a następnie HIL z rzeczywistymi opóźnieniami, nasyceniami i filtracją. Ponowna kwalifikacja dotyczy struktury interfejsu, kanału pomiarowego oraz zintegrowanej pętli; nie obejmuje automatycznie innych konfiguracji sensora ani innej wersji oprogramowania. Warunkiem zamknięcia jest $PM\geq35^\circ$ i $GM\geq6 \mathrm{dB}$ dla każdego przecięcia w całej zadeklarowanej domenie i w obwiedni $\pm2u$, a nie tylko brak niestateczności w jednej próbie. Ostatni symbol $\pm$ oznacza tę samą obwiednię rozszerzoną z $k=2$.

Przejściowa termika jako obciążenie geometrii

Temperatura nie jest jedną liczbą przypisaną do „gorącej konfiguracji”. Aerotermiczny transfer energii, promieniowanie, przewodzenie przez styki, pojemność cieplna i czas zdarzenia tworzą pole $T(\boldsymbol r,t)$. To pole zmienia własności materiału, napięcie wstępne połączeń, luzy, geometrię aerodynamiczną, oś sensora i czasem częstotliwości modalne. W pierwszym studium po cyklu cieplnym zmierzono siłę zacisku niższą od wartości nominalnej, bez pomiaru odniesienia przed cyklem; w trzecim ta sama domena zostanie połączona z błędem pomiarowym.

Dla ciągłego ciała równanie przewodzenia ma postać

$$\rho_m c_p\frac{\partial T}{\partial t} =\nabla\cdot(k\nabla T)+\dot q_v,$$

gdzie $\rho_m$ jest gęstością materiału w $\mathrm{kg/m^3}$, $c_p$ ciepłem właściwym w $\mathrm{J/(kg\,K)}$, $k$ przewodnością w $\mathrm{W/(m\,K)}$, a $\dot q_v$ objętościowym źródłem ciepła w $\mathrm{W/m^3}$. Na powierzchni bilans może obejmować zadany strumień $\dot q''$, konwekcję i promieniowanie. W modelu skupionym węzła $i$

$$C_i\dot T_i=Q_i(t)+\sum_j G_{ij}(T_j-T_i) -h_iA_i(T_i-T_\infty)-\varepsilon_i\sigma_{SB} A_i(T_i^4-T_{rad}^4),$$

gdzie $C_i=m_ic_{p,i}$ jest pojemnością cieplną w $\mathrm{J/K}$, $Q_i$ zewnętrzną mocą cieplną w watach, $G_{ij}$ przewodnością połączenia w $\mathrm{W/K}$, $h_i$ współczynnikiem konwekcji, $A_i$ polem wymiany, $T_\infty$ temperaturą płynu, $T_{rad}$ temperaturą radiacyjną otoczenia, $\varepsilon_i$ emisyjnością, a $\sigma_{SB}$ stałą Stefana–Boltzmanna. NASA SP-8105 opisuje taki bilans węzłowy, znaczenie warunków przejściowych, gradientów i przewodności interfejsów. Dokument dotyczy kontroli cieplnej statków kosmicznych, więc wspiera równanie i logikę korelacji, nie wartość strumienia dla pocisku.⁸

Opór kontaktowy $R_c=1/G_c$ obejmuje mikroskopowe pola styku, warstwę pośrednią i zależność od docisku. Wartość z katalogu materiału nie zastępuje pomiaru złącza. Jeżeli $R_c$ rośnie wskutek utraty docisku, zewnętrzna osłona nagrzewa się szybciej względem ramy; równocześnie maleje nośność cierna. Jest to fizyczne sprzężenie pierwszego i trzeciego studium, a nie przypadkowe użycie tej samej temperatury.

Odkształcenie swobodne izotropowego materiału wynosi $\varepsilon_{th}=\alpha_T\Delta T$, gdzie $\alpha_T$ jest współczynnikiem rozszerzalności cieplnej w $\mathrm{K^{-1}}$. Gdy dwa połączone elementy mają różne $\alpha_T$, temperaturę lub długość, powstaje różnicowe przemieszczenie. Jeśli więzy mu przeciwdziałają, rodzą się naprężenia termosprężyste; jeśli połączenie jest podatne, powstaje obrót i błąd osiowania. Przykłady z wytrzymałości materiałów pokazują, że reakcja cieplna zależy przede wszystkim od zgodności odkształceń i warunków podparcia, nie od samego przyrostu temperatury.⁹

Pełny model termosprężysty ma równanie równowagi ze związkiem materiałowym

$$\boldsymbol\sigma=\mathsf C: \left(\boldsymbol\varepsilon-\alpha_T\Delta T\,\boldsymbol I_2\right),$$

gdzie $\mathsf C$ jest tensorem sprężystości, dwukropek oznacza podwójną kontrakcję tensorów, $\boldsymbol\varepsilon$ jest tensorem odkształcenia, a $\boldsymbol I_2$ tensorem jednostkowym drugiego rzędu. Klej albo przekładka mogą wymagać lepkosprężystości, skurczu utwardzania i zależności własności od wilgotności. Model liniowy jest wystarczający tylko tak długo, jak długo historia temperatury nie uruchamia pełzania, poślizgu albo zmiany stanu styku.

Studium III: gradient osłony sensora i dryf osi

Wymaganie i konfiguracja. W A-07 całkowity błąd osi sensora wywołany środowiskiem nie może przekroczyć $0{,}25 \mathrm{mrad}$ w dowolnej chwili zdefiniowanego zdarzenia cieplnego. Osłona i rama sensora tworzą dwa główne węzły termiczne. Historia przekazana z przypadku ruchu ma postać $q''_{avail}(t)=25 \mathrm{kW/m^2}$ dla $0<t\leq30 \mathrm s$ i zero poza tym odcinkiem. Symbol $q''_{avail}$ oznacza syntetyczny potencjalny strumień entalpii dostępny przy ścianie, przed uwzględnieniem przejmowania przez jej rzeczywistą warstwę przyścienną; nie jest strumieniem radiacyjnym ani strumieniem netto do materiału. Bezwymiarowy $\eta_h=0{,}72$ jest syntetycznym efektywnym współczynnikiem przejmowania i akomodacji dla tego interfejsu, a nie absorpcyjnością powierzchni. Strumień netto wchodzący do węzła wynosi z definicji modelu $q''_{wall}=\eta_hq''_{avail}$. Pole wejścia ma $A_h=0{,}0120 \mathrm{m^2}$, masa nagrzewanej części osłony to $0{,}120 \mathrm{kg}$, a $c_p=800 \mathrm{J/(kg\,K)}$. Pomiar maksymalnego błędu dał $0{,}31\pm0{,}03 \mathrm{mrad}$; $0{,}03 \mathrm{mrad}$ jest niepewnością standardową.

Pięć hipotez wyjaśnia rozbieżność: H1 — rzeczywisty strumień ciepła jest większy; H2 — opór kontaktowy osłona–rama jest większy od modelowego; H3 — różnica współczynników rozszerzalności połączenia jest błędna; H4 — klej wykazuje dodatkowy skurcz albo pełzanie zależne od historii; H5 — własny dryf temperaturowy sensora został błędnie przypisany konstrukcji. H1 zmienia tempo nagrzewania powierzchni; H2 zwiększa różnicę temperatur przy poprawnym strumieniu; H3 skaluje przemieszczenie niemal proporcjonalnie do gradientu; H4 daje histerezę po wyrównaniu temperatur; H5 występuje także podczas izotermicznego nagrzewania sensora bez gradientu konstrukcji.

Moc przekazana do węzła osłony wynosi

$$Q_{in}(t)=\eta_hq''_{avail}(t)A_h =0{,}72(25\,000 \mathrm{W/m^2})(0{,}0120 \mathrm{m^2})=216 \mathrm W$$

dla $0<t\leq30 \mathrm s$ i zero później. Wartość $216 \mathrm W$ wynika wyłącznie z powyższej syntetycznej definicji interfejsu OR-D1 i iloczynu $0{,}72\cdot25\,000\cdot0{,}0120$; nie jest wartością ani definicją zaczerpniętą z NASA SP-8105. Gdyby wejściem był już konwekcyjny strumień netto $q''_{wall}$, żadnego dodatkowego mnożnika nie wolno byłoby stosować.

Pojemność cieplna węzła jest równa

$$C_h=(0{,}120 \mathrm{kg})(800 \mathrm{J/(kg\,K)})=96 \mathrm{J/K}.$$

Zmierzony opór kontaktowy ma wartość $R_c=0{,}55 \mathrm{K/W}$, czyli $G_c=1/R_c=1{,}818 \mathrm{W/K}$. Zlinearyzowane straty od osłony do otoczenia wynoszą $G_{he}=1{,}00 \mathrm{W/K}$. Rama ma pojemność $C_r=120 \mathrm{J/K}$ i straty $G_{re}=0{,}50 \mathrm{W/K}$. Przy przyrostach temperatury $\vartheta_h=T_h-T_0$ i $\vartheta_r=T_r-T_0$ kompletny model dwóch węzłów brzmi

$$\begin{aligned} C_h\dot\vartheta_h&=Q_{in}(t)-G_c(\vartheta_h-\vartheta_r)-G_{he}\vartheta_h,\\ C_r\dot\vartheta_r&=G_c(\vartheta_h-\vartheta_r)-G_{re}\vartheta_r, \end{aligned} \qquad \vartheta_h(0)=\vartheta_r(0)=0.$$

Jest to układ $\dot{\boldsymbol\vartheta}=A_T\boldsymbol\vartheta+\boldsymbol b_TQ_{in}$, gdzie

$$A_T= \begin{bmatrix} -(G_c+G_{he})/C_h&G_c/C_h\\ G_c/C_r&-(G_c+G_{re})/C_r \end{bmatrix},\qquad \boldsymbol b_T=\begin{bmatrix}1/C_h\\0\end{bmatrix}.$$

Dla stałej mocy na odcinku odpowiedź jest równa $\boldsymbol\vartheta(t)=A_T^{-1}[e^{A_Tt}-I]\boldsymbol b_TQ_{in}$. Podstawienie wszystkich powyższych parametrów daje w $t=30 \mathrm s$

$$\vartheta_h=46{,}536 \mathrm K,\qquad \vartheta_r=9{,}777 \mathrm K,\qquad \Delta T_{hr}=36{,}758 \mathrm K.$$

Wynik nie powstaje przez odjęcie arbitralnego ogrzania ramy od odpowiedzi jednowęzłowej. Uproszczenie nadal nie opisuje przestrzennego maksimum w skórze ani nieliniowego promieniowania; model szczegółowy musi zachować rozkład temperatury na interfejsie.

Różnica efektywnych współczynników rozszerzalności wynosi $\Delta\alpha=8{,}5\cdot10^{-6} \mathrm{K^{-1}}$, a efektywna długość drogi deformacji $L_T=0{,}420 \mathrm m$. Różnicowe przemieszczenie jest więc równe

$$\Delta u=\Delta\alpha\,L_T\,\Delta T_{hr} =(8{,}5\cdot10^{-6} \mathrm{K^{-1}})(0{,}420 \mathrm m)(36{,}758 \mathrm K) =1{,}312\cdot10^{-4} \mathrm m=0{,}1312 \mathrm{mm}.$$

Przy bazie kinematycznej $b=0{,}700 \mathrm m$ mały obrót wynosi

$$\theta_{th}=\frac{\Delta u}{b} =1{,}875\cdot10^{-4} \mathrm{rad}=0{,}1875 \mathrm{mrad}.$$

Wymiary $L_T$ i $b$ są parametrami zredukowanego modelu dydaktycznego, nie kompletną geometrią wykonawczą.

Plan pomiaru rozdziela wejście, drogę cieplną, deformację i sensor. Stan strugi i model kalibracyjny wyznaczają $q''_{avail}$, kalorymetr wyznacza $q''_{wall}$, a ich iloraz identyfikuje $\eta_h$; termopary powierzchniowe kontrolują odpowiedź ściany. Pary czujników po obu stronach styku wraz z mierzoną mocą wyznaczają $R_c$; bezkontaktowy pomiar przemieszczeń ramy daje obrót mechaniczny; niezależny wzorzec optyczny rejestruje oś sensora. Próbka materiałowa i świadek tego samego procesu klejenia dostarczają $\Delta\alpha$ oraz historię skurczu. Osobna próba izotermiczna sensora identyfikuje H5. Wszystkie kanały mają wspólny czas, a czujniki temperatury są skorygowane o bezwładność.

W budżecie pomiarowym otrzymano $q''_{avail}=24{,}6\pm1{,}2 \mathrm{kW/m^2}$, $\eta_h=0{,}72\pm0{,}04$, $R_c=0{,}57\pm0{,}06 \mathrm{K/W}$, $\Delta T_{hr}=36{,}5\pm2{,}0 \mathrm K$, $\Delta\alpha=(8{,}5\pm0{,}7)\cdot10^{-6} \mathrm{K^{-1}}$ i własny dryf sensora $0{,}010\pm0{,}010 \mathrm{mrad}$. Każde $\pm$ w tym zdaniu oznacza niepewność standardową. Iloczyn centralnych wartości wejścia cieplnego daje $q''_{wall}=17{,}7 \mathrm{kW/m^2}$, zgodny z modelowym $18{,}0 \mathrm{kW/m^2}$; w mechanicznym budżecie osiowania używamy jednak bezpośrednio wynikowego gradientu, więc nie liczymy ponownie niepewności $q''_{avail}$ i $\eta_h$ jako niezależnych wkładów.

Wszystkie wejścia relacji $\theta_{th}=\Delta\alpha L_T\Delta T_{hr}/b$ są jawne. Dla centralnej prognozy $\Delta T_{hr}=36{,}758 \mathrm K$, a jej standardową niepewność po korelacji termicznej przyjmujemy równą $u_T=2{,}0 \mathrm K$; pomiar $36{,}5\pm2{,}0 \mathrm K$ jest jej niezależnym sprawdzeniem. Dalej $\Delta\alpha=8{,}5\cdot10^{-6} \mathrm{K^{-1}}$ z $u_\alpha=0{,}7\cdot10^{-6} \mathrm{K^{-1}}$, $L_T=0{,}420 \mathrm m$ z $u_L=0{,}006 \mathrm m$ oraz $b=0{,}700 \mathrm m$ z $u_b=0{,}005 \mathrm m$. Są to syntetyczne, standardowe niepewności estymat traktowanych w tym budżecie jako wycentrowane, w przybliżeniu normalne i wzajemnie niezależne; wszystkie kowariancje przyjmujemy więc równe zeru. Osobny addytywny składnik $u_{red}=0{,}0234 \mathrm{mrad}$ jest standardową niepewnością błędu redukcji przestrzennego modelu termosprężystego do ilorazu $\Delta u/b$, oszacowaną syntetycznie z rozrzutu reszt walidacyjnych. Nie jest błędem czujnika ani arbitralnym współczynnikiem bezpieczeństwa i także jest traktowany jako niezależny.

Ponieważ czułości względne wynoszą $\partial\ln\theta/\partial\ln\Delta T=\partial\ln\theta/\partial\ln\Delta\alpha=\partial\ln\theta/\partial\ln L_T=1$ oraz $\partial\ln\theta/\partial\ln b=-1$, pełny rachunek RSS przy $\theta_{th}=0{,}187468 \mathrm{mrad}$ daje

$$\begin{aligned} u_{th}^2={}&\theta_{th}^2\left[ \left(\frac{2{,}0}{36{,}758}\right)^2+ \left(\frac{0{,}7}{8{,}5}\right)^2+ \left(\frac{0{,}006}{0{,}420}\right)^2+ \left(\frac{0{,}005}{0{,}700}\right)^2\right]+u_{red}^2\\ ={}&(0{,}01020)^2+(0{,}01544)^2+(0{,}00268)^2+(0{,}00134)^2+(0{,}02340)^2 \mathrm{mrad^2},\\ u_{th}={}&0{,}02998 \mathrm{mrad}\approx0{,}030 \mathrm{mrad}. \end{aligned}$$

Ani H1, ani H2 nie mają wymaganej sygnatury: zmierzone dostępne i przekazane do ściany strumienie oraz gradient pokrywają się z modelem. H3 także mieści się w budżecie, a H5 jest zbyt mała.

Zsynchronizowane przebiegi świadka kleju i próby izotermicznej pokazują, że podczas impulsu dodatkowy obrót kleju można odtworzyć funkcją. Po wyłączeniu wprowadzamy jawne warunki przełączenia: temperatury są ciągłe, $Q_{in}$ skacze z $216$ do $0 \mathrm W$, trwały składnik klejowy zachowuje wartość końcową, a własny dryf sensora zanika z syntetyczną stałą czasową $60 \mathrm s$:

$$\theta_{glue}(t)= \begin{cases} 0{,}120\dfrac{1-e^{-t/7}}{1-e^{-30/7}} \mathrm{mrad},&0\leq t\leq30 \mathrm s,\\ 0{,}120 \mathrm{mrad},&t>30 \mathrm s, \end{cases} \qquad \theta_s(t)= \begin{cases} 0{,}010\dfrac{t}{30} \mathrm{mrad},&0\leq t\leq30 \mathrm s,\\ 0{,}010e^{-(t-30)/60} \mathrm{mrad},&t>30 \mathrm s. \end{cases}$$

Składnik klejowy pozostaje po impulsie; jego maksimum ma wartość $0{,}120\pm0{,}015 \mathrm{mrad}$ z niepewnością standardową. Dla $0\leq t\leq30 \mathrm s$ gradient oraz oba pozostałe składniki rosną monotonicznie i mają zmierzony ten sam znak. Bezpośrednio przed przełączeniem $d\Delta T_{hr}/dt=+0{,}55286 \mathrm{K/s}$, a bezpośrednio po nim $-1{,}69714 \mathrm{K/s}$, więc gradient ma maksimum lokalne dokładnie w $30 \mathrm s$.

Przebieg chłodzenia wynika z tego samego modelu, bez nowego dopasowania. Dla $\tau=t-30\geq0$ zachodzi $\boldsymbol\vartheta(30+\tau)=e^{A_T\tau}\boldsymbol\vartheta(30)$, a po diagonalizacji

$$\Delta T_{hr}(30+\tau)= 41{,}0913e^{-0{,}0420049\tau}-4{,}33286e^{-0{,}00666935\tau} \mathrm K.$$

Daje to $29{,}116 \mathrm K$ w $35 \mathrm s$, $22{,}944 \mathrm K$ w $40 \mathrm s$ i $8{,}107 \mathrm K$ w $60 \mathrm s$. Jedyny późniejszy punkt stacjonarny gradientu jest minimum: $-1{,}684 \mathrm K$ dla $\tau=115{,}743 \mathrm s$; następnie gradient dąży od dołu do zera. Zatem dla całego chłodzenia $\Delta T_{hr}(t)\leq\Delta T_{hr}(30)$, $\theta_{glue}(t)\leq0{,}120 \mathrm{mrad}$ i $\theta_s(t)\leq0{,}010 \mathrm{mrad}$. Co więcej, najmniejszy wkład termomechaniczny to tylko $-0{,}00859 \mathrm{mrad}$, więc całkowity błąd pozostaje dodatni. Globalne maksimum także co do wartości bezwzględnej występuje zatem rzeczywiście w chwili przełączenia:

$$\theta_{tot}(30)=0{,}187468+0{,}120+0{,}010=0{,}317468 \mathrm{mrad}.$$

Standardowa niepewność tej sumy, przy niezależności trzech składników, wynosi $u_{tot}=\sqrt{0{,}02998^2+0{,}015^2+0{,}010^2}=0{,}03498 \mathrm{mrad}\approx0{,}0350 \mathrm{mrad}$. Predykcja $0{,}3175\pm0{,}0350 \mathrm{mrad}$ i pomiar $0{,}3100\pm0{,}0300 \mathrm{mrad}$, gdzie oba zapisy $\pm$ oznaczają niepewności standardowe, różnią się o $0{,}0075 \mathrm{mrad}$. Przy niezależności predykcji i pomiaru standardowa niepewność różnicy wynosi $u_\Delta=\sqrt{0{,}03498^2+0{,}0300^2}=0{,}04609 \mathrm{mrad}$. H4 zostaje przyjęta na podstawie bezpośredniej, czasowo zgodnej sygnatury: brakującą fizyką jest zależna od procesu odpowiedź warstwy klejowej.

Decyzja i zakres. Konfiguracja nie spełnia wymagania $0{,}25 \mathrm{mrad}$ i zostaje odrzucona do korekty procesu albo konstrukcji połączenia. Nie wolno zaliczyć jej przez podobieństwo do osłony o tym samym materiale zewnętrznym, ponieważ wynik kontroluje interfejs: grubość i stan kleju, docisk, utwardzanie, wilgotność oraz historia cieplna. Podobieństwo materiału osłony nie ustanawia podobieństwa drogi odkształcenia.

Regresja po zmianie obejmuje ponowny impuls termiczny z pomiarem $q''_{avail}(t)$, $q''_{wall}(t)$ i pełnego pola temperatury, identyfikację $\eta_h$ oraz oporu kontaktowego, szczelność przed i po cyklu, metrologię osi na gorąco i po ostygnięciu oraz próbę izotermiczną sensora. Do lokalnego modelu 6-DOF przekazuje się cały zsynchronizowany przebieg $y_{therm}(t)=\theta_{tot}(t)$, a nie samo maksimum; osobne przebiegi $q''_{avail}(t)$ i $\eta_h(t)$ tworzą wejście modelu cieplnego przez $q''_{wall}=\eta_hq''_{avail}$. Należy również powtórzyć kontrolę napięcia wstępnego z pierwszego studium i GVT z drugiego, ponieważ nowy interfejs zmienia jednocześnie przewodność, sztywność i tłumienie. Monografia pokazuje format tego zwrotnego przekazania, ale nie wykonuje ponownej całej trajektorii. Decyzja zachowuje ważność wyłącznie dla zakwalifikowanego procesu, zakresu temperatur, liczby cykli i partii materiałów objętych dowodem.

Historia egzemplarza: zmęczenie, wady i tolerancja uszkodzeń

Trzy studia analizowały krótkie zdarzenia, lecz płatowiec dziedziczy skutki wcześniejszej produkcji, transportu, przechowywania i prób. Naprężenie mniejsze od granicy plastyczności nie gwarantuje nieograniczonej trwałości. Zmienne obciążenie może inicjować pęknięcie w karbie, frettingu albo wadzie wykonania, a środowisko może przyspieszać propagację. Dlatego numer seryjny egzemplarza powinien wskazywać nie tylko rysunek wykonawczy, lecz również partię materiału, parametry procesów specjalnych, wyniki badań nieniszczących, liczbę i amplitudę cykli, przekroczenia oraz naprawy.

Cykl naprężenia opisują

$$\sigma_a=\frac{\sigma_{max}-\sigma_{min}}{2},\qquad \sigma_m=\frac{\sigma_{max}+\sigma_{min}}{2},\qquad R_\sigma=\frac{\sigma_{min}}{\sigma_{max}},$$

gdzie $\sigma_a$ jest amplitudą, $\sigma_m$ naprężeniem średnim, a $R_\sigma$ współczynnikiem asymetrii. Krzywa Wöhlera $S$–$N$ wiąże amplitudę z liczbą cykli do zniszczenia przy określonym materiale, wykończeniu, środowisku i stosunku $R_\sigma$. Przeniesienie krzywej z gładkiej próbki na spaw, otwór albo połączenie klejone wymaga odpowiedniej kategorii szczegółu, wpływu karbu i rozrzutu.

Dla widma o blokach amplitud prosty liniowy model kumulacji Palmgrena–Minera ma postać

$$D=\sum_i\frac{n_i}{N_i},$$

gdzie $n_i$ jest liczbą zrealizowanych cykli bloku $i$, a $N_i$ liczbą cykli do zniszczenia w warunkach tego bloku. $D=1$ nie jest prawem fizycznym. Model pomija kolejność obciążeń, przeciążenia, zamykanie pęknięcia i interakcje środowiskowe; jego próg decyzyjny oraz rozrzut muszą być wykazane dla danego szczegółu. Lokalne widmo nie może być liczone wyłącznie z komend sterowania. Musi pochodzić z reakcji konstrukcji, zatem zawierać przefiltrowanie modalne, drgania transportowe, cykle termiczne i zdarzenia obsługowe.

W OR-D1 jawne przekazanie do trwałości ma postać rekordu HIST-A07-03: historia reakcji przekrojowych z modelu obciążeń, w tym przypadek LC-A07-184, jest zamieniana przez skorelowany model lokalny na historię $\sigma(t)$, a zliczanie cykli tworzy trzy syntetyczne bloki. Dolne wartości $N_i^-$ są programowymi granicami pokrycia krzywych $S$--$N$ dla fikcyjnego szczegółu, nie danymi materiału rzeczywistego:

Blok historii bieżącej i planowanej	$\sigma_a$	$n_i$	$N_i^-$	$n_i/N_i^-$
transport i obsługa	$28 \mathrm{MPa}$	$12\,000$	$2{,}40\cdot10^6$	$0{,}0050$
odpowiedź manewrowa z uwzględnieniem postaci drgań	$74 \mathrm{MPa}$	$1\,600$	$1{,}60\cdot10^5$	$0{,}0100$
cykl termiczny interfejsu	$42 \mathrm{MPa}$	$400$	$8{,}00\cdot10^5$	$0{,}0005$

Stąd konserwatywna suma wynosi $D^+=0{,}0050+0{,}0100+0{,}0005=0{,}0155$. Dla ustanowionego wyłącznie w ćwiczeniu progu $D_{dop}=0{,}20$ warunkowy margines wynosi $MS_D=0{,}20/0{,}0155-1=11{,}90$. Liczba ta nie przywraca zgodności A-07. Studium I wykazało siłę zacisku po cyklu niższą od wartości nominalnej i deformację większą od predykcji; bez pomiaru odniesienia przed cyklem oraz ilościowej pętli histerezy nie dowodzi jednak spadku przed–po ani histerezy złącza. Nie zmierzyło też mikropoślizgu ani frettingu. Obliczony margines pokryciowy $MS_{slip,-2u}=-0{,}173$ wskazuje jednak możliwość poślizgu, a w dalszej eksploatacji potencjalnie frettingu. Kategoria szczegółu, na której oparto $N_i^-$, zakłada niezmienioną drogę obciążenia bez takiego mechanizmu, dlatego jej ważności nie wolno przyjąć dla obecnego złącza. Ilościowa decyzja trwałości brzmi więc: rachunek nominalnego widma ma dodatni margines, ale pozostaje warunkowo nieważny aż do zmiany połączenia i walidacji kategorii szczegółu w reprezentatywnym stanie styku. Po zmianie trzeba ponownie wyznaczyć $\sigma(t)$, widmo i granice $N_i^-$, sprawdzić mikroruch oraz ślady frettingu właściwą metodą i rozstrzygnąć mechanikę pękania. To jest ilościowy transfer obciążeń i historii egzemplarza, nie obietnica kompletnego programu trwałości.

Gdy dopuszcza się istnienie wady, właściwą wielkością jest współczynnik intensywności naprężeń

$$K_I=Y\sigma\sqrt{\pi a},$$

gdzie $a$ charakteryzuje wymiar pęknięcia, $Y$ zależy od geometrii i sposobu obciążenia, a $\sigma$ jest naprężeniem nominalnym. Kryterium kruchego pękania wymaga $K_I$ poniżej właściwej odporności na pękanie z odpowiednim marginesem i warunkami grubości, temperatury oraz środowiska. W obszarze stabilnego wzrostu popularna relacja Parisa

$$\frac{da}{dN}=C(\Delta K)^m$$

wiąże przyrost pęknięcia na cykl z zakresem $\Delta K$. Stałe $C,m$, próg propagacji i wpływ stosunku obciążeń muszą pochodzić z danych reprezentatywnych. Relacja nie opisuje poprawnie całego zakresu od inicjacji po niestabilne pękanie, dlatego integracja od wady wykrywalnej $a_0$ do krytycznej $a_c$ wymaga właściwego prawa i zachowawczych własności.

Dębski i Dębski rozróżniają podejścia bezpiecznej trwałości, uszkodzenia kontrolowanego i tolerancji uszkodzeń oraz podkreślają łańcuch od identyfikacji obciążeń przez analizę naprężeń do prób i monitorowania eksploatacji. Ich studia dotyczą konstrukcji lotniczych, więc są bezpośrednim źródłem metod zmęczeniowych, lecz tylko analogią dla specyficznych interfejsów OR-D1.¹⁰ Dla fikcyjnego obiektu przyjmujemy strategię mieszaną: elementy bez wiarygodnej inspekcji przed użyciem muszą wykazać bezpieczną trwałość z odpowiednim rozrzutem, natomiast wymienne elementy krytyczne i dostępne interfejsy mogą korzystać z programu tolerancji uszkodzeń, jeśli inspekcja ma udowodnione prawdopodobieństwo wykrycia.

Badania nieniszczące NDE (nondestructive evaluation) są częścią modelu ryzyka, nie rytuałem jakościowym. Stwierdzenie „brak wskazań” oznacza tylko brak wskazań przekraczających próg konkretnej metody, przy danym dostępie, orientacji wady, stanie powierzchni i kwalifikacji personelu. Krzywa POD (probability of detection) wiąże wymiar lub inną miarę wady z prawdopodobieństwem jej wykrycia. Wymiar początkowy do analizy propagacji nie może być arbitralnie mniejszy niż zdolność kwalifikowanej kontroli. NASA-STD-5009C ustanawia wymagania NDE dla systemów i komponentów NASA, lotnych lub naziemnych, gdy wymagana jest kontrola odporności na pękanie i ilościowa demonstracja POD. Opisuje demonstrację POD, procedury i kwalifikację personelu; tutaj wspiera strukturę procesu, nie wybór techniki ani dopuszczalnego rozmiaru wady dla OR-D1.¹²

NASA-STD-5019A ujmuje kontrolę odporności na pękanie (fracture control) jako program powiązany z planem, klasyfikacją części, analizą, próbami, NDE i identyfikowalnością. Zakres obowiązkowy standardu to zasadniczo załogowy sprzęt lotów kosmicznych (human-rated spaceflight hardware); programy niezałogowe mogą dobrowolnie narzucić wymagania. Zasada identyfikowalnego łańcucha od potencjalnej wady do ograniczenia eksploatacji jest tutaj tylko analogią procesową.¹¹ Podręcznik NASA-HDBK-5010, Volume 1, Revision A rozwija tę zasadę dla części metalowych, kompozytowych i klejonych oraz łączy widmo obciążeń, ocenę analityczną lub próbę, kontrolę wad i identyfikowalność. Jest wskazówką wdrożeniową dla sprzętu lotów kosmicznych, a nie gotowym planem OR-D1.⁷ Dla mocowania A-07 oznacza to co najmniej rozdzielenie trzech pytań: czy może powstać pęknięcie wskutek poślizgu i frettingu, czy dostępna metoda zobaczy je przed wymiarem krytycznym oraz czy pozostała liczba cykli przekracza udokumentowane przyszłe użycie z marginesem.

Jakość wykonania wchodzi do równań przez rozkład, a nie przez ogólne hasło. Mimośród otworu zmienia naprężenie lokalne; chropowatość i ślad narzędzia zmieniają inicjację; porowatość spoiny lub kleju zmienia efektywny przekrój; moment dokręcenia nie jest bezpośrednią miarą siły zacisku; grubość powłoki i owalizacja zmieniają wyboczenie; utwardzenie kleju zmienia sztywność i skurcz. Każda cecha krytyczna potrzebuje mierzalnej tolerancji, zdolności procesu oraz reakcji na przekroczenie. Statystyczna zgodność partii nie zastępuje oceny pojedynczej nieakceptowalnej wady w elemencie krytycznym.

Historia OR-D1 wiąże się teraz w jeden warunkowy ciąg hipotez do sprawdzenia. Jeżeli połączenie traci napięcie wstępne, rośnie jego mikroruch; mikroruch może powodować fretting i inicjację pęknięcia, a zmieniona podatność może obniżać częstotliwość oraz tłumienie postaci drgań. Jeżeli tej zmianie towarzyszy wzrost oporu kontaktowego, zmienia się również gradient cieplny. Studium wykazało wartość zacisku po cyklu niższą od nominalnej, lecz bez pomiaru przed cyklem nie dowiodło pierwszego ogniwa tego łańcucha. Dlatego cały ciąg pozostaje hipotezą wspólnej regresji, a nie opisem już zmierzonej historii. Cykl, który oddzielne zespoły nazwałyby „mechanicznym”, „dynamicznym” i „termicznym”, może być jednym zdarzeniem fizycznym; kontrola konfiguracji musi umożliwić sprawdzenie tej zależności.

Niepewność, weryfikacja i walidacja bez dopasowywania wyniku

Niepewność jest właściwością tezy i jej modelu, a nie jedną procentową liczbą dla całego pojazdu. Dla wielkości wyjściowej $y=f(x_1,\ldots,x_n)$ liniowa propagacja ma postać

$$u_y^2=\boldsymbol J\,\Sigma_x\,\boldsymbol J^T+u_{model}^2, \qquad J_i=\frac{\partial f}{\partial x_i},$$

gdzie $\Sigma_x$ jest macierzą kowariancji wejść, $\boldsymbol J$ wektorem wrażliwości, a $u_{model}$ składnikiem niedoskonałości postaci modelu. Zapis z pełną kowariancją jest ważny: $q_\infty$, $M$ i $Re$ korzystają z tych samych pomiarów atmosfery i prędkości, więc nie są niezależne. Podobnie sztywność i przewodność styku zależą od wspólnej siły docisku. Pierwiastkowe sumowanie wszystkich procentów jako niezależnych może zaniżyć lub zawyżyć wynik.

Losowanie Monte Carlo jest użyteczne dla nieliniowości, progów i rozkładów innych niż normalny, lecz liczba prób nie naprawia błędnej fizyki. Należy rozdzielić niepewność aleatoryczną, opisującą naturalny rozrzut, od epistemicznej, wynikającej z braku wiedzy. Pierwszą można opisać rozkładem na podstawie populacji; drugą często uczciwiej reprezentuje przedział modeli lub wariant hipotezy. W studium cieplnym przed pomiarem brak członu klejowego był niepewnością modelu, nie losowym odchyleniem czujnika.

Weryfikacja pyta, czy równania rozwiązano poprawnie. Obejmuje kontrolę jednostek, testy graniczne, zbieżność siatki i kroku, bilans energii i pędu, porównanie z rozwiązaniem analitycznym oraz niezależną implementację. Walidacja pyta, czy właściwe równania i parametry wystarczają do zamierzonego użycia. Obejmuje porównanie z pomiarem reprezentującym zjawisko, konfigurację i zakres. Korelacja jest kontrolowaną aktualizacją parametrów na podstawie danych; nie może naruszać fizycznych granic ani zużywać wszystkich danych, które miały sprawdzać zdolność predykcyjną. NASA-STD-7009B porządkuje wiarygodność modelu wokół zamierzonego użycia, założeń, weryfikacji, walidacji, niepewności i raportowania wyniku. Nie używamy terminu „akredytacja”, ponieważ samo wewnętrzne zatwierdzenie nie dodaje dowodu naukowego, a nie jest potrzebne do tej decyzji.¹³

Praktyczna hierarchia dowodu OR-D1 ma cztery poziomy. Najpierw kontrola rzędu wielkości wykrywa błąd jednostek lub znaku. Następnie model zredukowany pokazuje zależność przyczynową, jak w obliczeniach $N_a$, stałej czasowej i opóźnienia fazowego. Model szczegółowy rozwiązuje lokalne pola, geometrię i nieliniowość. Na końcu pomiar sprawdza sygnatury rozróżniające hipotezy. Zgodność na ostatnim poziomie jest przekonująca tylko wtedy, gdy trzy wcześniejsze pozostają spójne.

Miary korelacji muszą odpowiadać decyzji. Dla obciążeń porównuje się amplitudę, czas, fazę, rozkład i całkę impulsu; dla konstrukcji — reakcje, odkształcenia i przemieszczenia; dla dynamiki — częstotliwość, tłumienie, postać oraz zespolone FRF; dla termiki — tempo narastania, maksimum, gradient i chłodzenie. Jeden błąd RMS może ukryć przesunięcie fazowe albo lokalne maksimum. Kryteria ustala się przed obejrzeniem wyniku próby, aby uniknąć przesuwania bramki po fakcie.

NASA-STD-5002B wymaga udokumentowania konfiguracji, źródeł wymuszeń, matematycznego modelu, tłumienia, niepewności oraz walidacji i korelacji.¹ Dla OR-D1 przekładamy to na identyfikator pakietu danych: wersję geometrii, materiałów i połączeń; stan masowy; bazę aerodynamiczną z domeną; atmosferę; model MES i redukcję; model cieplny; sensor, aktuator i oprogramowanie; protokół próby; skrypt przetwarzania; budżet niepewności. Wynik bez takiej metryki nie jest powtarzalnym dowodem.

Aktualizacja po próbie musi zachować dane przed i po zmianie. Przykładowo skorygowana sztywność złącza poprawiła częstotliwość postaci drgań $42{,}3 \mathrm{Hz}$, ale powinna także przewidzieć wynik sztywnego imitatora $47{,}6 \mathrm{Hz}$ i nie pogorszyć innych postaci. Człon skurczu kleju powinien odtworzyć nie tylko maksimum $0{,}31 \mathrm{mrad}$, lecz także histerezę chłodzenia. Jeśli parametr dopasowany do jednego kanału niszczy zgodność innego, hipoteza jest niepełna.

Kontrola regresji działa po zależnościach. Zmiana połączenia wymusza ponowną ocenę nośności, poślizgu, zmęczenia, modalności, przewodzenia i osiowania; zmiana bazy aerodynamicznej — 6-DOF, obciążeń, aktuatora i aeroservoelastyczności; zmiana masy lub stacji środka masy — ruchu, momentów, reakcji i postaci drgań; zmiana filtru — zapasów, opóźnienia i odpowiedzi HIL. Nie oznacza to zawsze powtórzenia każdej próby. Oznacza udokumentowaną analizę wpływu, która rozstrzyga, co pozostaje ważne i dlaczego.

Od kwalifikacji do nadzoru nad egzemplarzem

Kwalifikacja (qualification) odpowiada na pytanie, czy projekt i proces mają wykazany zapas w określonej obwiedni. Odbiór (acceptance) odpowiada, czy konkretny egzemplarz został poprawnie wykonany i nie ma ujawnionych wad. Próba kwalifikacyjna może korzystać z poziomu lub czasu ostrzejszego od eksploatacyjnego i bywa zużywająca; próba odbiorcza powinna wykrywać wady wykonania bez niedopuszczalnego uszczuplenia trwałości. NASA-STD-7001C rozdziela cele kwalifikacji, podejścia protoflight, odbioru i kontroli jakości wykonania, wymaga również zdefiniowanej konfiguracji, tolerancji oraz danych. Standard dotyczy sprzętu statków kosmicznych, więc jego liczb i widm nie wolno przenosić do OR-D1; przydatna jest struktura rozumowania.¹⁵

NASA-STD-5001B podobnie rozróżnia kwalifikację projektu przez próbę prototypową od podejścia protoflight, w którym artykuł lotny jest badany na poziomie kwalifikacyjnym przy zmodyfikowanym czasie. To rozróżnienie ma konsekwencję: wynik egzemplarza kwalifikacyjnego nie obejmuje dowolnej późniejszej zmiany procesu, a zużyty artykuł nie staje się automatycznie egzemplarzem użytkowym.¹⁴ Dla OR-D1 wybór strategii musi uwzględnić możliwość inspekcji, rozrzut połączeń i kleju, zużycie zmęczeniowe oraz reprezentatywność masy i interfejsów.

Program kwalifikacyjny wywodzi przypadki z dowodu, nie z listy tradycyjnych badań. Próba obciążenia statycznego sprawdza drogę sił, nośność i stateczność. GVT sprawdza dynamikę modalną. Próba cieplna sprawdza pole temperatury, interfejs i osiowanie. HIL sprawdza zintegrowane opóźnienia oraz sprzężenie sensor–aktuator. Próby zmęczeniowe i NDE sprawdzają historię uszkodzenia. Każda ma zapisane wejście, mierzone sygnatury, niepewność, kryteria przerwania, kryteria zaliczenia i stan egzemplarza po badaniu. Próba bez wskazanej hipotezy może wytworzyć dużo danych i mało wiedzy.

Odbiór późniejszych egzemplarzy kontroluje cechy, które kwalifikacja uznała za krytyczne: masę i stację środka masy, napięcie wstępne oraz retencję połączenia, parametry procesu kleju, wybrane częstotliwości rezonansowe albo odpowiedź impulsową, szczelność, osiowanie, ciągłość torów pomiarowych i wskazania NDE. Granice odbioru wynikają z rozkładu kwalifikacyjnego oraz marginesu, a nie z wygody stanowiska. Odbiór nie powtarza próby niszczącej; korzysta z cech zastępczych tylko wtedy, gdy ich korelacja z mechanizmem uszkodzenia została wykazana.

Nadzór w okresie przechowywania (surveillance) sprawdza, czy dowód pozostaje aktualny. Starzenie kleju, relaksacja zacisku, korozja, nasiąkanie, degradacja uszczelnienia i dryf sensora zależą od czasu oraz środowiska. Plan próbkowania powinien być powiązany z modelem starzenia i krytycznością, a nie z arbitralnym kalendarzem. Wynik egzemplarza pobranego z partii aktualizuje ocenę partii w zakresie uzasadnionym wspólną historią; nie wolno ekstrapolować go na inną partię materiału albo odmienny magazyn bez analizy.

Po przekroczeniu, upadku, nieplanowanym cyklu termicznym albo zmianie oprogramowania status egzemplarza przechodzi w kwarantannę techniczną. Ocena obejmuje rekonstrukcję zdarzenia, porównanie z obwiednią, inspekcję mechanizmów wrażliwych oraz analizę wpływu. Sam brak widocznego uszkodzenia nie przywraca zgodności. Tak samo pozytywny test funkcjonalny nie wyklucza pęknięcia, utraty napięcia albo zmniejszenia zapasu fazy.

Bezpieczna granica badań

Opisane rachunki i badania są publiczną ilustracją metody. Nie stanowią instrukcji budowy, nastawiania sterowania ani prowadzenia próby rzeczywistego pocisku. Każdą próbę wykonuje wyłącznie uprawniony personel według zatwierdzonej analizy zagrożeń, w konfiguracji jednoznacznie zidentyfikowanej przed podaniem energii.

Artykuły obojętne i energetyczne muszą pozostawać organizacyjnie i fizycznie rozdzielone. Stan bezpieczny i stan uzbrojony podlegają formalnej kontroli konfiguracji, niezależnej weryfikacji oraz zasadom blokad; sam wpis w oprogramowaniu lub etykieta nie wystarcza. Próby stwarzające ryzyko prowadzi się zdalnie, w zabezpieczonym obszarze, z kontrolą dostępu, łącznością awaryjną i odpowiedzialnością służby bezpieczeństwa poligonu. Kryteria przerwania ustala się przed próbą i wiąże z mierzalnymi wielkościami — utratą telemetrii bezpieczeństwa, przekroczeniem obciążenia, temperatury, drgań, przemieszczenia albo granicy stanowiska — lecz konkretne progi należą do zatwierdzonej analizy danego obiektu, a nie do tej monografii. Po przerwaniu utrzymuje się stan bezpieczny do formalnej oceny; nie podejmuje się improwizowanej kontynuacji.

Jedna reguła decyzji

Dowód płatowca byłby zamknięty dopiero wtedy, gdy dla jednej kontrolowanej konfiguracji wspólny łańcuch atmosfera–aerodynamika–6-DOF–obciążenia–struktura–dynamika–termika–trwałość przewidywałby mierzalne sygnatury, próby rozróżniałyby konkurencyjne hipotezy, niepewność pozostawiałaby wymagany margines, a kwalifikacja, odbiór i nadzór chroniłyby te same cechy krytyczne. OR-D1 nie spełnia tego warunku: wykonano ilościowo tylko wskazane przekazania lokalne, a nie pełne obliczenie trajektorii i kwalifikacji. Dodatni margines w jednej domenie nie kompensuje braku dowodu w innej.

W OR-D1 ta reguła prowadzi do jednej decyzji końcowej: konfiguracja A-07 w opisanym stanie nie jest przyjęta. Powłoka ma dodatnie statyczne marginesy, lecz siła zacisku zmierzona po cyklu jest niższa od wartości nominalnej, retencji zacisku nie wykazano, a pokryciowy margines poślizgu jest ujemny. Skorelowany model wskazuje podatność złącza jako najlepiej wspartą, choć niewyłączną przyczynę obniżenia częstotliwości postaci drgań i zapasu wzmocnienia. Składnik odpowiedzi kleju sam nie przekracza limitu $0{,}25 \mathrm{mrad}$, lecz po dodaniu składnika termomechanicznego i dryfu sensora powoduje, że całkowity błąd osi osiąga $0{,}3175 \mathrm{mrad}$ i przekracza ten limit. Powrót do decyzji pozytywnej wymaga kontrolowanej zmiany interfejsu i wspólnej regresji trzech studiów, nie trzech niezależnych odstępstw.

Krótka kontrola zamknięcia ma zatem sens tylko jako indeks wykonanego dowodu:

• ta sama wersja geometrii, masy, stacji środka masy, bezwładności i interfejsów występuje w modelach i próbach;
• konwencje osi, znaków, punktów oraz pól odniesienia są jawne i sprawdzone testami granicznymi;
• każda wielkość graniczna ma fizycznie osiągalny przypadek, drogę obciążenia, budżet niepewności i właściwy mechanizm uszkodzenia;
• korelacja poprawia kilka niezależnych sygnatur, a nie tylko liczbę użytą do strojenia;
• każda zmiana ma analizę wpływu i regresję obejmującą wszystkie zależne domeny;
• egzemplarz ma śledzoną historię, pozytywny odbiór i aktualny status nadzoru;
• bezpieczeństwo próby oraz kryteria przerwania są zatwierdzone poza niniejszym tekstem.

Jeżeli choć jeden z tych warunków pozostaje nieudowodniony, wynik brzmi „brak podstaw do przyjęcia”, a nie „prawdopodobnie wystarczy”. Taka reguła jest bardziej wymagająca niż zbiór efektownych symulacji, ale dokładnie po to buduje się zintegrowany dowód: aby decyzja nie zależała od tego, która specjalność ostatnia oglądała własny wykres.