%--------------------------
 % @Author: Jingqiao Hu
 % @Date: 2022-05-05 21:23:30
 % @LastEditTime: 2022-05-06 19:53:43
%--------------------------
dbstop if error
feature accel on
addpath util util_simulation util_bezier util_micro util_mtimesx util_output
addpath util_cvt mexCVT util_tetmesh util_MMC util_intersection

%% init
datapath = 1;
[nodes, valid_tets] = read_tet(datapath);

minx = min(nodes(:,1));
maxx = max(nodes(:,1));
miny = min(nodes(:,2));
maxy = max(nodes(:,2));
minz = min(nodes(:,3));
maxz = min(nodes(:,3));

optDesign = 'halfsp';
optChanged = 0;
vF = 0.1;
[freedofs, fixeddofs, fext, alldofs] = design_domain(optDesign, ...
    maxx, minx, maxy, miny, maxz, minz, vF, nodes, optChanged)
U = zeros(alldofs, 1);

%% simulation
[c_sK, c_vol, vol0] = prepare_sK(c_nodesFine, c_tets, 1, nu0);

for ele = 1 : nele
    
end

%% utils
% tet file
% 1. Number of vertex, tet: n1, n2
% vertex:   n * 3
% tet:      m * 4
function [nodes, valid_tets] = read_tet(datapath)
    
    % datapath = [path0+'tet'+num2str(i)+'.txt']; 
    fileID = fopen(datapath, 'r');
    info   = fscanf(fileID, '%f');

    npnts  = info(1);
    ntets  = info(2); 

    nodes = reshape(info(3 : 3 + npnts*3 - 1), 3, [])'; % all nodes coodinates of cdt
    tets  = reshape(info(3 + npnts*3 : end), 4, [])' + 1;

    fclose(fileID);

    % [valid_tets] = checkValid_tetmesh(nodes, tets);
%     DT = triangulation(valid_tets, nodes(:,1), nodes(:,2), nodes(:,3));
%     figure; tetramesh(DT,'FaceAlpha',0.3);
end