Метод скінченних елементів — структурний аналіз від перших принципів

1. Лінійна пружність та визначальні рівняння

Лінійна пружність припускає малі деформації, лінійний відгук матеріалу (закон Гука) та відсутність пластичної деформації. Рівняння рівноваги в 3D має вигляд:

2. Слабка форма та метод Гальоркіна

3. Функції форми та ізопараметричні елементи

4. Матриця жорсткості елемента K_e

5. Глобальне складання та граничні умови

6. Відновлення напружень та критерій фон Мізеса

7. Розв'язувач CST-елементів на JavaScript

// Мінімальний 2D МСЕ для плоского напруженого стану з CST-трикутниками
class FEM2D {
  constructor(nodes, elements, E, nu, thickness = 1) {
    this.nodes    = nodes;    // [[x,y], ...] — вузли
    this.elements = elements; // [[n0,n1,n2], ...] — елементи
    this.t = thickness;
    this.D = planeStressD(E, nu);
    const ndof = 2 * nodes.length;
    this.K = Array.from({length: ndof}, () => new Float64Array(ndof));
    this.f = new Float64Array(ndof);
  }

  assemble() {
    for (const [n0, n1, n2] of this.elements) {
      const [Ke, Be, A] = cstElement(this.nodes[n0], this.nodes[n1],
                                         this.nodes[n2], this.D, this.t);
      const dofs = [2*n0, 2*n0+1, 2*n1, 2*n1+1, 2*n2, 2*n2+1];
      for (let i = 0; i < 6; i++)
        for (let j = 0; j < 6; j++)
          this.K[dofs[i]][dofs[j]] += Ke[i][j];
    }
  }

  applyFixed(nodeIndex, axis) { // axis: 0=x, 1=y // закріпити
    const d = 2 * nodeIndex + axis;
    const n = this.K.length;
    for (let j = 0; j < n; j++) this.K[d][j] = this.K[j][d] = 0;
    this.K[d][d] = 1;
    this.f[d] = 0;
  }

  applyForce(nodeIndex, fx, fy) {
    this.f[2*nodeIndex]   += fx;
    this.f[2*nodeIndex+1] += fy;
  }

  solve() { // Наївне виключення Гаусса (у продакшені використовуйте розріджений розв'язувач)
    const A = this.K.map(r => Array.from(r));
    const b = Array.from(this.f);
    const n = b.length;
    for (let p = 0; p < n; p++) {
      for (let i = p+1; i < n; i++) {
        const f = A[i][p] / A[p][p];
        for (let j = p; j < n; j++) A[i][j] -= f*A[p][j];
        b[i] -= f*b[p];
      }
    }
    const u = new Float64Array(n);
    for (let i = n-1; i >= 0; i--) {
      u[i] = b[i];
      for (let j = i+1; j < n; j++) u[i] -= A[i][j]*u[j];
      u[i] /= A[i][i];
    }
    return u;
  }
}

function planeStressD(E, nu) {
  const c = E / (1 - nu*nu);
  return [[c, c*nu, 0], [c*nu, c, 0], [0, 0, c*(1-nu)/2]];
}

function cstElement([x1,y1], [x2,y2], [x3,y3], D, t) {
  const A = 0.5 * Math.abs((x2-x1)*(y3-y1) - (x3-x1)*(y2-y1));
  const b = [y2-y3, y3-y1, y1-y2];
  const c = [x3-x2, x1-x3, x2-x1];
  // B має розмір 3×6
  const B = [
    [b[0],0,  b[1],0,  b[2],0  ],
    [0,  c[0],0,  c[1],0,  c[2]],
    [c[0],b[0],c[1],b[1],c[2],b[2]]
  ].map(row => row.map(v => v / (2*A)));
  // Ke = B^T D B A t  (6×6)
  const DB = matMul(D, B);  // 3×6
  const Ke = matMul(transpose(B), DB).map(row => row.map(v => v*A*t));
  return [Ke, B, A];
}

8. Типи елементів та програмне забезпечення

Комерційні пакети МСЕ — Ansys, Abaqus, NASTRAN, LS-DYNA — додають контактну механіку, нелінійні матеріали, геометричну нелінійність, динамічний (модальний та перехідний) аналіз і мультифізичне поєднання. Серед відкритих альтернатив — FEniCS (Python, мова форм FEniCS), Code_Aster (загальний) та deal.II (шаблони C++ для адаптивних методів високого порядку).