DEMO_febio_0050_foot_insole_01

Below is a demonstration for:

Contents

Keywords

clear; close all; clc;

Plot settings

fontSize=15;
faceAlpha1=1;
faceAlpha2=0.3;
markerSize1=15;
markerSize2=10;
lineWidth=2;

Control parameters

% Path names
defaultFolder = fileparts(fileparts(mfilename('fullpath')));
savePath=fullfile(defaultFolder,'data','temp');
loadPathSurfaces=fullfile(defaultFolder,'data','STL','leg','post');

% Defining file names
febioFebFileNamePart='tempModel';
febioFebFileName=fullfile(savePath,[febioFebFileNamePart,'.feb']); %FEB file name
febioLogFileName=fullfile(savePath,[febioFebFileNamePart,'.txt']); %FEBio log file name
febioLogFileName_disp=[febioFebFileNamePart,'_disp_out.txt']; %Log file name for exporting displacement
febioLogFileName_force=[febioFebFileNamePart,'_force_out.txt']; %Log file name for exporting force
febioLogFileName_strainEnergy=[febioFebFileNamePart,'_energy_out.txt']; %Log file name for exporting strain energy density

% Surface model file names
fileNameBones=fullfile(loadPathSurfaces,'Foot_bulk.stl');
fileNameSkin=fullfile(loadPathSurfaces,'Skin_coarse.stl');

%Geometric parameters
soleOffsetOutward=3; %How much the insole protrudes outward from the foot
numSmoothIterations_sole_Z=15; %Number of smoothing iterations for the sole Z value
soleMinThickness=6; %Sole thickness
volumeFactor=1; %Volume factor used in tetgen, larger means larger internal elements
footSize=251; %Foot size in mm size 6=241 mm, size 7=251 mm

% Material parameters
%-> Tissue (Ogden, uncoupled)
% Data source: Petre et al. 2013, Optimization of Nonlinear Hyperelastic
% Coefficients for Foot Tissues Using a Magnetic Resonance Imaging
% Deformation Experiment. https://dx.doi.org/10.1115%2F1.4023695
c1_1=0.02652; %Shear-modulus-like parameter in MPa
m1_1=17.71; %Material parameter setting degree of non-linearity
k_1=c1_1*100; %Bulk modulus

%-> Sole
%Sole material parameters (Ogden unconstrained, with m=2 i.e. Neo-Hookean)
c1_min=0.05; %Shear-modulus-like parameter in MPa
c1_max=0.1; %Shear-modulus-like parameter in MPa
numMaterialLevels=50; %Max number of materials to use
initialMaterialLevel=numMaterialLevels;
c1_range=linspace(c1_min,c1_max,numMaterialLevels); %Shear-modulus-like parameter
bulkModulusFactor=10;

testCase=1; % Change to test different sole cases, e.g. stiff, soft, spatially varying.
maxLevelColorbar_SED=0.0005;

% FEA control settings
numTimeSteps=10; %Number of time steps desired
max_refs=25; %Max reforms
max_ups=0; %Set to zero to use full-Newton iterations
opt_iter=10; %Optimum number of iterations
max_retries=8; %Maximum number of retires
dtmin=(1/numTimeSteps)/100; %Minimum time step size
dtmax=1/numTimeSteps; %Maximum time step size
symmetric_stiffness=0;
min_residual=1e-20;
runMode='internal';

initialSoleSpacing=0.1;

%Boundary condition parameters
displacementMagnitude=-2.63; %Displacement applied to the bones
bodyWeight=65/2;% Body weight in Kg
forceBody=bodyWeight.*-9.81; %Force due to body weight
forceDifferenceToleranceFraction=0.05;
forceDifferenceTolerance=abs(forceBody.*forceDifferenceToleranceFraction);
optimizeForceOption=1; %set to 0 to turn off and use input displacement
maxProposedDisplacementDifference=1;

%Contact parameters
contactPenalty=10;
laugon=0;
minaug=1;
maxaug=10;
fric_coeff=0.25;

Import surface models

%Import STL file for the bones
[stlStruct] = import_STL(fileNameBones);
F1=stlStruct.solidFaces{1}; %Faces
V1=stlStruct.solidVertices{1}; %Vertices
[F1,V1]=mergeVertices(F1,V1); %Merge nodes

%Import STL file for the skin
[stlStruct] = import_STL(fileNameSkin);
F2=stlStruct.solidFaces{1}; %Faces
V2=stlStruct.solidVertices{1}; %Vertices
[F2,V2]=mergeVertices(F2,V2); %Merge nodes
V2=V2*100; %Scale

Subdivide skin mesh and smooth

% %Split elements once
% [F2,V2]=subtri(F2,V2,1);
%
% %Smooth
% clear cParSmooth;
% cParSmooth.n=5;
% cParSmooth.Method='HC';
% [V2]=patchSmooth(F2,V2,[],cParSmooth);

Visualize imported surfaces

cFigure; hold on;
gpatch(F1,V1,'bw','none',1);
gpatch(F2,V2,'rw','none',0.5);
camlight('headlight');
axisGeom(gca,fontSize);
drawnow;

Compute mesh derived parameters

Compute point spacings for the surfaces meshes. These are useful to relate other mesh sizes to.

pointSpacing1=mean(patchEdgeLengths(F1,V1));
pointSpacing2=mean(patchEdgeLengths(F2,V2));
pointSpacing=mean([pointSpacing1 pointSpacing2]);
pointSpacingSole=pointSpacing;

Reorient surfaces

Reorient so that the leg points up in a forward view with the toes pointing towards the viewer.

R=euler2DCM([-0.5*pi 0 -0.5*pi]);
V1=V1*R;
V2=V2*R;

Visualize reoriented surfaces

cFigure; hold on;
gpatch(F1,V1,'bw','none',1);
gpatch(F2,V2,'rw','none',0.5);
camlight('headlight');
axisGeom(gca,fontSize);
drawnow;

Cut skin surface

The skin surface is cut such that it stops where the bones of the ankle end in the z direction.

%Create a logic for cutting away faces
max_Z1=max(V1(:,3))+2*pointSpacing; %Max Z-level used for cutting
logicVertices=V2(:,3)<max_Z1; %Logic for the points below this level
logicFaces=all(logicVertices(F2),2); %Logic for the faces
logicFaces=triSurfLogicSharpFix(F2,logicFaces,3); %Altered logic so it is smoother

Visualize

cFigure; hold on;
gpatch(F1,V1,'bw','none',1);
gpatch(F2,V2,logicFaces,'none',0.5);
camlight('headlight');
axisGeom(gca,fontSize);
drawnow;

Cut away faces using logic

F2=F2(logicFaces,:); %The faces to keep
[F2,V2]=patchCleanUnused(F2,V2); %Remove unused points

%Attempt to self triangulate potentially jagged edge
Eb=patchBoundary(F2,V2); %Get boundary edges
indBoundary=edgeListToCurve(Eb); %Convert boundary edges to a curve list
indBoundary=indBoundary(1:end-1); %Trim off last point since it is equal to first on a closed loop
angleThreshold=pi*(120/180); %threshold for self triangulation
[F2,V2,indBoundaryTop]=triSurfSelfTriangulateBoundary(F2,V2,indBoundary,angleThreshold,1);

%Force boundary to have the max Z level chosen
V2(indBoundaryTop,3)=max_Z1;

Visualize

cFigure; hold on;
gpatch(F1,V1,'bw','none',1);
gpatch(F2,V2,'rw','k',0.5);
plotV(V2(indBoundaryTop,:),'r-','LineWidth',lineWidth);
camlight('headlight');
axisGeom(gca,fontSize);
drawnow;

Scale foot to desired size

currentFootSize=abs(max(V2(:,1))-min(V2(:,1)));
V2=V2./currentFootSize; %Scale to unit foot
V2=V2.*footSize; %Scale to desired foot size
V1=V1./currentFootSize; %Scale to unit foot
V1=V1.*footSize; %Scale to desired foot size

Close over top of skin

The top boundary curve of the cut surface is filled with triangles. This is a 2D method. The z-coordinate is added after.

[F2t,V2t]=regionTriMesh2D({V2(indBoundaryTop,[1 2])},pointSpacing2,0,0);
V2t(:,3)=mean(V2(indBoundaryTop,3)); %Add/set z-level

Visualize

cFigure; hold on;
gpatch(F1,V1,'bw','none',1);
gpatch(F2,V2,'rw','k',0.5);
gpatch(F2t,V2t,'gw','k',1);
plotV(V2(indBoundaryTop,:),'r-','LineWidth',lineWidth);
camlight('headlight');
axisGeom(gca,fontSize);
drawnow;

Joining surface features

Add all surface sets together in joint list of faces, vertices

