Przejdź do treści

Praca wirtualna sił zewnętrznych i składanie równań globalnych

W Dyskretyzacji pracy wirtualnej sił wewnętrznych lewa strona postaci słabej została sprowadzona do elementowego wektora sił wewnętrznych \(\boldsymbol{q}^e\) (metoda UL) lub \(\boldsymbol{Q}^e\) (metoda TL). W tym rozdziale wprowadzono elementowy węzłowy wektor sił zewnętrznych \(\boldsymbol{F}^e\) na podstawie pracy wirtualnej sił zewnętrznych, a następnie zastosowano operację składania, która porządkuje i sumuje wielkości węzłowe elementów według globalnych numerów węzłów, ostatecznie prowadząc do nieliniowego układu równań względem przemieszczeń węzłowych rozwiązywanego w nieliniowej analizie strukturalnej FrontISTR.

Dekompozycja elementowa pracy wirtualnej sił zewnętrznych

Prawą stronę zasady pracy wirtualnej można rozłożyć element po elemencie na pracę wirtualną sił zewnętrznych, składającą się z sił objętościowych i zadanych obciążeń powierzchniowych na mechanicznych brzegach. Aby zapisać interpolację przemieszczenia wprowadzoną w Funkcjach kształtu i aproksymacji metodą elementów skończonych w postaci macierzowej, dla węzła \(\alpha\) i jego funkcji kształtu \(N_\alpha^e\) używa się bloku \(d \times d\) \(\boldsymbol{N}_\alpha\), w którym funkcja ta znajduje się na przekątnej, oraz macierzy \(\boldsymbol{N} = [\boldsymbol{N}_1, \ldots, \boldsymbol{N}_{n_e}]\) utworzonej przez poziome zestawienie tych bloków, tak aby \(\delta\boldsymbol{u} = \boldsymbol{N}\, \delta\boldsymbol{u}^e\). Po podstawieniu tego wyrażenia do pracy wirtualnej sił zewnętrznych zapisanej w konfiguracji odniesienia otrzymuje się

\[ \delta W^{\mathrm{ext}} = \sum_e \delta\boldsymbol{u}^{eT} \boldsymbol{F}^e, \qquad \boldsymbol{F}^e_\alpha = \int_{\Omega^e_0} \boldsymbol{N}_\alpha^T \rho_0 \boldsymbol{g}\, dV + \int_{\Gamma^e_{0t}} \boldsymbol{N}_\alpha^T \bar{\boldsymbol{t}}_0\, d\Gamma_0 \]

gdzie elementowy węzłowy wektor sił zewnętrznych jest uporządkowany jako \(\boldsymbol{F}^e = (\boldsymbol{F}^{eT}_1, \ldots, \boldsymbol{F}^{eT}_{n_e})^T\). W ten sposób praca wirtualna sił zewnętrznych zostaje sprowadzona do tej samej postaci „elementowy wektor węzłowy × funkcja testowa” co strona sił wewnętrznych (przy zapisie w konfiguracji bieżącej tę samą postać otrzymuje się przez podstawienia \(dV \to dv\), \(\rho_0 \to \rho\) oraz \(\bar{\boldsymbol{t}}_0 \to \bar{\boldsymbol{t}}\)).

Składanie wielkości węzłowych elementów

Wielkości węzłowe \(\boldsymbol{Q}^e_\alpha, \boldsymbol{F}^e_\alpha\) wyznaczone dla każdego elementu są sumowane do wektorów globalnych uporządkowanych według globalnych numerów węzłów. Niech globalny numer węzła odpowiadający lokalnemu numerowi węzła \(\Omega^e\) elementu \(\alpha\) będzie równy

\[ \mathrm{gdx}(e, \alpha) = i_g \]

Wtedy elementowa wielkość węzłowa jest identyczna z odpowiadającą jej składową globalnej wielkości węzłowej (na przykład \(\boldsymbol{u}^e_\alpha = \boldsymbol{u}_{i_g}\)). Ponieważ węzeł \(i_g\) jest na ogół współdzielony przez wiele elementów, zdefiniuj zbiór par \(i_g\), których globalny numer węzła wynosi \((e, \alpha)\), jako

\[ \mathcal{E}(i_g) = \{ (e, \alpha) \mid \mathrm{gdx}(e, \alpha) = i_g \} \]

