Przejdź do treści

Translation update pending

The Polish translation may be out of date because the Japanese source has changed. The existing translation remains available while it is being updated.

Wysokowydajne sformułowania elementów

Standardowa formulacja przemieszczeniowa przedstawiona w Formulacji trójwymiarowych elementów bryłowych wykazuje nadmierną sztuczną sztywność, zwaną blokowaniem (blokowaniem objętościowym lub blokowaniem ścinania), gdy jest stosowana do materiałów prawie nieściśliwych lub cienkościennych konstrukcji zdominowanych przez zginanie. Aby temu zapobiec, FrontISTR udostępnia metody B-bar i F-bar, które zastępują jedynie część objętościową macierzy B lub gradientu deformacji; elementy niezgodne z wewnętrznymi stopniami swobody; mieszane elementy u-p, w których ciśnienie jest traktowane jako niezależne pole niewiadome; oraz elementy powłokowe MITC i elementy belkowe przeznaczone dla konstrukcji płytowych i belkowych.

W tym rozdziale podsumowano formulacje tych zaawansowanych elementów i elementów konstrukcyjnych na poziomie poszczególnych elementów.

Metoda B-bar

Gdy 8-węzłowy liniowy element sześciościenny jest używany dla materiału prawie nieściśliwego, pole odkształceń wewnątrz elementu może pozostawać w konflikcie z warunkiem stałej objętości, powodując nadmierną sztywność zwaną blokowaniem objętościowym. Metoda B-bar łagodzi to nadmierne ograniczenie przez zastąpienie składowych macierzy B odpowiedzialnych za rozszerzenie objętościowe wartościami wyznaczonymi w środku elementu [Hughes1980].

W środku elementu \(\boldsymbol{r} = \boldsymbol{0}\) niech macierz B obliczona z pochodnych przestrzennych funkcji kształtu będzie oznaczona przez \(\bar{\boldsymbol{B}}\); w punkcie całkowania \(\boldsymbol{r}\) niech zwykła macierz B będzie oznaczona przez \(\boldsymbol{B}(\boldsymbol{r})\). W zależności przemieszczenie–odkształcenie dla węzła \(\alpha\) i stopnia swobody \(i\) składowe odkształcenia objętościowego \((\varepsilon_{11}, \varepsilon_{22}, \varepsilon_{33})\) są modyfikowane przez dodanie

\[ \Delta B_{i\alpha} = \tfrac{1}{3}\bigl(\bar{B}_{i\alpha}(\boldsymbol{0}) - B_{i\alpha}(\boldsymbol{r})\bigr) \]

Dla składowych ścinania \((\varepsilon_{12}, \varepsilon_{23}, \varepsilon_{31})\) używana jest zwykła \(\boldsymbol{B}\). Macierz sztywności elementu i wektor sił wewnętrznych są składane przy użyciu otrzymanej macierzy B-bar.

FrontISTR udostępnia tę formulację specjalnie dla 8-węzłowego liniowego elementu sześciościennego (ID elementu 361; zob. System numeracji elementów) i można ją stosować w formulacjach małych deformacji, Total Lagrange oraz Updated Lagrange.

Metoda F-bar

Przy dużych deformacjach zmiana objętości wchodzi nieliniowo przez gradient deformacji \(\boldsymbol{F}\). Metoda F-bar [deSouzaNeto1996] stosuje na poziomie gradientu deformacji zabieg przeciw blokowaniu objętościowemu analogiczny do metody B-bar.

Niech stosunek objętości gradientu deformacji wyznaczony w środku elementu \(\boldsymbol{r} = \boldsymbol{0}\) będzie równy \(J_0 = \det \boldsymbol{F}(\boldsymbol{0})\), a stosunek objętości w punkcie całkowania niech będzie równy \(J = \det \boldsymbol{F}(\boldsymbol{r})\). Gradient deformacji w punkcie całkowania zastępuje się przez

\[ \bar{\boldsymbol{F}} = \left(\frac{J_0}{J}\right)^{1/3} \boldsymbol{F} \]

