%---------------------------------------------------------
% B_BUILD_MODEL.M
% building the openFEM computer model
% Fiala Peter
% 2008 BME, Laboratory of Acoustics
%---------------------------------------------------------
clear all;
close all;
clc;
%% loading parameters
load data/parameters.mat
%% element materials and properties
% materials
pl = [
%   MatID type                       E    nu    rho    G    eta
    1     fe_mat('m_elastic','SI',1) Este nuste rhoste Gste 0
];
midste = 1; % material id for steel
% properties
il = [
%   ProID type                    J    I1    I2    A
    1     fe_mat('p_beam','SI',1) Jbea Ibea  Ibea  Abea % column property
%   ProID type                     f d 0 w
    2     fe_mat('p_shell','SI',1) 0 0 0 wpla % plate property
];
pidbea = 1; % property id for beam element
pidpla = 2; % property id for plate element
%% Initialisation
FEnode = [];
FEelt = [];
%% x beams
nodes = [
    0  0  0
    Lx 0  0
];
[FEnode, nodeind] = femesh('AddNode', FEnode, nodes);
FEel0 = [
    Inf abs('beam1') 0
%   n1 n2     MatID  ProID  vx  vy  vz
    nodeind.' midste pidbea 0.0 1.0 0.0
];
femesh(['Divide ' num2str(Nex)]);
femesh(['RepeatSel 2 ' num2str([0 Ly 0])]);
femesh(['RepeatSel 2 ' num2str([0 0  Lz])]);
femesh('AddSel');
%% y beams
nodes = [
    0  0  0
    0  Ly 0
];
[FEnode, nodeind] = femesh('AddNode', FEnode, nodes);
FEel0 = [
    Inf abs('beam1') 0
%   n1 n2     MatID  ProID  vx  vy  vz
    nodeind.' midste pidbea 1.0 0.0 0.0
];
femesh(['Divide ' num2str(Ney)]);
femesh(['RepeatSel 2 ' num2str([Lx 0 0])]);
femesh(['RepeatSel 2 ' num2str([0 0  Lz])]);
femesh('AddSel');
%% z beams 
nodes = [
    0 0 0
    0 0 Lz
];
[FEnode, nodeind] = femesh('AddNode', FEnode, nodes);
FEel0 = [
    Inf abs('beam1') 0
%   n1 n2     MatID  ProID  vx  vy  vz
    nodeind.' midste pidbea 1.0 0.0 0.0
];
femesh(['Divide ' num2str(Nez)]);
femesh(['RepeatSel 2 ' num2str([Lx 0  0])]);
femesh(['RepeatSel 2 ' num2str([0  Ly 0])]);
femesh('AddSel');
%% XY walls 
nodes = [
    0 0  0
    Lx 0  0
    Lx Ly 0
    0 Ly 0
];
[FEnode, nodeind] = femesh('AddNode', FEnode, nodes);
FEel0 = [
    Inf abs('quad4')
%   n1 n2 n3 n4 MatID  ProID
    nodeind.'   midste pidpla
];
femesh(['Divide ' num2str([Ney Nex])]);
femesh(['RepeatSel 2 ' num2str([0 0 Lz])]);
femesh('AddSel');
%% YZ walls 
nodes = [
    0 0  0
    0 0  Lz
    0 Ly Lz
    0 Ly 0
];
[FEnode, nodeind] = femesh('AddNode', FEnode, nodes);
FEel0 = [
    Inf abs('quad4')
%   n1 n2 n3 n4 MatID  ProID
    nodeind.'   midste pidpla
];
femesh(['Divide ' num2str([Ney Nez])]);
femesh(['RepeatSel 2 ' num2str([Lx 0 0])]);
femesh('AddSel');
%% XZ walls 
nodes = [
    0 0  0
    Lx 0  0
    Lx 0 Lz
    0 0 Lz
];
[FEnode, nodeind] = femesh('AddNode', FEnode, nodes);
FEel0 = [
    Inf abs('quad4')
%   n1 n2 n3 n4 MatID  ProID
    nodeind.'   midste pidpla
];
femesh(['Divide ' num2str([Nez Nex])]);
femesh(['RepeatSel 2 ' num2str([0 Ly 0])]);
femesh('AddSel');
%% Definition of the SDT model
box = femesh('Model');
box.il = il;
box.pl = pl;
DOF = feutil('GetDOF', box);
%
save data/boxmesh.mat box DOF