[FT,VT,CT]=joinElementSets({F1,F2,F2t},{V1,V2,V2t});

Merge shared nodes

The join operation only adds the sets together. Nodes with the same coordinates are not seen as the same yet and need to be merged.

[FT,VT]=mergeVertices(FT,VT); %Merge nodes

Visualize

cFigure; hold on;
gpatch(FT,VT,CT,'k',0.5);
patchNormPlot(FT,VT);
camlight('headlight');
axisGeom(gca,fontSize);
colormap gjet; icolorbar;
drawnow;

Mesh foot with tetrahedral elements

Tet meshing is based on tetgen. TetGen requires a interior points for regions to be meshed, as well as intertior points for holes.

Define region points

[V_region]=getInnerPoint({FT(CT==2 | CT==3,:),FT(CT==1,:)},{VT,VT});
[V_hole]=getInnerPoint(FT(CT==1,:),VT);

Visualize interior points

cFigure; hold on;
gpatch(FT,VT,'kw','none',0.2);
hp1=plotV(V_region,'r.','markerSize',markerSize1);
hp2=plotV(V_hole,'b.','markerSize',markerSize1);
legend([hp1 hp2],{'Region point','Hole point'});
camlight('headlight');
axisGeom(gca,fontSize);
drawnow;

Mesh using tetgen

inputStruct.stringOpt='-pq1.2AaY'; %TetGen option string
inputStruct.Faces=FT; %The faces
inputStruct.Nodes=VT; %The vertices
inputStruct.holePoints=V_hole; %The hole interior points
inputStruct.faceBoundaryMarker=CT; %Face boundary markers
inputStruct.regionPoints=V_region; %The region interior points
inputStruct.regionA=tetVolMeanEst(F2,V2)*volumeFactor; %Volume for regular tets

Mesh model using tetrahedral elements using tetGen

[meshOutput]=runTetGen(inputStruct); %Run tetGen
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
--- TETGEN Tetrahedral meshing --- 13-Mar-2020 12:22:39
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
--- Writing SMESH file --- 13-Mar-2020 12:22:39
----> Adding node field
----> Adding facet field
----> Adding holes specification
----> Adding region specification
--- Done --- 13-Mar-2020 12:22:39
--- Running TetGen to mesh input boundary--- 13-Mar-2020 12:22:39
Opening /mnt/data/MATLAB/GIBBON/data/temp/temp.smesh.
Delaunizing vertices...
Delaunay seconds:  0.039103
Creating surface mesh ...
Surface mesh seconds:  0.011248
Recovering boundaries...
Boundary recovery seconds:  0.027423
Removing exterior tetrahedra ...
Spreading region attributes.
Exterior tets removal seconds:  0.017467
Recovering Delaunayness...
Delaunay recovery seconds:  0.012961
Refining mesh...
Refinement seconds:  0.153073
Optimizing mesh...
Optimization seconds:  0.009416

Writing /mnt/data/MATLAB/GIBBON/data/temp/temp.1.node.
Writing /mnt/data/MATLAB/GIBBON/data/temp/temp.1.ele.
Writing /mnt/data/MATLAB/GIBBON/data/temp/temp.1.face.
Writing /mnt/data/MATLAB/GIBBON/data/temp/temp.1.edge.

Output seconds:  0.141619
Total running seconds:  0.412759

Statistics:

  Input points: 6540
  Input facets: 13072
  Input segments: 19608
  Input holes: 1
  Input regions: 1

  Mesh points: 10963
  Mesh tetrahedra: 49003
  Mesh faces: 104542
  Mesh faces on exterior boundary: 13072
  Mesh faces on input facets: 13072
  Mesh edges on input segments: 19608
  Steiner points inside domain: 4423

--- Done --- 13-Mar-2020 12:22:40
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
--- Importing TetGen files --- 13-Mar-2020 12:22:40
--- Done --- 13-Mar-2020 12:22:40

Access model element and patch data

Fb_foot=meshOutput.facesBoundary; %Boundary faces of the foot
Cb_foot=meshOutput.boundaryMarker; %Boundary marker/color data for the foot
V_foot=meshOutput.nodes; %The vertices/nodes
E_foot=meshOutput.elements; %The tet4 elements

Visualizing mesh using meshView, see also anim8

meshView(meshOutput);

Build insole

The insole is created by taking the 2D convex hull of the foot (ignoring the z-direction). Next this convex hull is resampled and filled with triangular elements. The z-coordinates are then based on the nearest foot nodes. The z-coordinate data is next smoothed to create a smooth surface. This surface is next thickened to form hexahedral elements.

Create sole boundary curve in 2D

% Offset surface outward to thicken so sole is enlarged outward
[~,~,Nv]=patchNormal(FT,VT);
VT_sole=VT+soleOffsetOutward.*Nv;

P=VT_sole(:,[1 2]); %Point set flattened to 2D
DT=delaunayTriangulation(P); %Delaunay triangulation of 2D set
VD=DT.Points; %Delaunay point set
VD(:,3)=min(VT(:,3)); %Set z-coord to minimum for now
indChull=DT.convexHull; %Ordered point list for conhex hull
indChull=indChull(1:end-1); %Trim away last (=start) point to avoid double
V_chull=VD(indChull,:); %Vertices for convex hull
D=max(pathLength(V_chull)); %Get length of sole curve for resampling
numResample=ceil(D./pointSpacingSole);
V_sole_curve=evenlySampleCurve(V_chull,numResample,'pchip',1);

Visualize

cFigure; hold on;
gpatch(FT,VT,'kw','none',0.25);
plotV(P,'k.','MarkerSize',markerSize2);
plotV(V_chull,'r.','MarkerSize',markerSize1*2);
plotV(V_sole_curve,'b.-','MarkerSize',markerSize1,'LineWidth',lineWidth);
camlight('headlight');
axisGeom(gca,fontSize);
colormap gjet; icolorbar;
drawnow;

Build sole top surface

[F_sole_top,V_sole_top]=regionTriMesh2D({V_sole_curve(:,[1 2])},pointSpacingSole,0,0);
V_sole_top(:,3)=mean(V_sole_curve(:,3));

Eb=patchBoundary(F_sole_top,V_sole_top);
indBoundary=unique(Eb(:));

%Get z-coordinate
[~,indMin]=minDist(V_sole_top,VT);
V_sole_top(:,3)=VT(indMin,3);

%Free smoothing of boundary
clear cParSmooth;
cParSmooth.n=numSmoothIterations_sole_Z;
cParSmooth.Method='LAP';
[p]=patchSmooth(Eb,V_sole_top,[],cParSmooth);
V_sole_top(indBoundary,3)=p(indBoundary,3);

% Constrained general smoothing
cParSmooth.n=numSmoothIterations_sole_Z;
cParSmooth.Method='LAP';
cParSmooth.RigidConstraints=indBoundary;
[V_sole_top(:,3)]=patchSmooth(F_sole_top,V_sole_top(:,3),[],cParSmooth);

%Shift a small amount to there is no initial contact
[~,indMin]=minDist(V_sole_top,VT);
dz=VT(indMin,3)-V_sole_top(:,3);
V_sole_top(:,3)=V_sole_top(:,3)-abs(min(dz));

V_sole_top(:,3)=V_sole_top(:,3)-initialSoleSpacing;

Visualize

cFigure; hold on;
gpatch(FT,VT,'kw','none',0.25);
gpatch(F_sole_top,V_sole_top,'gw','k',1);
patchNormPlot(F_sole_top,V_sole_top);
camlight('headlight');
axisGeom(gca,fontSize);
colormap gjet; icolorbar;
drawnow;

Create bottom surface of sole

dirSet=[0 0 -1];
numElementsSoleThickness=ceil(soleMinThickness./pointSpacingSole);
[E_sole,V_sole,F_sole_top,F_sole_bottom]=patchThick(fliplr(F_sole_top),V_sole_top,dirSet,soleMinThickness,numElementsSoleThickness);
F_sole_top=fliplr(F_sole_top);

% Use element2patch to get patch data
F_E_sole=element2patch(E_sole,[],'penta6');
% indBoundaryFaces=tesBoundary(F_E_sole,V_sole);
% Fb_sole=F_E_sole(indBoundaryFaces,:);
cFigure; hold on;
gpatch(FT,VT,'kw','none',0.25);
gpatch(F_sole_top,V_sole,'gw','k',1);
gpatch(F_sole_bottom,V_sole,'bw','k',1);
patchNormPlot(F_sole_top,V_sole);
camlight('headlight');
axisGeom(gca,fontSize);
colormap gjet; icolorbar;
drawnow;

Visualize

