function [xdot, cm, grf, mforce, mforceCE, mtorqueDOF, stick, outputFEA, Premind] = legs3d(t, x, u, biofile);

% Ordinary Differential Equation for 3D musculoskeletal model : dx/dt = f(x,u,t);
%
% The model has 21 kinematic degrees of freedom and 86 muscles.  The model therefore has 214 state
% variables: 2*21 + 2*86.
%
% Source code:
%	legs3d.sd		Multibody model
%	legs3d*.bio		Muscles and ground contact points
%                   make sure the right grf comes before left grf
%   legs3d*_i       information file for the model, from .sd file
%	legs3d.c		MEX function & ground force-deformation model
%
% Input:
%	t		Times at which the ODE must be evaluated (1 x Nt)
%	x		States (214 x Nt) : q(21), qdot(21), Lce(86), act(86)
%	u		Controls (86 x Nt): stimulation level for muscles
%   biofile    (optional) Name of .bio file and _i file.  File must have 86 muscles.
%   The default one is legs3d.bio and legs3d_i in the same directory if
%   there is no this input argument.
%
% Output
%	xdot    State derivatives (214 x Nt)
%   cm    Center of mass positions and velocities (6 x Nt)
%   grf    right foot GRF vector and torque at COP, then left foot (12xNt)
%   mforce    muscle force (86 x Nt)
%   mforceCE    contract elements of muscle force (86 x Nt)
%   mtorqueDOF    the total muscle torque on each degree of freedom (21 x Nt)
%   stick    global XYZ coordinates of all joints (15) and all points (91*2) of each GRF elements: 3*(15+2*91)
%            (591 x Nt)
%   outputFEA    FEA outputs (force and moment) on right foot, then left foot (12 x Nt)
%
%
% Notes:
% (1) When used with an ODE solver, legs3d is usually evaluated at a single time (Nt=1).
%     The capability to evaluate at multiple times is needed for solving optimal control
%     problems or boundary value problems using collocation.
% (2) The "biofile" parameter can be used to select a model with simple contact (legs3d.bio) or
%     a model with complex contact implemented as a FE model interfaced to Matlab (legs3dFE.bio).
%     When legs3d is called for the first time, the default "legs3d" is used if there is no 
%     this input argument. Subsequent function calls may omit the biofile parameters.
