Облако своими руками для расчета пространственных стержней методом конечных элементов на Node js, React js и Three js / Habr

Моделирование пространственных стержневых систем

ammonit3d в каталоге RuStore

www.rustore.ru

В данной статье (а возможно цикле статей) речь пойдет о собственной разработке облачного SPA приложения по моделированию пространственных стержневых систем методом конечных элементов с численно-аналитическим решением для инженеров-проектировщиков в основе которого математическая модель Эйлера-Бернулли, вариационные принципы и итерационный метод сопряжённых градиентов применяемый для большеразмерных СЛАУ с разреженной матрицей жёсткости с одной стороны, и JavaScripts экосистема облака, выполненного в стеке Node js, Express js бэкенд части, и React js, MobX, Three js, glsl shaders фронтенд части с другой стороны. Отображение эпюр усилий в пространственных стержневых элементах реализовано на шейдерах vertexShader и fragmentShader. Это позволяет вычислять эпюры для каждого стержня на лету и выполнять отображение графиков (в общем случае полиномов 5 степени) в пространстве мгновенно.

Клиентская часть представляет собой одностраничное SPA приложение с возможностью регистрации новых пользователей, профиля пользователя для управления приватными моделями, списком публичных моделей и основного интерфейса приложения, который состоит из навигационной панели и окна 3Д сцены. Веб приложение построено на адаптивной вёрстке. Это позволяет организовать удобное отображение элементов интерфейса на мониторах, планшетах и смартфонах.

Математическая модель Эйлера-Бернулли

Тех кого мало интересует математическая или механическая основа модели, обыкновенные дифференциальные уравнения и различные формульные выкладки советую пропустить этот раздел. Как было упомянуто выше, в основе численно-аналитического решения пространственных стержней (балок, ферм, рам, пространственных связей) методом конечных элементов лежит механическая модель изогнутой балки или стержня длиной L с заданной изгибной жёсткостью EJ на которую действуют сосредоточенная сила F или момент M, а также распределённая сила q(x) или момент m(x) по длине стержня, продольная ось которого x1, вертикальная x2 проходит через начало стержня и x3 направлена на нас из точки пересечения x1, x2. Полагая сечения плоскими до и после изгиба, а также одноосное деформированное состояние, то уравнение изогнутой оси стержня записывается в виде обыкновенного неоднородного дифференциального уравнения 4-го порядка:

$EJ_3 \dfrac{d^4v}{dx^4} = q(x) - \dfrac{dm}{dx}, \quad q(x) = kx + q_0, \quad m(x) = rx + m_0,\quad x\in[0,L],$

общее решение которого состоит из частного решения правой части (первые два слагаемых) и решения однородного дифференциального уравнения:

$v(x)=\dfrac{kx^5}{120EJ_3}+\dfrac{(q_0-r)x^4}{24EJ_3}+a_4x^3+a_3x^2+a_2x+a_1,$

Если записать граничные условия

$\begin{equation}\begin{cases}\begin{array} \, v(0)=v^1,\quad v,_1(0)=\varphi^1,\\v(L)=v^2, \quad v,_1(L)= \varphi^2,\end{array}\end{cases}\end{equation}$

то неизвестные коэффициенты можно выразить через свёртку (сумму произведений) прогибов вдоль x2 и поворотов вокруг x3 на концах изгибаемой балки с эрмитовыми полиномами N:

$\begin{array} \, v(x) = \dfrac{k}{EJ_3}(\dfrac{x^5}{120}-\dfrac{x^3L^2}{40}+\dfrac{x^2L^2}{60})+\dfrac{q_0-r}{EJ_3}(\dfrac{x^4}{24}-\dfrac{x^3L}{12}+\dfrac{x^2L^2}{24})+N_p^i\alpha_i^p,\\N_p^1\alpha_1^p=(\dfrac{2x^3}{L^3}-\dfrac{3x^2}{L^2}+1)\,v^1 + (-\dfrac{2x^3}{L^3}+\dfrac{3x^2}{L^2})\,v^2,\quad i=1,2, \quad p = 1,2,\\ N_p^2\alpha_2^p=(\dfrac{x^3}{L^2}-\dfrac{2x^2}{L}+x)\,\varphi^1 + (\dfrac{x^3}{L^2}-\dfrac{x^2}{L})\,\varphi^2,\quad \alpha_1^p=v^p, \quad \alpha_2^p=\varphi^p. \end{array}$