cFigure; hold on;
title('Hexahedral mesh');
gpatch(FT,VT,'kw','none',0.25);
gpatch(F_E_sole,V_sole,'bw','k',1);
% patchNormPlot(Fb_sole,V_sole);
axisGeom;
camlight headlight;
drawnow;

Joining node sets

V=[V_foot;V_sole;]; %Combined node sets
E_sole=E_sole+size(V_foot,1); %Fixed element indices
F_sole_top=F_sole_top+size(V_foot,1); %Fixed element indices
F_sole_bottom=F_sole_bottom+size(V_foot,1); %Fixed indices
F_E_sole=element2patch(E_sole,[],'penta6');

Visualize

cFigure; hold on;
title('Hexahedral mesh');
gpatch(Fb_foot,V,Cb_foot,'k',1);
gpatch(F_sole_top,V,'kw','k',1);
gpatch(F_sole_bottom,V,'kw','k',1);
% patchNormPlot(FEs,V);
colormap gjet; icolorbar;
axisGeom;
camlight headlight;
drawnow;

Define contact surfaces

% The rigid master surface of the sphere
F_contact_master=F_sole_top;

% The deformable slave surface of the slab
F_contact_slave=fliplr(Fb_foot(Cb_foot==2,:));

Visualize contact surfaces

cFigure; hold on;
title('Contact sets and normal directions','FontSize',fontSize);

gpatch(Fb_foot,V,'kw','none',faceAlpha2);
% gpatch(Fb_sole,V,'kw','none',faceAlpha2);

hl(1)=gpatch(F_contact_master,V,'gw','k',1);
patchNormPlot(F_contact_master,V);
hl(2)=gpatch(F_contact_slave,V,'bw','k',1);
patchNormPlot(F_contact_slave,V);

legend(hl,{'Master','Slave'});

axisGeom(gca,fontSize);
camlight headlight;
drawnow;

Define boundary conditions

%Supported nodes
bcSupportList=unique(F_sole_bottom);

%Prescribed displacement nodes
bcPrescribeList=unique(Fb_foot(Cb_foot==1,:));

Visualize BC's

hf=cFigure; hold on;
title('Boundary conditions model','FontSize',fontSize);
gpatch(Fb_foot,V,'kw','none',faceAlpha2);

hl2(1)=plotV(V(bcPrescribeList,:),'r.','MarkerSize',markerSize2);
hl2(2)=plotV(V(bcSupportList,:),'k.','MarkerSize',markerSize2);
legend(hl2,{'BC prescribe','BC support'});
axisGeom(gca,fontSize);
camlight headlight;
drawnow;

Visualize rigid body bone surface elements

F_bone=Fb_foot(Cb_foot==1,:);

center_of_mass_bone=mean(V(bcPrescribeList,:),1);

cFigure; hold on;
% title('Boundary conditions model','FontSize',fontSize);
gpatch(Fb_foot,V,'kw','none',faceAlpha2);
gpatch(F_bone,V,'bw','b',1);
axisGeom(gca,fontSize);
camlight headlight;
drawnow;

Set-up materials

switch testCase
    case 0 %Spatially varying based on bone distance
        VE=patchCentre(E_sole,V); %Element centres
        V_vulnerable=VT(bcPrescribeList,:);
        W=minDist(VE,V_vulnerable); %Distance to bone points
        W=W-min(W(:));
        W=W./max(W(:));
        W=W.*(c1_max-c1_min);
        W=W+c1_min;
        c1_desired=W;
    case 1 %
        c1_desired=c1_max*ones(size(E_sole,1),1);
    case 2
        c1_desired=(c1_max+c1_min)/2*ones(size(E_sole,1),1);
    case 3
        c1_desired=c1_min*ones(size(E_sole,1),1);
end

%Snap to nearest available materials
[~,elementMaterialID]=minDist(c1_desired(:),c1_range(:));

Visualize material levels

%Use element2patch to get patch and color data
[F_E_sole,C_F_E_sole]=element2patch(E_sole,elementMaterialID,'penta6');

hf=cFigure; hold on;
title('Material levels','FontSize',fontSize);
gpatch(Fb_foot,V,'w','none',0.1);
for q=1:1:numel(F_E_sole)
    gpatch(F_E_sole{q},V,C_F_E_sole{q},'k',1);
end
axisGeom(gca,fontSize);
camlight headlight;
colormap parula;
icolorbar;
%caxis([0 maxLevelColorbar_SED]);
drawnow;

Convert list of materials to

%Sorting elements according to material label
[elementMaterialID,indSortElements]=sort(elementMaterialID);
E_sole=E_sole(indSortElements,:);

%Removing unused materials
[indUsed,ind1,ind2]=unique(elementMaterialID(:));
elementMaterialID=ind2;
c1=c1_range(indUsed);

numMaterials=numel(indUsed);

Defining the FEBio input structure

See also febioStructTemplate and febioStruct2xml and the FEBio user manual.

%Get a template with default settings
[febio_spec]=febioStructTemplate;

%febio_spec version
febio_spec.ATTR.version='2.5';

%Module section
febio_spec.Module.ATTR.type='solid';

%Control section
febio_spec.Control.analysis.ATTR.type='static';
febio_spec.Control.time_steps=numTimeSteps;
febio_spec.Control.step_size=1/numTimeSteps;
febio_spec.Control.time_stepper.dtmin=dtmin;
febio_spec.Control.time_stepper.dtmax=dtmax;
febio_spec.Control.time_stepper.max_retries=max_retries;
febio_spec.Control.time_stepper.opt_iter=opt_iter;
febio_spec.Control.max_refs=max_refs;
febio_spec.Control.max_ups=max_ups;
febio_spec.Control.symmetric_stiffness=symmetric_stiffness;
febio_spec.Control.min_residual=min_residual;

%Material section
febio_spec.Material.material{1}.ATTR.type='Ogden';
febio_spec.Material.material{1}.ATTR.id=1;
febio_spec.Material.material{1}.c1=c1_1;
febio_spec.Material.material{1}.m1=m1_1;
febio_spec.Material.material{1}.k=k_1;

febio_spec.Material.material{2}.ATTR.type='rigid body';
febio_spec.Material.material{2}.ATTR.id=2;
febio_spec.Material.material{2}.density=1;
febio_spec.Material.material{2}.center_of_mass=center_of_mass_bone;

%Material section
for q=1:1:numMaterials
    febio_spec.Material.material{q+2}.ATTR.type='Ogden unconstrained';
    febio_spec.Material.material{q+2}.ATTR.id=q+2;
    febio_spec.Material.material{q+2}.c1=c1(q);
    febio_spec.Material.material{q+2}.m1=2; %m1=2 to make it Neo-Hookean
    febio_spec.Material.material{q+2}.cp=c1(q)*bulkModulusFactor;
end

%Geometry section
% -> Nodes
febio_spec.Geometry.Nodes{1}.ATTR.name='nodeSet_all'; %The node set name
febio_spec.Geometry.Nodes{1}.node.ATTR.id=(1:size(V,1))'; %The node id's
febio_spec.Geometry.Nodes{1}.node.VAL=V; %The nodel coordinates

% -> Elements
febio_spec.Geometry.Elements{1}.ATTR.type='tet4'; %Element type of this set
febio_spec.Geometry.Elements{1}.ATTR.mat=1; %material index for this set
febio_spec.Geometry.Elements{1}.ATTR.name='Foot'; %Name of the element set
febio_spec.Geometry.Elements{1}.elem.ATTR.id=(1:1:size(E_foot,1))'; %Element id's
febio_spec.Geometry.Elements{1}.elem.VAL=E_foot;

febio_spec.Geometry.Elements{2}.ATTR.type='tri3'; %Element type of this set
febio_spec.Geometry.Elements{2}.ATTR.mat=2; %material index for this set
febio_spec.Geometry.Elements{2}.ATTR.name='Bone'; %Name of the element set
febio_spec.Geometry.Elements{2}.elem.ATTR.id=(size(E_foot,1)+1:1:size(E_foot,1)+size(F_bone,1))'; %Element id's
febio_spec.Geometry.Elements{2}.elem.VAL=F_bone;

n=size(E_foot,1)+size(F_bone,1)+1;
for q=1:1:numMaterials
    logicMaterialNow=(elementMaterialID==q);
    febio_spec.Geometry.Elements{q+2}.ATTR.type='penta6'; %Element type of this set
    febio_spec.Geometry.Elements{q+2}.ATTR.mat=q+2; %material index for this set
    febio_spec.Geometry.Elements{q+2}.ATTR.name=['Sole_mat',num2str(q)]; %Name of the element set
    febio_spec.Geometry.Elements{q+2}.elem.ATTR.id=(n:1:(n-1+nnz(logicMaterialNow)))'; %Element id's
    febio_spec.Geometry.Elements{q+2}.elem.VAL=E_sole(logicMaterialNow,:);
    n=n+nnz(logicMaterialNow);