Daje to \(\det \bar{\boldsymbol{F}} = J_0\), dzięki czemu stosunek objętości w całym elemencie zostaje dopasowany do wartości w środku elementu. Zastąpiona \(\bar{\boldsymbol{F}}\) jest używana do wyznaczania naprężeń i konstrukcji macierzy odkształcenie–przemieszczenie, a styczna sztywność elementu obejmuje dodatkowe składniki związane z tym zastąpieniem.

FrontISTR implementuje metodę F-bar specjalnie dla 8-węzłowego liniowego elementu sześciościennego i może ona być stosowana zarówno dla małych deformacji, jak i nieliniowej geometrii Total Lagrange / Updated Lagrange.

Elementy niezgodne

8-węzłowy liniowy element sześciościenny nie ma modów odkształcenia wymaganych do zginania i wykazuje blokowanie zginania w problemach zdominowanych przez zginanie. Elementy niezgodne [Taylor1976] wprowadzają dodatkowe mody przemieszczeń wewnątrz elementu w celu skompensowania tego braku.

Oprócz przemieszczeń węzłowych elementu \(\boldsymbol{u}^e\) wprowadza się stopnie swobody modów niezgodnych \(\boldsymbol{\alpha} \in \mathbb{R}^{9}\), istniejące wyłącznie wewnątrz elementu, jako 3 kierunki × 3 mody na element, a pole przemieszczeń aproksymuje się jako

\[ \boldsymbol{u}(\boldsymbol{r}) = \sum_{\alpha=1}^{8} N_\alpha^e(\boldsymbol{r})\, \boldsymbol{u}^e_\alpha + \sum_{k=1}^{3} M_k(\boldsymbol{r})\, \boldsymbol{\alpha}_k \]

Dla współrzędnych naturalnych \(\boldsymbol{r} = (\xi, \eta, \zeta)\) niezgodne funkcje kształtu przyjmuje się jako \(M_1 = 1 - \xi^2\), \(M_2 = 1 - \eta^2\) oraz \(M_3 = 1 - \zeta^2\). Nie gwarantują one ciągłości na granicach elementów, ale dodają przestrzeń wewnętrzną zdolną odtwarzać odkształcenia związane z modami zginania.

Sztywność elementu jest najpierw składana w postaci blokowej zewnętrzne DOF/wewnętrzne DOF

\[ \begin{bmatrix} \boldsymbol{K}_{dd} & \boldsymbol{K}_{d\alpha} \\ \boldsymbol{K}_{\alpha d} & \boldsymbol{K}_{\alpha\alpha} \end{bmatrix} \begin{bmatrix} d\boldsymbol{u}^e \\ d\boldsymbol{\alpha} \end{bmatrix} = \begin{bmatrix} \boldsymbol{F}^e_{\text{ext}} \\ \boldsymbol{0} \end{bmatrix} \]

Po złożeniu w tej postaci kondensacja statyczna eliminuje wewnętrzne stopnie swobody za pomocą \(d\boldsymbol{\alpha} = -\boldsymbol{K}_{\alpha\alpha}^{-1}\boldsymbol{K}_{\alpha d}\,d\boldsymbol{u}^e\), dając sztywność elementu obejmującą wyłącznie zewnętrzne stopnie swobody

\[ \boldsymbol{K}^e = \boldsymbol{K}_{dd} - \boldsymbol{K}_{d\alpha}\,\boldsymbol{K}_{\alpha\alpha}^{-1}\,\boldsymbol{K}_{\alpha d} \]

która jest przekazywana do składania globalnego.

FrontISTR implementuje element niezgodny specjalnie dla 8-węzłowego liniowego elementu sześciościennego (C3D8IC) i można go stosować w formulacjach małych deformacji, Total Lagrange oraz Updated Lagrange.

Mieszane elementy U-P

Podczas gdy metody B-bar i F-bar korygują składową objętościową w ramach formulacji przemieszczeniowej, mieszany element u-p (U-P) wykorzystuje formulację mieszaną [Bathe1996], w której ciśnienie \(\lambda\) jest wprowadzane jako pole niewiadome niezależne od przemieszczenia. Dla materiałów prawie nieściśliwych, takich jak materiały gumopodobne ze współczynnikiem Poissona bardzo bliskim 0.5 lub metale po odkształceniu plastycznym, wymuszanie warunku stałej objętości wyłącznie za pomocą pola przemieszczeń powoduje blokowanie objętościowe; traktowanie ciśnienia jako niezależnej zmiennej łagodzi to ograniczenie.

