menu
Anh-Thi DINH

NXFEM - idea and description (coding)

Posted on 19/09/2017, in PhD.

Cái này là phiên bản khác của cái ghi chú file nxfem_matlab_algorithm.pdf trong thư mục \matlab\readme.

Ý tưởng của việc code matlab với FEM

Làm sao từ bilinear form (dạng $a(u,v)=(f,v)$) cho đến AU=F? Có thể đọc trang 20 (1D) sách của Claes. Bên dưới là ghi ngắn gọn.

Mà $u=\sum_iu_i\varphi_i$ nên

Trong đó $A_{ij}=\int_{\Omega}\varphi_i\varphi_j = \Sigma_K\int_K\varphi_i\varphi_j$ và $F_j=\int_{\Omega}f\varphi_j = \Sigma_K\int_Kf\varphi_j$.

The main idea to implement FEM is that we will find the value of matrix A at nodes of each triangle and then take a sum of all triangles in the domain. Values at nodes with the same number will be added together.

Mesh trong matlab

  • Tam giác đánh số i-j-k theo thứ tự ngược chiều kim đồng hồ. Xem thêm mesh data.
  • Cái cut triangles CTs giống cấu trúc của msh.t nhưng thêm hàng thứ 5 chứa index trong msh.t

Ý tưởng code của NXFEM

If we are working on the standard FEM, we can consider at the same time all of triangles in the mesh because, basically, they are “he same”. However, it’s more difficult in NXFEM because it’s different in the cut triangles. We need to add “more nodes” (more basic functions) to the standard mesh.

Ý tưởng code $L^2$ norm

Có đoạn nói là áp dụng ý tưởng khi tính load vector. Xem lại ý tưởng tìm load vector. Để tính cái này chủ yếu cần tính từng lượng

Với NCTs (not cut triangles), ta tính bình thường, chỉ cẩn trọng tại những node around interface vì những node đó có thể là k(i).

Với CTs thì khó hơn do mỗi đỉnh có tới 2 basic functin (ik(i)) nên trong file getLoadCTs có tính FCT1, FCT2 tại mỗi đỉnh. Chủ yếu là dùng quadrature rule cho tam giác (getLoadPartTri), tứ giác thì lấy cái lớn trừ cái nhỏ (getLoadWholeTri-getLoadPartTri).

Còn $L^2$ có chút khác biệt ở CTs. Chúng ta chỉ có dạng $\varphi_i\varphi_j$ hoặc $\varphi_{k(i)}\varphi_{k(j)}$. Đúng là ý tưởng cũng tương tự cái load vector, có điều nếu load vector tính cho 1 cái j thì cái này tính cho 2 cái i,j.

(cần bổ sung thêm trong file pdf)

Normal vector & GradnPhi ($\nabla_n \varphi$)

Normal vector của đoạn $\overrightarrow{AB}$ là vector luôn nằm bên trái!

Mong muốn code 1 hàm GradnPhi để khi nhập cái đỉnh, tam giác, đoạn là nó cho ra kết quả!

Cần lưu ý là có sự khác nhau lớn giữa hai cái sau đây

Ý tưởng tính 3 cái $\eta_K, \eta_S, \zeta_S$

Ba cái này là dùng để dự đoán giá trị của $\gamma$ trong thesis của Barrau, trang 33. Cái ý tưởng này có nói ở trong bài NXFEM rồi, search để xem.

  • Hai cái $\eta_K, \eta_S$ đều phải xây dựng một ma trận riêng để tính giá trị của chúng.

  • Riêng cái $\zeta_S$ do trùng dạng với chuẩn

    chỉ khác cái hệ số $\gamma$ thôi nên có thể viết về cùng 1 form tính chuẩn của jump được. Còn cái hệ số đi kèm thì tính riêng giống như file getLambda vậy, cái đó cũng tính riêng.

  • Ý tưởng tính cái $\zeta_S$ là tính cái ma trận, tính cái này cũng giống như tính trong hàm getMatrixOnGamCTs, đó chính là file getMatNormJumpU.

Cái quadrature

Hàm số intGradnPhi, nếu xét $P^1$ thì có ngay công thức sau (từ thầy và cô), không cần dùng quadrature

Trong đó $M$ là trung điểm của đoạn $AB$ và $\varphi_i(x_0,y_0)=N_i(\hat{x}_0,\hat{y}_0)$, chứng minh như bên dưới

Theo đó $Q_x,Q_y$ được xác định theo sự tương ứng $i-(0,0),j-(1,0),k-(0,1)$. Ta sẽ đi chứng minh sự tương đồng giữa $\varphi_i$ và $N_i$ tại các node này, sau đó tương tự cho $j,k$. Có ghi chú trong viết bảng. Đã kiểm chứng bằng matlab là hai cách tính này GIỐNG NHAU!

Top