end

% -> NodeSets
febio_spec.Geometry.NodeSet{1}.ATTR.name='bcSupportList';
febio_spec.Geometry.NodeSet{1}.node.ATTR.id=bcSupportList(:);

febio_spec.Geometry.NodeSet{2}.ATTR.name='bcPrescribeList';
febio_spec.Geometry.NodeSet{2}.node.ATTR.id=bcPrescribeList(:);

% -> Surfaces
febio_spec.Geometry.Surface{1}.ATTR.name='contact_master';
febio_spec.Geometry.Surface{1}.tri3.ATTR.lid=(1:1:size(F_contact_master,1))';
febio_spec.Geometry.Surface{1}.tri3.VAL=F_contact_master;

febio_spec.Geometry.Surface{2}.ATTR.name='contact_slave';
febio_spec.Geometry.Surface{2}.tri3.ATTR.lid=(1:1:size(F_contact_slave,1))';
febio_spec.Geometry.Surface{2}.tri3.VAL=F_contact_slave;

% -> Surface pairs
febio_spec.Geometry.SurfacePair{1}.ATTR.name='Contact1';
febio_spec.Geometry.SurfacePair{1}.master.ATTR.surface=febio_spec.Geometry.Surface{1}.ATTR.name;
febio_spec.Geometry.SurfacePair{1}.slave.ATTR.surface=febio_spec.Geometry.Surface{2}.ATTR.name;

%Boundary condition section
% -> Fix boundary conditions
febio_spec.Boundary.fix{1}.ATTR.bc='x';
febio_spec.Boundary.fix{1}.ATTR.node_set=febio_spec.Geometry.NodeSet{1}.ATTR.name;
febio_spec.Boundary.fix{2}.ATTR.bc='y';
febio_spec.Boundary.fix{2}.ATTR.node_set=febio_spec.Geometry.NodeSet{1}.ATTR.name;
febio_spec.Boundary.fix{3}.ATTR.bc='z';
febio_spec.Boundary.fix{3}.ATTR.node_set=febio_spec.Geometry.NodeSet{1}.ATTR.name;

febio_spec.Boundary.fix{4}.ATTR.bc='x';
febio_spec.Boundary.fix{4}.ATTR.node_set=febio_spec.Geometry.NodeSet{2}.ATTR.name;
febio_spec.Boundary.fix{5}.ATTR.bc='y';
febio_spec.Boundary.fix{5}.ATTR.node_set=febio_spec.Geometry.NodeSet{2}.ATTR.name;

% -> Prescribed boundary conditions on the rigid body
febio_spec.Boundary.rigid_body{1}.ATTR.mat=2;
febio_spec.Boundary.rigid_body{1}.fixed{1}.ATTR.bc='x';
febio_spec.Boundary.rigid_body{1}.fixed{2}.ATTR.bc='y';
febio_spec.Boundary.rigid_body{1}.fixed{3}.ATTR.bc='Rx';
febio_spec.Boundary.rigid_body{1}.fixed{4}.ATTR.bc='Ry';
febio_spec.Boundary.rigid_body{1}.fixed{5}.ATTR.bc='Rz';

febio_spec.Boundary.rigid_body{1}.prescribed.ATTR.bc='z';
febio_spec.Boundary.rigid_body{1}.prescribed.ATTR.lc=1;
febio_spec.Boundary.rigid_body{1}.prescribed.VAL=displacementMagnitude;

% febio_spec.Boundary.rigid_body{1}.force.ATTR.bc='z';
% febio_spec.Boundary.rigid_body{1}.force.ATTR.lc=1;
% febio_spec.Boundary.rigid_body{1}.force.VAL=forceBody;

%Contact section
febio_spec.Contact.contact{1}.ATTR.surface_pair=febio_spec.Geometry.SurfacePair{1}.ATTR.name;
febio_spec.Contact.contact{1}.ATTR.type='sliding-elastic';
febio_spec.Contact.contact{1}.two_pass=1;
febio_spec.Contact.contact{1}.laugon=laugon;
febio_spec.Contact.contact{1}.tolerance=0.2;
febio_spec.Contact.contact{1}.gaptol=0;
febio_spec.Contact.contact{1}.minaug=minaug;
febio_spec.Contact.contact{1}.maxaug=maxaug;
febio_spec.Contact.contact{1}.search_tol=0.01;
febio_spec.Contact.contact{1}.search_radius=0.1;
febio_spec.Contact.contact{1}.symmetric_stiffness=0;
febio_spec.Contact.contact{1}.auto_penalty=1;
febio_spec.Contact.contact{1}.penalty=contactPenalty;
febio_spec.Contact.contact{1}.fric_coeff=fric_coeff;

%Output section
% -> log file
febio_spec.Output.logfile.ATTR.file=febioLogFileName;
febio_spec.Output.logfile.node_data{1}.ATTR.file=febioLogFileName_disp;
febio_spec.Output.logfile.node_data{1}.ATTR.data='ux;uy;uz';
febio_spec.Output.logfile.node_data{1}.ATTR.delim=',';
febio_spec.Output.logfile.node_data{1}.VAL=1:size(V,1);

febio_spec.Output.logfile.rigid_body_data{1}.ATTR.file=febioLogFileName_force;
febio_spec.Output.logfile.rigid_body_data{1}.ATTR.data='Fx;Fy;Fz';
febio_spec.Output.logfile.rigid_body_data{1}.ATTR.delim=',';
febio_spec.Output.logfile.rigid_body_data{1}.VAL=2; %Rigid body material id

febio_spec.Output.logfile.element_data{1}.ATTR.file=febioLogFileName_strainEnergy;
febio_spec.Output.logfile.element_data{1}.ATTR.data='sed';
febio_spec.Output.logfile.element_data{1}.ATTR.delim=',';
febio_spec.Output.logfile.element_data{1}.VAL=1:size(E_foot,1);

Running the FEBio analysis

To run the analysis defined by the created FEBio input file the runMonitorFEBio function is used. The input for this function is a structure defining job settings e.g. the FEBio input file name. The optional output runFlag informs the user if the analysis was run succesfully.

febioAnalysis.run_filename=febioFebFileName; %The input file name
febioAnalysis.run_logname=febioLogFileName; %The name for the log file
febioAnalysis.disp_on=1; %Display information on the command window
febioAnalysis.disp_log_on=1; %Display convergence information in the command window
febioAnalysis.runMode=runMode;%'internal';
febioAnalysis.t_check=0.25; %Time for checking log file (dont set too small)
febioAnalysis.maxtpi=1e99; %Max analysis time
febioAnalysis.maxLogCheckTime=10; %Max log file checking time

if optimizeForceOption==0

Export input file

Exporting the febio_spec structure to an FEBio input file is done using the febioStruct2xml function.

    febioStruct2xml(febio_spec,febioFebFileName); %Exporting to file and domNode

Run febio

    [runFlag]=runMonitorFEBio(febioAnalysis);%START FEBio NOW!!!!!!!!


    if runFlag~=1
        error('FEBio error');
    end
elseif optimizeForceOption==1

    forceDifference=1e6;
    while 1

Set displacement magnitude in input file structure

        febio_spec.Boundary.rigid_body{1}.prescribed.VAL=displacementMagnitude;

Export input file

Exporting the febio_spec structure to an FEBio input file is done using the febioStruct2xml function.

        febioStruct2xml(febio_spec,febioFebFileName); %Exporting to file and domNode

Run febio

        [runFlag]=runMonitorFEBio(febioAnalysis);%START FEBio NOW!!!!!!!!

        if runFlag~=1
            error('FEBio error');
        end
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
--- STARTING FEBIO JOB --- 13-Mar-2020 12:22:50
===========================================================================
         ________    _________   _________     __     _________            
        |        |\ |        |\ |        |\   |  |\  /         \\          
        |    ____|| |    ____|| |    __  ||   |__|| |    ___    ||         
        |   |\___\| |   |\___\| |   |\_| ||    \_\| |   //  \   ||         
        |   ||      |   ||      |   || | ||    __   |  ||    |  ||         
        |   ||__    |   ||__    |   ||_| ||   |  |\ |  ||    |  ||         
        |       |\  |       |\  |         \\  |  || |  ||    |  ||         
        |    ___||  |    ___||  |    ___   || |  || |  ||    |  ||         
        |   |\__\|  |   |\__\|  |   |\__|  || |  || |  ||    |  ||         
        |   ||      |   ||      |   ||  |  || |  || |  ||    |  ||         
        |   ||      |   ||___   |   ||__|  || |  || |   \\__/   ||         
        |   ||      |        |\ |          || |  || |           ||         
        |___||      |________|| |__________|| |__||  \_________//          
                                                                           
      F I N I T E   E L E M E N T S   F O R   B I O M E C H A N I C S      
                                                                           
                 --- v e r s i o n - 2 . 9 . 1 ---                 
                                                                           
                                                                           
  http://febio.org                                                         
  FEBio is a registered trademark.                                         
  copyright (c) 2006-2019 - All rights reserved                            
                                                                           