В случае одноосного деформированного состояния для изотропного материала из функционала Лагранжа можно получить СЛАУ относительно прогибов и углов поворота на концах для данного стержневого элемента:

$K_{pq}^{ij}\alpha_j^q=F_p^i,\quad K_{pq}^{ij}=EJ_3 \int_0^L N_{p,_{11}}^iN_{q,_{11}}^j dx, \quad F_p^i = \int_0^L [q(x) N_p^i + m(x) N_{p,_{1}}^i] dx + f_p^i,$ $N_{p,_1}=\dfrac{dN_p}{dx},\quad N_{p,_{11}}=\dfrac{d^2N_p}{dx^2},\quad F_p^1 = F_p, \quad F_p^2=M_p,$ $\dfrac{EJ_3}{L^3}\begin{bmatrix} 12&6L&-12&6L\\ 6L&4L^2&-6L&2L^2\\-12&-6L&12&-6L\\6L&2L^2&-6L&4L^2\end{bmatrix} \left\lbrace \begin{matrix}v^1\\ \varphi^1\\v^2\\ \varphi^2 \end{matrix}\right\rbrace=\left\lbrace \begin{matrix} F_1\\M_1\\F_2\\M_2 \end{matrix}\right\rbrace$

Здесь Fp, Mp - узловые силы и моменты, заданные изначально. При этом СЛАУ записаны в локальной системе координат в плоскости 0xy. Аналогичные уравнения можно выписать в плоскости 0zx. Для связи локальной и глобальной области обычно используется матрица ортогональных преобразований (для которой обратная матрица равна транспонированной), компоненты которой являются координатами i,j,k локального базиса в глобальной системе координат. Матрицу жёсткости К в глобальной системе координат можно получить из локальной, если выписать квадратичную часть функционала Лагранжа в терминах глобальных прогибов и поворотов:

$\left\lbrace \begin{matrix} \alpha_i^p = a_{im} \beta_m^p, \\ \alpha_j^q = a_{jn}\beta_j^q,\end{matrix}\right. \, \Rightarrow \,\alpha_i^p K_{pq}^{ij} \alpha_j^q = a_{im}\beta_m^p K_{pq}^{ij} a_{jn}\beta_n^q, \, \Rightarrow \mathbb{K}_{pq}^{mn}\beta_n^q = F_p^m$ $\mathbb{K}_{pq}^{mn} = a_{im}^T K_{pq}^{ij} a_{jn}, \quad a_{ij}=\begin{bmatrix}ix&iy&iz\\jx&jy&jz\\kx&ky&kz\end{bmatrix}.$

Таким образом можно составить глобальную СЛАУ и найти неизвестные узловые прогибы (перемещения) и повороты. В этом и заключается основная идея метода конечных элементов. Далее, зная узловые перемещения и повороты в глобальной системе координат можно с помощью матрицы (a) транспонированной перейти к локальным перемещениям (прогибам) и поворотам в узлах и восстановить аналитическую функцию прогиба, функцию момента и поперечной силы:

$M_3(x)=-EJ_3\dfrac{d^2v}{dx^2}, \quad Q(x) = -\dfrac{dM}{dx} - m(x).$

Кто более подробно интересуется построением данной математической модели может оставить комментарий.

В качестве примера представлена иллюстрация одной из публичной конечно-элементной модели с эпюрами прогибов, моментов и поперечных сил, найденных численно-аналитическим методом в рамках гипотезы Эйлера-Бернулли.

Эпюры прогибов Uy для пространственной рамы под действием крутящего момента M3, приложенного сверху. Если приложена линейная распределённая нагрузка q(x), то функция прогиба вдоль элемента в общем случае является полиномом 5-й степени.

Эпюры моментов для 3-элементной рамы. В общем случае вдоль стержня моменты представлены полиномом 3-й степени.

В следующем разделе речь пойдёт об особенности построения архитектуры SPA приложения, база которой представлена стеком React js, MobX, Three js...