% A script that contains a bunch of things to be computed.
% ReWT Aug. 2001 
% Modified 2005
% 
reload=1; % set outside
initialdeform=0;
Nloop = 1;
littleh=0.01;

fibequal=0;

if reload==1,
  
  filename = 'Total.FEM';
  Fe = readFEMfile(filename);

   fprintf('fiberdensities overwritten ... \n');

  if fibequal==1,
    for k=1:length(Fe.fiberdensity),
      Fe.fiberdensity{k} = ones(size(Fe.fiberdensity{k}));
    end;
  end;
  
 fprintf('FEMMatrices ... takes loooooong\n');

  Fe = FEMMatrices(Fe);
  
  Fe.U = zeros(Fe.nnodes,3);

  con = Fe.constr;
  con = (con==0) | (con>2);
end;

% [Fe.cindxbar,Fe.dirbar,Fe.multibar]=fourthindicies;
% [Fe.cindxlag,Fe.dirlag,Fe.multilag]=lagrangeindicies;

% deformation in x dir: 
initialdeform=1;
if initialdeform==1,
Fd = [1 0.0   0.0
      0 1.0   0.0;
      0 0.0   1.1];

Fd = [    0.9400   -0.0800   -0.5900;
	 -0.3800    1.2400   -0.0600;
	 -0.2100   -0.8200    1.4200];

uu = Fe.xyz;
uu = (Fd*uu')'-Fe.xyz;
Fe.U = uu; % uu.*con;
end;

El = 10     ;
pois = 0.49;
Fe.lamda = pois*El/((1+pois)*(1-2*pois));
Fe.mu = El/(2*(1+pois));
Fe.kappa = 10 ;

activation=-1;

Nloop = 1