===========================================================================


Reading file /mnt/data/MATLAB/GIBBON/data/temp/tempModel.feb ...SUCCESS!
 ]0;(0%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 1 : 0.1 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1005190

 1
 Nonlinear solution status: time= 0.1
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787056e+04    1.796322e-22    0.000000e+00 
	   energy              1.239875e+03    8.598305e-12    1.239875e+01 
	   displacement        4.124547e+02    4.124547e+02    4.124547e-04 
 *************************************************************************
 *                               WARNING                                 *
 *                                                                       *
 * No force acting on the system.                                        *
 *************************************************************************

------- converged at time : 0.1

 ]0;(10%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 2 : 0.2 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1005190

 1
 Nonlinear solution status: time= 0.2
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787056e+04    2.097654e+00    0.000000e+00 
	   energy              1.239875e+03    5.605370e-01    1.239875e+01 
	   displacement        4.124547e+02    4.124547e+02    4.124547e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1006558

 2
 Nonlinear solution status: time= 0.2
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787056e+04    2.609836e-02    0.000000e+00 
	   energy              1.239875e+03    1.085745e-02    1.239875e+01 
	   displacement        4.124547e+02    4.816696e-02    4.117548e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1006216

 3
 Nonlinear solution status: time= 0.2
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787056e+04    1.816699e-04    0.000000e+00 
	   energy              1.239875e+03    6.768522e-06    1.239875e+01 
	   displacement        4.124547e+02    2.223353e-03    4.115790e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1006216

 4
 Nonlinear solution status: time= 0.2
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787056e+04    2.087088e-08    0.000000e+00 
	   energy              1.239875e+03    4.469180e-07    1.239875e+01 
	   displacement        4.124547e+02    7.165760e-05    4.115767e-04 

------- converged at time : 0.2

 ]0;(20%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 3 : 0.3 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1006216

 1
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    1.722103e+01    0.000000e+00 
	   energy              1.239561e+03    2.390009e+00    1.239561e+01 
	   displacement        4.096040e+02    4.096040e+02    4.096040e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1009510

 2
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    2.104562e-01    0.000000e+00 
	   energy              1.239561e+03    6.199884e-02    1.239561e+01 
	   displacement        4.096040e+02    1.463780e-01    4.086101e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1008916

 3
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    9.370039e-03    0.000000e+00 
	   energy              1.239561e+03    2.760157e-03    1.239561e+01 
	   displacement        4.096040e+02    4.789680e-02    4.072461e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1008754

 4
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    1.114058e-03    0.000000e+00 
	   energy              1.239561e+03    2.938874e-04    1.239561e+01 
	   displacement        4.096040e+02    2.838250e-03    4.069964e-04 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1008574

 5
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    1.088413e-03    0.000000e+00 
	   energy              1.239561e+03    3.332791e-05    1.239561e+01 
	   displacement        4.096040e+02    6.424351e-04    4.068862e-04 
Reforming stiffness matrix: reformation #6

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1008574

 6
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 7
	stiffness matrix reformations = 6
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    1.459179e-03    0.000000e+00 
	   energy              1.239561e+03    1.019832e-04    1.239561e+01 
	   displacement        4.096040e+02    1.106113e-03    4.068260e-04 
Reforming stiffness matrix: reformation #7

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1008574

 7
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 9
	stiffness matrix reformations = 7
	step from line search         = 0.586716
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    8.930748e-04    0.000000e+00 
	   energy              1.239561e+03    5.631585e-05    1.239561e+01 
	   displacement        4.096040e+02    4.927070e-04    4.068242e-04 
Reforming stiffness matrix: reformation #8

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1008574

 8
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 11
	stiffness matrix reformations = 8
	step from line search         = 0.591175
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    2.720942e-04    0.000000e+00 
	   energy              1.239561e+03    1.446273e-04    1.239561e+01 
	   displacement        4.096040e+02    5.490122e-04    4.068333e-04 
Reforming stiffness matrix: reformation #9

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1008574

 9
 Nonlinear solution status: time= 0.3
	stiffness updates             = 0
	right hand side evaluations   = 13
	stiffness matrix reformations = 9
	step from line search         = 0.500803
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.787348e+04    3.000237e-04    0.000000e+00 
	   energy              1.239561e+03    2.713233e-05    1.239561e+01 
	   displacement        4.096040e+02    3.906708e-05    4.068287e-04 

------- converged at time : 0.3

 ]0;(30%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 4 : 0.4 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1008574

 1
 Nonlinear solution status: time= 0.4
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.789113e+04    2.250555e+01    0.000000e+00 
	   energy              1.239151e+03    3.103029e+00    1.239151e+01 
	   displacement        4.044517e+02    4.044517e+02    4.044517e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1013092

 2
 Nonlinear solution status: time= 0.4
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.789113e+04    2.083536e-01    0.000000e+00 
	   energy              1.239151e+03    9.373163e-02    1.239151e+01 
	   displacement        4.044517e+02    3.606027e-01    4.018699e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1011922

 3
 Nonlinear solution status: time= 0.4
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.789113e+04    4.392346e-03    0.000000e+00 
	   energy              1.239151e+03    1.888481e-03    1.239151e+01 
	   displacement        4.044517e+02    7.273390e-02    4.004870e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1011490

 4
 Nonlinear solution status: time= 0.4
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.789113e+04    2.718162e-04    0.000000e+00 
	   energy              1.239151e+03    4.739171e-05    1.239151e+01 
	   displacement        4.044517e+02    2.268288e-03    4.002549e-04 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1011490

 5
 Nonlinear solution status: time= 0.4
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.789113e+04    1.870046e-04    0.000000e+00 
	   energy              1.239151e+03    1.374579e-05    1.239151e+01 
	   displacement        4.044517e+02    3.830920e-05    4.002507e-04 

------- converged at time : 0.4

 ]0;(40%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 5 : 0.5 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1011490

 1
 Nonlinear solution status: time= 0.5
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.793398e+04    1.365974e+02    0.000000e+00 
	   energy              1.238889e+03    9.492616e+00    1.238889e+01 
	   displacement        3.969408e+02    3.969408e+02    3.969408e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1020904

 2
 Nonlinear solution status: time= 0.5
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.793398e+04    9.601382e-01    0.000000e+00 
	   energy              1.238889e+03    2.350116e-01    1.238889e+01 
	   displacement        3.969408e+02    5.473637e-01    3.944169e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1019230

 3
 Nonlinear solution status: time= 0.5
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.793398e+04    2.801279e-02    0.000000e+00 
	   energy              1.238889e+03    8.768358e-03    1.238889e+01 
	   displacement        3.969408e+02    1.415678e-01    3.924549e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1018222

 4
 Nonlinear solution status: time= 0.5
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.793398e+04    2.549438e-03    0.000000e+00 
	   energy              1.238889e+03    6.214277e-04    1.238889e+01 
	   displacement        3.969408e+02    1.032958e-02    3.919769e-04 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1017862

 5
 Nonlinear solution status: time= 0.5
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.793398e+04    1.342158e-03    0.000000e+00 
	   energy              1.238889e+03    1.274436e-03    1.238889e+01 
	   displacement        3.969408e+02    1.934621e-02    3.916384e-04 
Reforming stiffness matrix: reformation #6

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1017862

 6
 Nonlinear solution status: time= 0.5
	stiffness updates             = 0
	right hand side evaluations   = 8
	stiffness matrix reformations = 6
	step from line search         = 0.603135
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.793398e+04    2.172798e-04    0.000000e+00 
	   energy              1.238889e+03    3.811654e-04    1.238889e+01 
	   displacement        3.969408e+02    8.634988e-03    3.918352e-04 
Reforming stiffness matrix: reformation #7

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1017862

 7
 Nonlinear solution status: time= 0.5
	stiffness updates             = 0
	right hand side evaluations   = 10
	stiffness matrix reformations = 7
	step from line search         = 0.501006
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.793398e+04    1.171206e-04    0.000000e+00 
	   energy              1.238889e+03    4.693792e-05    1.238889e+01 
	   displacement        3.969408e+02    2.438081e-04    3.918270e-04 

------- converged at time : 0.5

 ]0;(50%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 6 : 0.6 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1017862

 1
 Nonlinear solution status: time= 0.6
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.800767e+04    2.763142e+02    0.000000e+00 
	   energy              1.238453e+03    1.871757e+01    1.238453e+01 
	   displacement        3.880372e+02    3.880372e+02    3.880372e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1037662

 2
 Nonlinear solution status: time= 0.6
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.800767e+04    2.260090e+00    0.000000e+00 
	   energy              1.238453e+03    4.880951e-01    1.238453e+01 
	   displacement        3.880372e+02    1.203140e+00    3.852797e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1034872

 3
 Nonlinear solution status: time= 0.6
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.800767e+04    9.906327e-02    0.000000e+00 
	   energy              1.238453e+03    2.228380e-02    1.238453e+01 
	   displacement        3.880372e+02    2.637873e-01    3.826161e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1032712

 4
 Nonlinear solution status: time= 0.6
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.800767e+04    7.451678e-03    0.000000e+00 
	   energy              1.238453e+03    2.131071e-03    1.238453e+01 
	   displacement        3.880372e+02    4.337167e-02    3.820211e-04 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1031794

 5
 Nonlinear solution status: time= 0.6
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.800767e+04    8.264332e-04    0.000000e+00 
	   energy              1.238453e+03    2.517899e-04    1.238453e+01 
	   displacement        3.880372e+02    6.937007e-03    3.819552e-04 
Reforming stiffness matrix: reformation #6

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1031578

 6
 Nonlinear solution status: time= 0.6
	stiffness updates             = 0
	right hand side evaluations   = 7
	stiffness matrix reformations = 6
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.800767e+04    1.420511e-04    0.000000e+00 
	   energy              1.238453e+03    4.034721e-05    1.238453e+01 
	   displacement        3.880372e+02    2.509630e-03    3.819175e-04 
Reforming stiffness matrix: reformation #7

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1031578

 7
 Nonlinear solution status: time= 0.6
	stiffness updates             = 0
	right hand side evaluations   = 8
	stiffness matrix reformations = 7
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.800767e+04    3.787369e-07    0.000000e+00 
	   energy              1.238453e+03    4.867610e-07    1.238453e+01 
	   displacement        3.880372e+02    7.388610e-05    3.819164e-04 

------- converged at time : 0.6

 ]0;(60%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 7 : 0.7 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1031578

 1
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    5.584429e+02    0.000000e+00 
	   energy              1.237178e+03    2.891242e+01    1.237178e+01 
	   displacement        3.737553e+02    3.737553e+02    3.737553e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1053934

 2
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    3.017690e+00    0.000000e+00 
	   energy              1.237178e+03    7.763288e-01    1.237178e+01 
	   displacement        3.737553e+02    2.153642e+00    3.740582e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1050316

 3
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    1.336491e-01    0.000000e+00 
	   energy              1.237178e+03    3.244602e-02    1.237178e+01 
	   displacement        3.737553e+02    3.676329e-01    3.729038e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1048588

 4
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    1.954225e-02    0.000000e+00 
	   energy              1.237178e+03    7.116544e-03    1.237178e+01 
	   displacement        3.737553e+02    2.175623e-01    3.727095e-04 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1047652

 5
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    4.411175e-03    0.000000e+00 
	   energy              1.237178e+03    1.859939e-03    1.237178e+01 
	   displacement        3.737553e+02    7.942120e-02    3.730796e-04 
Reforming stiffness matrix: reformation #6

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1047364

 6
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 7
	stiffness matrix reformations = 6
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    1.848239e-03    0.000000e+00 
	   energy              1.237178e+03    4.077019e-04    1.237178e+01 
	   displacement        3.737553e+02    6.239163e-03    3.732067e-04 
Reforming stiffness matrix: reformation #7

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1047256

 7
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 8
	stiffness matrix reformations = 7
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    9.016213e-04    0.000000e+00 
	   energy              1.237178e+03    1.747055e-04    1.237178e+01 
	   displacement        3.737553e+02    1.548487e-03    3.732262e-04 
Reforming stiffness matrix: reformation #8

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1047238

 8
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 9
	stiffness matrix reformations = 8
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    1.365365e-05    0.000000e+00 
	   energy              1.237178e+03    4.574865e-05    1.237178e+01 
	   displacement        3.737553e+02    1.100289e-03    3.732093e-04 
Reforming stiffness matrix: reformation #9

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1047238

 9
 Nonlinear solution status: time= 0.7
	stiffness updates             = 0
	right hand side evaluations   = 10
	stiffness matrix reformations = 9
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.812714e+04    3.126822e-10    0.000000e+00 
	   energy              1.237178e+03    2.307643e-08    1.237178e+01 
	   displacement        3.737553e+02    4.246053e-06    3.732100e-04 

------- converged at time : 0.7

 ]0;(70%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 8 : 0.8 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1047202

 1
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    7.199969e+02    0.000000e+00 
	   energy              1.236419e+03    3.807546e+01    1.236419e+01 
	   displacement        3.739304e+02    3.739304e+02    3.739304e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1070188

 2
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    2.598768e+00    0.000000e+00 
	   energy              1.236419e+03    6.458493e-01    1.236419e+01 
	   displacement        3.739304e+02    2.484371e+00    3.748844e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1066768

 3
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    1.819910e-01    0.000000e+00 
	   energy              1.236419e+03    3.939338e-02    1.236419e+01 
	   displacement        3.739304e+02    4.209526e-01    3.735262e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1064644

 4
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    6.588946e-02    0.000000e+00 
	   energy              1.236419e+03    1.913890e-02    1.236419e+01 
	   displacement        3.739304e+02    3.241463e-01    3.745058e-04 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1063348

 5
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    3.295739e-02    0.000000e+00 
	   energy              1.236419e+03    4.408418e-03    1.236419e+01 
	   displacement        3.739304e+02    4.416815e-01    3.769015e-04 
Reforming stiffness matrix: reformation #6

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062862

 6
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 7
	stiffness matrix reformations = 6
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    8.928106e-03    0.000000e+00 
	   energy              1.236419e+03    1.236309e-03    1.236419e+01 
	   displacement        3.739304e+02    5.730830e-02    3.776230e-04 
Reforming stiffness matrix: reformation #7

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062646

 7
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 8
	stiffness matrix reformations = 7
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    2.330891e-02    0.000000e+00 
	   energy              1.236419e+03    3.325264e-03    1.236419e+01 
	   displacement        3.739304e+02    9.329633e-03    3.778358e-04 
Reforming stiffness matrix: reformation #8

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062646

 8
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 9
	stiffness matrix reformations = 8
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    1.187611e-02    0.000000e+00 
	   energy              1.236419e+03    2.430642e-03    1.236419e+01 
	   displacement        3.739304e+02    4.247735e-03    3.778612e-04 
Reforming stiffness matrix: reformation #9

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062646

 9
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 10
	stiffness matrix reformations = 9
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    1.646694e-02    0.000000e+00 
	   energy              1.236419e+03    2.948592e-03    1.236419e+01 
	   displacement        3.739304e+02    3.571670e-03    3.778315e-04 
Reforming stiffness matrix: reformation #10

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062646

 10
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 11
	stiffness matrix reformations = 10
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    5.754030e-03    0.000000e+00 
	   energy              1.236419e+03    1.690045e-03    1.236419e+01 
	   displacement        3.739304e+02    3.146832e-03    3.778046e-04 
Reforming stiffness matrix: reformation #11

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062646

 11
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 12
	stiffness matrix reformations = 11
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    6.419186e-03    0.000000e+00 
	   energy              1.236419e+03    1.549591e-03    1.236419e+01 
	   displacement        3.739304e+02    1.796040e-03    3.777765e-04 
Reforming stiffness matrix: reformation #12

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062700

 12
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 13
	stiffness matrix reformations = 12
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    7.633837e-03    0.000000e+00 
	   energy              1.236419e+03    3.678916e-04    1.236419e+01 
	   displacement        3.739304e+02    1.178137e-03    3.777824e-04 
Reforming stiffness matrix: reformation #13

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062700

 13
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 15
	stiffness matrix reformations = 13
	step from line search         = 0.556838
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    1.046965e-02    0.000000e+00 
	   energy              1.236419e+03    9.402932e-04    1.236419e+01 
	   displacement        3.739304e+02    7.069564e-04    3.777679e-04 
Reforming stiffness matrix: reformation #14

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062700

 14
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 16
	stiffness matrix reformations = 14
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    1.451043e-02    0.000000e+00 
	   energy              1.236419e+03    1.211603e-03    1.236419e+01 
	   displacement        3.739304e+02    6.298333e-04    3.777757e-04 
Reforming stiffness matrix: reformation #15

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062700

 15
 Nonlinear solution status: time= 0.8
	stiffness updates             = 0
	right hand side evaluations   = 19
	stiffness matrix reformations = 15
	step from line search         = 0.315941
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            2.831646e+04    7.408438e-03    0.000000e+00 
	   energy              1.236419e+03    8.640738e-04    1.236419e+01 
	   displacement        3.739304e+02    3.148576e-04    3.777645e-04 

------- converged at time : 0.8

 ]0;(80%) tempModel.feb - FEBio 2.9.1.0  
