Вектор нагрузки 1D

Теперь вычислим вектор нагрузки для одномерного случая. Вектор нагрузки связан с интегралом произведения функции источника на пробную функцию. Рассмотрим интеграл для одного отрезка с вершинами xi,xi+1x_i, x_{i+1}

∫xixi+1f(x)⋅υ dx.\int_{x_i}^{x_{i+1}} f(x) \cdot \upsilon \,dx.
(6.36)

Пробная функция на отрезке имеет вид υ(x)=qi⋅ϕi(x)+qi+1⋅ϕi+1(x)\upsilon(x) = q_i \cdot \phi_i(x) + q_{i+1} \cdot \phi_{i+1}(x). Для вычисления интеграла необходимо интерполировать функцию источника f(x)f(x) на отрезке [xi,xi+1][x_i, x_{i+1}] линейной функцией. Представим f(x)f(x) в виде f~(x)=e1+e2⋅x\widetilde{f}(x) = e_1 + e_2 \cdot x. Коэффициенты e1e_1 и e2e_2 определяются из условий интерполяции в узлах

{e1+e2⋅xi=fie1+e2⋅xi+1=fi+1\begin{cases} e_1 + e_2 \cdot x_i = f_i\\ e_1 + e_2 \cdot x_{i+1} = f_{i+1} \end{cases}
(6.37)

Решая эту систему по аналогии с (5.8), получаем

{e1=fi⋅xi+1−fi+1⋅xixi+1−xie2=fi+1−fixi+1−xi\begin{cases} e_1 = \frac{\displaystyle f_i \cdot x_{i+1} - f_{i+1} \cdot x_i}{\displaystyle x_{i+1} - x_i}\\ e_2 = \frac{\displaystyle f_{i+1} - f_i}{\displaystyle x_{i+1} - x_i} \end{cases}
(6.38)

Подставим интерполированную функцию источника и пробную функцию в (6.36). Примем во внимание (5.7) и (5.8) из раздела о функциях <<крышек>>

∫xixi+1f(x)⋅υ dx=∫xixi+1(e1+e2⋅x)⋅[qi⋅(ai+bi⋅x)+qi+1⋅(ai+1+bi+1⋅x)] dx.\int_{x_i}^{x_{i+1}} f(x) \cdot \upsilon \,dx = \int_{x_i}^{x_{i+1}} (e_1 + e_2 \cdot x) \cdot \left[ q_i \cdot (a_i + b_i \cdot x) + q_{i+1} \cdot (a_{i+1} + b_{i+1} \cdot x) \right] \,dx.

Раскроем скобки и разделим интеграл на части

∫xixi+1f(x)⋅υ dx=∫xixi+1[e1⋅(qi⋅ai+qi+1⋅ai+1)+e1⋅(qi⋅bi+qi+1⋅bi+1)⋅x+e2⋅(qi⋅ai+qi+1⋅ai+1)⋅x+e2⋅(qi⋅bi+qi+1⋅bi+1)⋅x2] dx\begin{aligned} &\int_{x_i}^{x_{i+1}} f(x) \cdot \upsilon \,dx = \int_{x_i}^{x_{i+1}} \Big[ e_1 \cdot (q_i \cdot a_i + q_{i+1} \cdot a_{i+1}) + e_1 \cdot (q_i \cdot b_i + q_{i+1} \cdot b_{i+1}) \cdot x\\ &+ e_2 \cdot (q_i \cdot a_i + q_{i+1} \cdot a_{i+1}) \cdot x + e_2 \cdot (q_i \cdot b_i + q_{i+1} \cdot b_{i+1}) \cdot x^2 \Big] \,dx \end{aligned}

Вычислим каждый из интегралов

∫xixi+1f(x)⋅υ dx=e1⋅(qi⋅ai+qi+1⋅ai+1)⋅x∣xixi+1+[e1⋅(qi⋅bi+qi+1⋅bi+1)+e2⋅(qi⋅ai+qi+1⋅ai+1)]⋅x22∣xixi+1+e2⋅(qi⋅bi+qi+1⋅bi+1)⋅x33∣xixi+1\begin{aligned} &\int_{x_i}^{x_{i+1}} f(x) \cdot \upsilon \,dx = e_1 \cdot (q_i \cdot a_i + q_{i+1} \cdot a_{i+1}) \cdot x \bigg|_{x_i}^{x_{i+1}}\\ &+ \left[ e_1 \cdot (q_i \cdot b_i + q_{i+1} \cdot b_{i+1}) + e_2 \cdot (q_i \cdot a_i + q_{i+1} \cdot a_{i+1}) \right] \cdot \frac{\displaystyle x^2}{\displaystyle 2} \bigg|_{x_i}^{x_{i+1}}\\ &+ e_2 \cdot (q_i \cdot b_i + q_{i+1} \cdot b_{i+1}) \cdot \frac{\displaystyle x^3}{\displaystyle 3} \bigg|_{x_i}^{x_{i+1}} \end{aligned}

Подставляя выражения для e1,e2e_1, e_2 из (6.38) и используя соотношения (5.8), после упрощения получаем

∫xixi+1f(x)⋅υ dx=(fi+fi+1)⋅(xi+1−xi)2⋅(qi⋅ai+qi+1⋅ai+1)+xi+1−xi6⋅[fi⋅xi+fi+1⋅xi+1+(fi+fi+1)⋅(xi+xi+1)]⋅(qi⋅bi+qi+1⋅bi+1)\begin{aligned} &\int_{x_i}^{x_{i+1}} f(x) \cdot \upsilon \,dx = \frac{\displaystyle (f_i + f_{i+1}) \cdot (x_{i+1} - x_i)}{\displaystyle 2} \cdot (q_i \cdot a_i + q_{i+1} \cdot a_{i+1})\\ &+ \frac{\displaystyle x_{i+1} - x_i}{\displaystyle 6} \cdot \left[ f_i \cdot x_i + f_{i+1} \cdot x_{i+1} + (f_i + f_{i+1}) \cdot (x_i + x_{i+1}) \right] \cdot (q_i \cdot b_i + q_{i+1} \cdot b_{i+1}) \end{aligned}

Дальнейшие преобразования с учётом (5.8) приводят к окончательному результату

∫xixi+1f(x)⋅υ dx=xi+1−xi6⋅[qi⋅(2⋅fi+fi+1)+qi+1⋅(fi+2⋅fi+1)].\int_{x_i}^{x_{i+1}} f(x) \cdot \upsilon \,dx = \frac{\displaystyle x_{i+1} - x_i}{\displaystyle 6} \cdot \left[ q_i \cdot (2 \cdot f_i + f_{i+1}) + q_{i+1} \cdot (f_i + 2 \cdot f_{i+1}) \right].

Раскроем и перегруппируем члены

∫xixi+1f(x)⋅υ dx=L⋅[2⋅fi+fi+16⋅qi+fi+2⋅fi+16⋅qi+1],\int_{x_i}^{x_{i+1}} f(x) \cdot \upsilon \,dx = L \cdot \left[ \frac{\displaystyle 2 \cdot f_i + f_{i+1}}{\displaystyle 6} \cdot q_i + \frac{\displaystyle f_i + 2 \cdot f_{i+1}}{\displaystyle 6} \cdot q_{i+1} \right],

где LL — длина отрезка, которая вычисляется по формуле

L=xi+1−xi.L = x_{i+1} - x_i.
(6.39)

Введём обозначения для элементов локального вектора нагрузки отрезка

ri=L6⋅(2⋅fi+fi+1)ri+1=L6⋅(fi+2⋅fi+1)\begin{split} &r_i = \frac{\displaystyle L}{\displaystyle 6} \cdot (2 \cdot f_i + f_{i+1})\\ &r_{i+1} = \frac{\displaystyle L}{\displaystyle 6} \cdot (f_i + 2 \cdot f_{i+1}) \end{split}
(6.40)

Таким образом, локальный вектор нагрузки для одномерного элемента имеет вид

R=[riri+1]=L6[2⋅fi+fi+1fi+2⋅fi+1].\begin{gathered}\mathbf{R} = \begin{bmatrix} r_i\\ r_{i+1} \end{bmatrix} = \frac{\displaystyle L}{\displaystyle 6} \begin{bmatrix} 2 \cdot f_i + f_{i+1}\\ f_i + 2 \cdot f_{i+1} \end{bmatrix}.\end{gathered}
(6.41)

Глобальный вектор нагрузки R\mathbf{R} получается путём суммирования вкладов от всех отрезков сетки методом сборки: элементы локальных векторов добавляются к соответствующим элементам глобального вектора согласно глобальной нумерации узлов. Размерность глобального вектора нагрузки равна NN, где NN — общее количество узлов сетки. Важно отметить, что для внутренних узлов происходит суммирование вкладов от двух смежных элементов, в результате чего элемент глобального вектора нагрузки для внутреннего узла jj принимает вид

Rj=Lj−16⋅(fj−1+2⋅fj)+Lj6⋅(2⋅fj+fj+1),R_j = \frac{\displaystyle L_{j-1}}{\displaystyle 6} \cdot (f_{j-1} + 2 \cdot f_j) + \frac{\displaystyle L_j}{\displaystyle 6} \cdot (2 \cdot f_j + f_{j+1}),
(6.42)

где Lj−1=xj−xj−1L_{j-1} = x_j - x_{j-1} и Lj=xj+1−xjL_j = x_{j+1} - x_j — длины смежных отрезков. Для граничных узлов вектор нагрузки определяется граничными условиями задачи.

Матрица демпфирования 3DВектор нагрузки 2D