crack.m

来自「C++编译指导」· M 代码 · 共 31 行

function crack
clear all
% Initial Mesh
node = [1,0; 0,1; -1,0; 0,-1; 0,0; 1,0];
elem = [5,1,2; 5,2,3; 5,3,4; 5,4,6];
Dirichlet = [1,2; 2,3; 3,4; 4,6; 1,5; 5,6];
Neumann = [];
% Adaptive Finite Element Method
for step = 1:28
    % Step1: Solver
    u = Poisson(node, elem, Dirichlet, Neumann, @f, @u_D, @g);
    % Step2: Estimate
    eta = estimate(node,elem,u);
    % Step3: Refine
    [node,elem,marker,d2p] = bisection(node,elem,eta.^2,0.4,0);
    % refine boundary edges
    Dirichlet = updatebd(Dirichlet,marker,d2p);
    Neumann = updatebd(Neumann,marker,d2p);
end
u = Poisson(node, elem, Dirichlet, Neumann, @f, @u_D, @g);

% -------------------------------------------------------------------------
function z = f(p) % data
z = 1;

function z = u_D(p) % Dirichlet boundary condition
r = sqrt(sum(p.^2,2));
z = sqrt(0.5*(r-p(:,1)))-0.25*r.^2;

function z = g(p) % Neumann boundary condition
z = 0;