AUTO STEPPER: decreasing time step, dt = 0.0818332


===== beginning time step 9 : 0.881833 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1062628

 1
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    1.061215e+03    0.000000e+00 
	   energy              8.277608e+02    4.028039e+01    8.277608e+00 
	   displacement        2.492718e+02    2.492718e+02    2.492718e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1082302

 2
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    1.682896e+00    0.000000e+00 
	   energy              8.277608e+02    4.443259e-01    8.277608e+00 
	   displacement        2.492718e+02    2.558685e+00    2.516055e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1079098

 3
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    2.224807e-01    0.000000e+00 
	   energy              8.277608e+02    4.157959e-02    8.277608e+00 
	   displacement        2.492718e+02    3.341720e-01    2.526942e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1077370

 4
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    1.076953e-01    0.000000e+00 
	   energy              8.277608e+02    2.957869e-02    8.277608e+00 
	   displacement        2.492718e+02    5.809898e-01    2.554200e-04 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1076578

 5
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    5.069453e-02    0.000000e+00 
	   energy              8.277608e+02    9.503665e-03    8.277608e+00 
	   displacement        2.492718e+02    8.047864e-01    2.608152e-04 
Reforming stiffness matrix: reformation #6

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1075912

 6
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 7
	stiffness matrix reformations = 6
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    1.266298e-02    0.000000e+00 
	   energy              8.277608e+02    2.393130e-03    8.277608e+00 
	   displacement        2.492718e+02    1.135219e-01    2.622602e-04 