Użycie tego zbioru do przepisania sumy w postaci \(\sum_e \sum_\alpha = \sum_{i_g} \sum_{(e,\alpha) \in \mathcal{E}(i_g)}\) daje węzłową siłę wewnętrzną oraz globalny wektor sił wewnętrznych dla wszystkich \(n_g\) węzłów:

\[ \boldsymbol{Q}_{i_g} = \sum_{(e,\alpha) \in \mathcal{E}(i_g)} \boldsymbol{Q}^e_\alpha, \qquad \boldsymbol{Q} = (\boldsymbol{Q}^T_1, \ldots, \boldsymbol{Q}^T_{n_g})^T \]

Tutaj \(\boldsymbol{Q}_{i_g}\) jest wypadkową elementowych węzłowych sił wewnętrznych działających w węźle \(i_g\) i jest równe \(\boldsymbol{0}\), gdy nie działa żadna siła zewnętrzna i spełniona jest równowaga. Ta sama procedura w metodzie UL daje \(\boldsymbol{q}_{i_g}, \boldsymbol{q}\); ponieważ ich wartości liczbowe spełniają \(\boldsymbol{q} = \boldsymbol{Q}\), poniżej używa się oznaczenia \(\boldsymbol{Q}\), z wyjątkiem miejsc, w których konieczne jest rozróżnienie. Globalny wektor sił zewnętrznych \(\boldsymbol{F}\) otrzymuje się przez takie samo sumowanie.

W implementacji zbiór \(\mathcal{E}(i_g)\) nie jest tworzony jawnie; zamiast tego wkłady są dodawane do odpowiednich składowych wewnątrz pętli po elementach.

Zainicjuj globalny wektor sił wewnętrznych Q wartością 0: Q_{i_g} = 0  (i_g = 1, ..., n_g)
for e = 1 to (liczba elementów)
    for α = 1 to n_e
        i_g = gdx(e, α)
        Q_{i_g} += Q^e_α
    end for
end for

Globalny wektor sił zewnętrznych \(\boldsymbol{F}\) jest konstruowany tą samą procedurą. Operację dodawania i zapisywania elementowych wielkości węzłowych do wektorów i macierzy numerowanych według globalnych numerów węzłów nazywa się składaniem. Dla tensorów drugiego rzędu związanych z dwoma numerami węzłów (takich jak macierze sztywności) analogiczne składanie otrzymuje się przy użyciu zbioru \(\mathcal{E}^2(i_g, i_h) = \{ (e, \alpha, \beta) \mid \mathrm{gdx}(e, \alpha) = i_g\ \mathrm{and}\ \mathrm{gdx}(e, \beta) = i_h \}\) (konkretną konstrukcję opisano w Stycznej macierzy sztywności).

Nieliniowe równania do rozwiązania

Po podstawieniu złożonych sił wewnętrznych i zewnętrznych do zasady pracy wirtualnej oraz wykorzystaniu faktu, że obowiązuje ona dla dowolnej funkcji testowej \(\delta\boldsymbol{u}^n\) spełniającej geometryczne warunki brzegowe, otrzymuje się

\[ \boldsymbol{Q}(\boldsymbol{u}^n) - \boldsymbol{F}(\boldsymbol{u}^n) = \boldsymbol{0} \]

W kontekście analizy przyrostowej (Schemat analizy przyrostowej) przywróć indeks czasu \(_{n+1}\) i pomiń indeks górny \(^n\) oznaczający globalny wektor węzłowy. Równanie do rozwiązania ma wtedy postać

\[ \boldsymbol{Q}(\boldsymbol{u}_{n+1}) - \boldsymbol{F}(\boldsymbol{u}_{n+1}) = \boldsymbol{0} \]

W ten sposób zdyskretyzowane zagadnienie brzegowe polegające na wyznaczeniu przemieszczenia węzłowego \(t_{n+1}\) w czasie \(\boldsymbol{u}_{n+1}\) sprowadza się do rozwiązania tego nieliniowego równania względem przemieszczenia wraz z geometrycznymi warunkami brzegowymi. Linearyzację równania i konstrukcję stycznej macierzy sztywności opisano w Stycznej macierzy sztywności, a iteracyjną metodę rozwiązania w Metodzie Newtona-Raphsona.

Powiązane tematy