Added comments in test_contactforces.py for the supposed comparison.
Fixed the imported functions from bicycle.py.
Built the module of simulation.py

Simulation of bicycle model in reference and steady turning equilibrium
condition.
Calculate the contact force with the states values from the simulation.
Simulation and plotting the states value varying over time.

Arrange the matrix of dynamic, f_a and kindiffs for mass and forcing matrix.
Built function, constructed the initial conditions and time vector for
simulation of odeint and final plotting.
Simulated the dynamic performance of the bicycle model in reference
configuration. The model was imported from BicycleModel in model.py.

However, I found the forcing_full matrix has uas and uads of auxiliary
speeds in it, which mean a problem of substitution may exist in KaneMethod
of mecha, although the forcing matrix is right after replacing uas and
uads with zeros. Overall, the forcing matrix issue should go to Gilbert for
an answer.
Added ua_zero for substitution into the points accelerations in order to
eliminate the ua and uads in forcing matrix. However, it did not work; I
should go to Gilbert for an answer. Simultaneously, in order to simulate
the motion of biycle model, I introduced the ua_zero and uad_zero in
simulation.py.
Renamed the initial values in reference configuration with "_ref",
in order to start new initial values in a new configuration.
Fix q2, q4 with q2_ref and q4_ref in pitch_from_roll_and_steer.
Built a function of plotting for plotting states values.
Updated the way to check the equilibrium values of steady turning.

Starting from u5_value(yaw rate), we check the sqare of it to be positive
value which means the steady turning is feasible. In fact, the stability
of the steady turning still needs to be checked later on.