Reforming stiffness matrix: reformation #7

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1075678

 7
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 8
	stiffness matrix reformations = 7
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    5.014789e-03    0.000000e+00 
	   energy              8.277608e+02    8.774460e-05    8.277608e+00 
	   displacement        2.492718e+02    1.727181e-02    2.627329e-04 
Reforming stiffness matrix: reformation #8

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1075660

 8
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 9
	stiffness matrix reformations = 8
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    5.720579e-03    0.000000e+00 
	   energy              8.277608e+02    6.785223e-04    8.277608e+00 
	   displacement        2.492718e+02    8.844292e-04    2.628132e-04 
Reforming stiffness matrix: reformation #9

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1075588

 9
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 10
	stiffness matrix reformations = 9
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    2.551387e-03    0.000000e+00 
	   energy              8.277608e+02    4.402403e-04    8.277608e+00 
	   displacement        2.492718e+02    8.772928e-04    2.627603e-04 
Reforming stiffness matrix: reformation #10

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1075588

 10
 Nonlinear solution status: time= 0.881833
	stiffness updates             = 0
	right hand side evaluations   = 12
	stiffness matrix reformations = 10
	step from line search         = 0.559850
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.918121e+04    1.257270e-04    0.000000e+00 
	   energy              8.277608e+02    1.580275e-05    8.277608e+00 
	   displacement        2.492718e+02    1.067426e-04    2.627590e-04 

------- converged at time : 0.881833

 ]0;(88%) tempModel.feb - FEBio 2.9.1.0  
===== beginning time step 10 : 0.963666 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1075588

 1
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    1.205023e+03    0.000000e+00 
	   energy              8.273898e+02    4.099237e+01    8.273898e+00 
	   displacement        2.516370e+02    2.516370e+02    2.516370e-04 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1090222

 2
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    1.435996e+00    0.000000e+00 
	   energy              8.273898e+02    4.127033e-01    8.273898e+00 
	   displacement        2.516370e+02    3.455140e+00    2.565281e-04 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1088440

 3
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    3.532258e-01    0.000000e+00 
	   energy              8.273898e+02    6.184469e-02    8.273898e+00 
	   displacement        2.516370e+02    3.143129e-01    2.586696e-04 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1087738

 4
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    1.760236e-01    0.000000e+00 
	   energy              8.273898e+02    5.391190e-02    8.273898e+00 
	   displacement        2.516370e+02    1.363059e+00    2.653211e-04 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1087216

 5
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    5.275117e-02    0.000000e+00 
	   energy              8.273898e+02    1.670607e-02    8.273898e+00 
	   displacement        2.516370e+02    1.102399e+00    2.733665e-04 
Reforming stiffness matrix: reformation #6

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1086118

 6
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 7
	stiffness matrix reformations = 6
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    1.788545e-02    0.000000e+00 
	   energy              8.273898e+02    3.172355e-03    8.273898e+00 
	   displacement        2.516370e+02    8.939029e-02    2.754797e-04 
Reforming stiffness matrix: reformation #7

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1086064

 7
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 8
	stiffness matrix reformations = 7
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    7.359259e-03    0.000000e+00 
	   energy              8.273898e+02    1.345060e-03    8.273898e+00 
	   displacement        2.516370e+02    1.752152e-02    2.760877e-04 
Reforming stiffness matrix: reformation #8

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1086082

 8
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 9
	stiffness matrix reformations = 8
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    2.952413e-03    0.000000e+00 
	   energy              8.273898e+02    4.525459e-04    8.273898e+00 
	   displacement        2.516370e+02    5.663173e-03    2.764006e-04 
Reforming stiffness matrix: reformation #9

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1086100

 9
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 10
	stiffness matrix reformations = 9
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    1.882748e-03    0.000000e+00 
	   energy              8.273898e+02    2.730614e-04    8.273898e+00 
	   displacement        2.516370e+02    1.404900e-03    2.765395e-04 
Reforming stiffness matrix: reformation #10

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1086082

 10
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 11
	stiffness matrix reformations = 10
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    9.154775e-04    0.000000e+00 
	   energy              8.273898e+02    5.897824e-05    8.273898e+00 
	   displacement        2.516370e+02    5.023820e-04    2.765667e-04 
Reforming stiffness matrix: reformation #11

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1086082

 11
 Nonlinear solution status: time= 0.963666
	stiffness updates             = 0
	right hand side evaluations   = 12
	stiffness matrix reformations = 11
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            1.948911e+04    5.248040e-05    0.000000e+00 
	   energy              8.273898e+02    2.875790e-05    8.273898e+00 
	   displacement        2.516370e+02    4.838424e-05    2.765716e-04 

------- converged at time : 0.963666

 ]0;(96%) tempModel.feb - FEBio 2.9.1.0  
AUTO STEPPER: decreasing time step, dt = 0.0780714

MUST POINT CONTROLLER: adjusting time step. dt = 0.0363337


===== beginning time step 11 : 1 =====
Reforming stiffness matrix: reformation #1

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1086064

 1
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 2
	stiffness matrix reformations = 1
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    1.412823e+01    0.000000e+00 
	   energy              1.636941e+02    9.494760e-02    1.636941e+00 
	   displacement        5.103260e+01    5.103260e+01    5.103260e-05 
Reforming stiffness matrix: reformation #2

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1090492

 2
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 3
	stiffness matrix reformations = 2
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    4.038932e-01    0.000000e+00 
	   energy              1.636941e+02    5.308038e-02    1.636941e+00 
	   displacement        5.103260e+01    5.268280e-02    5.139238e-05 
Reforming stiffness matrix: reformation #3

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1090078

 3
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 4
	stiffness matrix reformations = 3
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    6.941697e-02    0.000000e+00 
	   energy              1.636941e+02    1.530452e-02    1.636941e+00 
	   displacement        5.103260e+01    4.497812e-02    5.180234e-05 
Reforming stiffness matrix: reformation #4

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1089736

 4
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 5
	stiffness matrix reformations = 4
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    4.368358e-02    0.000000e+00 
	   energy              1.636941e+02    1.012663e-02    1.636941e+00 
	   displacement        5.103260e+01    3.411521e-01    5.351097e-05 
Reforming stiffness matrix: reformation #5

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1089646

 5
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 6
	stiffness matrix reformations = 5
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    1.308180e-02    0.000000e+00 
	   energy              1.636941e+02    3.182855e-03    1.636941e+00 
	   displacement        5.103260e+01    1.947570e-01    5.500068e-05 
Reforming stiffness matrix: reformation #6

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1089430

 6
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 7
	stiffness matrix reformations = 6
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    3.498226e-03    0.000000e+00 
	   energy              1.636941e+02    8.378021e-04    1.636941e+00 
	   displacement        5.103260e+01    2.906466e-02    5.556552e-05 