Naprężenie rozdziela się na składową dewiatorową i ciśnieniową jako

\[ \boldsymbol{\sigma} = \boldsymbol{\sigma}_{\mathrm{dev}} + \lambda\,\boldsymbol{I}, \qquad \boldsymbol{\sigma}_{\mathrm{dev}} = \mathbf{D}_{\mathrm{dev}}\,\boldsymbol{\varepsilon} \]

Tutaj \(\mathbf{D}_{\mathrm{dev}}\) jest dewiatorową macierzą sprężystości uzyskaną przez usunięcie z macierzy sprężystości części objętościowej proporcjonalnej do modułu objętościowego \(K\). Ciśnienie \(\lambda\) i odkształcenie objętościowe \(g = \mathrm{tr}\,\boldsymbol{\varepsilon}\) są powiązane warunkiem przez ściśliwość \(\alpha^{-1} = 1/K\)

\[ g - \alpha^{-1}\lambda = 0 \]

Dyskretyzacja z przemieszczeniem \(\boldsymbol{u}\) i ciśnieniem \(\lambda\) jako niewiadomymi prowadzi do sprzężonego układu elementu

\[ \begin{bmatrix} \mathbf{K}_{uu} & \mathbf{K}_{up} \\ \mathbf{K}_{up}^{T} & \mathbf{K}_{pp} \end{bmatrix} \begin{bmatrix} d\boldsymbol{u} \\ d\lambda \end{bmatrix} = \begin{bmatrix} \boldsymbol{f}_{u} \\ \boldsymbol{f}_{p} \end{bmatrix} \]

gdzie \(\mathbf{K}_{uu}\) zawiera dewiatorowy wkład sprężysty oraz, dla dużych deformacji, sztywność geometryczną; \(\mathbf{K}_{up}\) sprzęga odkształcenie objętościowe i ciśnienie; a \(\mathbf{K}_{pp} = -\int \alpha^{-1}\,\boldsymbol{N}_p \boldsymbol{N}_p^{T}\,dV\) jest składnikiem stabilizacji ciśnienia (\(\boldsymbol{N}_p\) jest funkcją kształtu ciśnienia). Ponieważ stopnie swobody ciśnienia są wewnętrzne dla elementu,

\[ \mathbf{K}_{\mathrm{eff}} = \mathbf{K}_{uu} - \mathbf{K}_{up}\,\mathbf{K}_{pp}^{-1}\,\mathbf{K}_{up}^{T} \]

w tej postaci wykonuje się kondensację statyczną, a efektywna sztywność obejmująca wyłącznie zewnętrzne (przemieszczeniowe) stopnie swobody jest przekazywana do składania globalnego.

FrontISTR implementuje element U-P specjalnie dla 8-węzłowego liniowego elementu sześciościennego, z jednym stopniem swobody ciśnienia na element (stałym w obrębie elementu). Można go stosować w formulacjach małych deformacji, Total Lagrange oraz Updated Lagrange. W Updated Lagrange naprężenie dewiatorowe jest aktualizowane z użyciem obiektywnej prędkości naprężenia (typu Jaumann/Hughes-Winget), a następnie ciśnienie \(\lambda\,\boldsymbol{I}\) jest narzucane przy użyciu wartości uzyskanej przez kondensację statyczną.

Elementy powłokowe

Elementy powłokowe oparte na teorii płyt/powłok Reissnera-Mindlina są stosowane do cienkościennych konstrukcji płytowych i powłokowych. Wraz ze zmniejszaniem grubości niskiego rzędu przemieszczeniowe elementy płytowe/powłokowe zawyżają poprzeczne odkształcenie ścinające przez mechanizm analogiczny do blokowania objętościowego (blokowanie ścinania), powodując nadmierny wzrost sztywności dla modów zginania. Metoda MITC (Mixed Interpolation of Tensorial Components) [Dvorkin1984] [Bathe1986] unika tego problemu przez ponowne próbkowanie wyłącznie składowych odkształcenia ścinającego w zdefiniowanych punktach wiążących elementu i interpolację tych wartości z powrotem w całym elemencie.