Reforming stiffness matrix: reformation #7

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1089376

 7
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 8
	stiffness matrix reformations = 7
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    9.317263e-04    0.000000e+00 
	   energy              1.636941e+02    3.014558e-04    1.636941e+00 
	   displacement        5.103260e+01    4.511123e-03    5.576635e-05 
Reforming stiffness matrix: reformation #8

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1089376

 8
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 9
	stiffness matrix reformations = 8
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    1.062968e-04    0.000000e+00 
	   energy              1.636941e+02    4.940105e-05    1.636941e+00 
	   displacement        5.103260e+01    1.486396e-03    5.587206e-05 
Reforming stiffness matrix: reformation #9

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1089376

 9
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 10
	stiffness matrix reformations = 9
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    7.139726e-07    0.000000e+00 
	   energy              1.636941e+02    5.424882e-07    1.636941e+00 
	   displacement        5.103260e+01    8.153484e-05    5.589392e-05 
Reforming stiffness matrix: reformation #10

===== reforming stiffness matrix:
	Nr of equations ........................... : 25522
	Nr of nonzeroes in stiffness matrix ....... : 1089376

 10
 Nonlinear solution status: time= 1
	stiffness updates             = 0
	right hand side evaluations   = 11
	stiffness matrix reformations = 10
	step from line search         = 1.000000
	convergence norms :     INITIAL         CURRENT         REQUIRED
	   residual            3.940237e+03    1.690901e-09    0.000000e+00 
	   energy              1.636941e+02    1.268048e-09    1.636941e+00 
	   displacement        5.103260e+01    1.046078e-07    5.589427e-05 

------- converged at time : 1

 ]0;(100%) tempModel.feb - FEBio 2.9.1.0  

N O N L I N E A R   I T E R A T I O N   I N F O R M A T I O N

	Number of time steps completed .................... : 11

	Total number of equilibrium iterations ............ : 88

	Average number of equilibrium iterations .......... : 8

	Total number of right hand evaluations ............ : 108

	Total number of stiffness reformations ............ : 88

	Time in linear solver: 0:00:32


 Elapsed time : 0:01:21


 N O R M A L   T E R M I N A T I O N

Waiting for log file...
Proceeding to check log file...13-Mar-2020 12:24:12
------- converged at time : 0.1
------- converged at time : 0.2
------- converged at time : 0.3
------- converged at time : 0.4
------- converged at time : 0.5
------- converged at time : 0.6
------- converged at time : 0.7
------- converged at time : 0.8
------- converged at time : 0.881833
------- converged at time : 0.963666
------- converged at time : 1
--- Done --- 13-Mar-2020 12:24:13

Importing rigid body reaction forces from a log file

        [time_mat, R_force_mat,~]=importFEBio_logfile(fullfile(savePath,febioLogFileName_force)); %Nodal forces
        time_mat=[0; time_mat(:)]; %Time

        F_reaction=[0 0 0; R_force_mat(:,2:end)]; %Get force, add zeros, remove rigidbody id column

        Fz_final=F_reaction(end,3);

        forceDifference=abs(Fz_final-forceBody);

Use inter/extra-polation to guess displacement

        u=time_mat.*displacementMagnitude;
        f=F_reaction(:,3);
        displacementProposed=interp1(f,u,forceBody,'linear','extrap');
        cFigure; hold on;
        hp(1)=plot(f,u,'r-','LineWidth',lineWidth);
        hp(2)=plot(forceBody,displacementProposed,'b.','MarkerSize',25);
        hl=legend(hp,{'displacement - force FEA','Proposed displacement'},'FontSize',fontSize);
        axis square; axis tight; grid on; box on;
        drawnow;
     displacementDifference=displacementProposed-displacementMagnitude;
     if displacementDifference<-maxProposedDisplacementDifference
         displacementMagnitude=displacementMagnitude-maxProposedDisplacementDifference;
     else
         displacementMagnitude=displacementProposed;
     end
        disp(['Force is: ',num2str(Fz_final),' N. Setting displacement to: ',num2str(displacementMagnitude)]);

        if forceDifference<forceDifferenceTolerance
            break
        end
Force is: -321.048 N. Setting displacement to: -2.6247
    end
end

Visualize force

cFigure; hold on;
hp(1)=plot(time_mat,F_reaction(:,1),'r-','LineWidth',lineWidth);
hp(2)=plot(time_mat,F_reaction(:,2),'g-','LineWidth',lineWidth);
hp(3)=plot(time_mat,F_reaction(:,3),'b-','LineWidth',lineWidth);
hl=legend(hp,{'$F_x$','$F_y$','$F_z$'},'Interpreter','Latex','FontSize',fontSize);
axis square; axis tight; grid on; box on;
drawnow;

Importing nodal displacements from a log file

[time_mat, N_disp_mat,~]=importFEBio_logfile(fullfile(savePath,febioLogFileName_disp)); %Nodal displacements
time_mat=[0; time_mat(:)]; %Time

N_disp_mat=N_disp_mat(:,2:end,:);
sizImport=size(N_disp_mat);
sizImport(3)=sizImport(3)+1;
N_disp_mat_n=zeros(sizImport);
N_disp_mat_n(:,:,2:end)=N_disp_mat;
N_disp_mat=N_disp_mat_n;
DN=N_disp_mat(:,:,end);

V_def=V+DN;
V_DEF=N_disp_mat+repmat(V,[1 1 size(N_disp_mat,3)]);
X_DEF=V_DEF(:,1,:);
Y_DEF=V_DEF(:,2,:);
Z_DEF=V_DEF(:,3,:);

C=sqrt(sum(DN(:,3).^2,2));
[CF]=vertexToFaceMeasure(Fb_foot,C);

Importing element strain energies from a log file

[~,E_energy,~]=importFEBio_logfile(fullfile(savePath,febioLogFileName_strainEnergy)); %Element strain energy

%Remove nodal index column
E_energy=E_energy(:,2:end,:);

%Add initial state i.e. zero energy
sizImport=size(E_energy);
sizImport(3)=sizImport(3)+1;
E_energy_mat_n=zeros(sizImport);
E_energy_mat_n(:,:,2:end)=E_energy;
E_energy=E_energy_mat_n;

[FE_foot,C_energy_foot]=element2patch(E_foot,E_energy(:,:,end));
%     [FE_foot,C_energy_foot]=element2patch(E_foot,E_energy(1:size(E_foot,1),:,1));
indBoundaryFacesFoot=tesBoundary(FE_foot,V);

Plot animation

Plotting the simulated results using anim8 to visualize and animate deformations

% Create basic view and store graphics handle to initiate animation
hf=cFigure; hold on; %Open figure
gtitle([febioFebFileNamePart,': Press play to animate']);
CV=faceToVertexMeasure(FE_foot(indBoundaryFacesFoot,:),V_def,C_energy_foot(indBoundaryFacesFoot));
hp1=gpatch(FE_foot(indBoundaryFacesFoot,:),V_def,CV,'k',1); %Add graphics object to animate
hp1.FaceColor='Interp';
hp2=gpatch(F_E_sole{1},V_def,'kw','none',0.25); %Add graphics object to animate

axisGeom(gca,fontSize);
colormap(gjet(250)); colorbar;
% caxis([0 max(C_energy_foot)/100]);
caxis([0 maxLevelColorbar_SED]);
axis([min(X_DEF(:)) max(X_DEF(:)) min(Y_DEF(:)) max(Y_DEF(:)) min(Z_DEF(:)) max(Z_DEF(:))]);
view([-40 -40]);
camlight headlight;

% Set up animation features
animStruct.Time=time_mat; %The time vector
for qt=1:1:size(N_disp_mat,3) %Loop over time increments
    DN=N_disp_mat(:,:,qt); %Current displacement
    V_def=V+DN; %Current nodal coordinates

    [FE_foot,C_energy_foot]=element2patch(E_foot,E_energy(:,:,qt));
    CV=faceToVertexMeasure(FE_foot(indBoundaryFacesFoot,:),V_def,C_energy_foot(indBoundaryFacesFoot));

    %Set entries in animation structure
    animStruct.Handles{qt}=[hp1 hp1 hp2]; %Handles of objects to animate
    animStruct.Props{qt}={'Vertices','CData','Vertices'}; %Properties of objects to animate
    animStruct.Set{qt}={V_def,CV,V_def}; %Property values for to set in order to animate
end
anim8(hf,animStruct); %Initiate animation feature
drawnow;

GIBBON www.gibboncode.org

Kevin Mattheus Moerman, [email protected]

GIBBON footer text

License: https://github.com/gibbonCode/GIBBON/blob/master/LICENSE

GIBBON: The Geometry and Image-based Bioengineering add-On. A toolbox for image segmentation, image-based modeling, meshing, and finite element analysis.

Copyright (C) 2019 Kevin Mattheus Moerman

This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/.