Węzły elementu powłokowego MITC leżą na powierzchni środkowej, a każdy węzeł ma sześć stopni swobody: trzy składowe translacyjne i trzy składowe obrotowe względem układu współrzędnych normalnego do powierzchni środkowej. Sztywność elementu wyznacza się przez trójwymiarowe całkowanie Gaussa po naturalnych współrzędnych powierzchni środkowej i kierunku grubości, a grubość \(h\) jest podawana jako właściwość elementu podczas wyznaczania odpowiedzi konstytutywnej.

FrontISTR udostępnia MITC3 (ID elementu 731), MITC4 (741) i MITC9 (743), reprezentowane przez pojedynczą warstwę powierzchni środkowej, a także warstwowe elementy powłokowe MITC3-shell361 (761, 3\(\times\)2 węzłów, 3 DOF na węzeł) i MITC4-shell361 (781, 4\(\times\)2 węzłów, 3 DOF na węzeł), których węzły są rozmieszczone w dwóch warstwach przez grubość. W warstwowych elementach powłokowych stopnie swobody węzłów obejmują tylko trzy składowe translacyjne, natomiast mody zginania odpowiadające obrotowym stopniom swobody są reprezentowane przez układ dwóch warstw.

Elementy belkowe

Elementy liniowe, takie jak dźwigary i konstrukcje ramowe, są dyskretyzowane elementami belkowymi. Elementy belkowe opisują wydłużenie, zginanie i skręcanie wyłącznie za pomocą przemieszczeń i obrotów węzłów leżących na osi belki, natomiast deformację wewnątrz przekroju reprezentują stałe przekrojowe.

Każdy węzeł belki ma sześć stopni swobody: trzy składowe translacyjne i trzy składowe obrotowe. Sztywność elementu otrzymuje się, tworząc w lokalnym układzie współrzędnych, wyznaczonym przez oś elementu i oś odniesienia zadaną w danych przekroju, macierze sztywności w postaci zamkniętej dla wydłużenia, zginania w dwóch kierunkach i skręcania, a następnie transformując je do globalnego układu współrzędnych. Pole przekroju \(A\), drugie momenty pola \(I_y\) i \(I_z\) względem lokalnych osi \(y\) i \(z\) oraz stała skręcania są zadawane jako dane przekroju; wraz z modułem Younga materiału \(E\) i modułem ścinania \(G\) wyznaczanym z modułu Younga i współczynnika Poissona określają odpowiednie współczynniki sztywności.

Dla zginania dostępne są dwie formulacje: belka Bernoulliego-Eulera, w której zakłada się, że przekrój po deformacji pozostaje prostopadły do osi belki i pomija się deformację ścinającą, oraz belka Timoshenki, która uwzględnia deformację ścinającą. Dla dwuwęzłowego elementu belkowego 611 można wybrać obie formulacje, przy czym domyślna jest belka Bernoulliego-Eulera. Dla belki Timoshenki, przy długości elementu \(L\), współczynniki opisujące udział deformacji ścinającej

\[ \Phi_y = \frac{12 E I_z}{G A_s L^2}, \qquad \Phi_z = \frac{12 E I_y}{G A_s L^2} \]

są uwzględniane we współczynnikach sztywności dla odpowiednich płaszczyzn zginania. Na przykład współczynnik sztywności dla ugięcia w lokalnym kierunku \(y\) wynosi \(12 E I_z / \{ L^3 (1 + \Phi_y) \}\). Wielkość \(A_s\) jest polem przekroju na ścinanie. Ponieważ dane przekroju w FrontISTR nie zawierają osobnej wartości \(A_s\), przyjmuje się \(A_s = A\). Odpowiada to przyjęciu współczynnika korekcyjnego ścinania równego 1. Dla \(\Phi_y = \Phi_z = 0\) otrzymuje się sztywność belki Bernoulliego-Eulera.

Element 641 jest elementem o mieszanych stopniach swobody, reprezentującym belkę za pomocą dwóch węzłów z translacyjnymi stopniami swobody oraz dwóch węzłów z obrotowymi stopniami swobody, i jest sformułowany jako belka Bernoulliego-Eulera.

Tematy powiązane