Hey,
I am currently working with multibodies and manipulating them during simulation. I noticed that i had different results for the exact same parameters. I manipulate the revolute joints by setting a desired angular velocity (max torque is set). However when i use more complex creatures they start to diverge mainly positioning wise!
I wanted to ask if this behaviour is known and if so if there is some way to manipulate the library to be deterministic?
Best regards
Here is a screenshot:
Hi!
If the two simulations you are running use the exact same initial state, and are running on the same machine, you should indeed get the same results (otherwise it’s considered a bug).
Is there some code I could run to reproduce your problem ?
What exactly do you mean by initial state? I initiate the world by default, building a ground and then add two exact same creatures with the exact same parameters at the same time and let then the world run.
I will try to simplify my code and send the relevant parts to you asap. It is currently a larger project and therefore not too much of use to post it all.
extern crate nalgebra as na;
use na::{Vector2, Isometry2, Point2};
use core::f32::consts::PI;
use ncollide2d::pipeline::CollisionGroups;
use ncollide2d::shape::{Cuboid, Capsule, ShapeHandle};
use nphysics2d::force_generator::DefaultForceGeneratorSet;
use nphysics2d::joint::{DefaultJointConstraintSet, FreeJoint, RevoluteJoint};
use nphysics2d::object::{
BodyPartHandle, ColliderDesc, DefaultBodySet, DefaultColliderSet, Ground, MultibodyDesc,
};
use nphysics2d::world::{DefaultGeometricalWorld, DefaultMechanicalWorld};
use nphysics_testbed2d::Testbed;
static GROUND_GROUP_ID: usize = 0;
static CREATURE_GROUP_ID: usize = 1;
pub fn main() {
let mut testbed = Testbed::new_empty();
init_world(&mut testbed);
testbed.run()
}
pub fn init_world(testbed: &mut Testbed) {
//Generate World
let mechanical_world = DefaultMechanicalWorld::new(Vector2::new(0.0, -10.0)); //Gravity
let geometrical_world = DefaultGeometricalWorld::new();
let mut bodies = DefaultBodySet::new();
let mut colliders = DefaultColliderSet::new();
let joint_constraints = DefaultJointConstraintSet::new();
let force_generators = DefaultForceGeneratorSet::new();
//Add Ground
let ground_group = create_collision_group_ground();
let ground_shape = ShapeHandle::new(Cuboid::new(Vector2::new(1000.0, 1.0)));
let ground_handle = bodies.insert(Ground::new());
let co = ColliderDesc::new(ground_shape)
.translation(-Vector2::y() * 1.0)
.collision_groups(ground_group)
.build(BodyPartHandle(ground_handle, 0));
colliders.insert(co);
//Build the body
let mut shapes: Vec<ColliderDesc<f32>> = Vec::new();
let body = build_body(&mut shapes);
//Build to creatures form the same body
for _ in 0..2 {
//build creature and get handle
let creature = body.build();
let creature_handle = bodies.insert(creature);
//add all colliders for the current creature
for i in 0..shapes.len() {
colliders.insert(shapes[i].build(BodyPartHandle(creature_handle, i)));
}
}
testbed.set_world(
mechanical_world,
geometrical_world,
bodies,
colliders,
joint_constraints,
force_generators,
);
testbed.look_at(Point2::new(0.0, 0.0), 50.0);
}
fn build_body(colliders: &mut Vec<ColliderDesc<f32>>) -> MultibodyDesc<f32> {
//each four consecutive numbers correspond to a line, which will be used for the multibody creation
//he first number is the first x, the second number the first y, the third the second x and the last the second y coordinate
let bones_flat_temp = vec![
48, 216, 120, 208, 120, 208, 112, 328, 112, 328, 88, 352, 112, 328, 144, 344, 120, 208,
168, 152, 168, 152, 112, 88, 168, 152, 224, 144, 224, 144, 256, 72, 224, 144, 304, 216,
304, 216, 304, 320, 304, 320, 272, 344, 304, 320, 336, 344, 304, 216, 416, 144];
//scale the coordinates down
let mut bones_flat = Vec::new();
for i in bones_flat_temp{
bones_flat.push(i as f32/100.0);
}
//create the collider shape for each 4 pairs that correspond to a bone
let mut i = 0;
while i < bones_flat.len()-1 {
let bone_length = (Vector2::new(bones_flat[i+2],bones_flat[i+3]) - Vector2::new(bones_flat[i+0],bones_flat[i+1])).norm();
let shape = ShapeHandle::new(Capsule::new(bone_length/2.0, bone_length/2.0/7.0));//thin them
let collider = ColliderDesc::new(shape.clone())
.collision_groups(create_collision_group_creature(false))
.density(1.0);
colliders.push(collider);
i = i+4;
}
//build the first free body
let x = bones_flat[0];
let y = bones_flat[1];
let x2 = bones_flat[2];
let y2 = bones_flat[3];
let free = FreeJoint::new(Isometry2::new(Vector2::new(0.0, y),(y2-y).atan2(x2-x)+PI/2.0));
let mut body = MultibodyDesc::new(free);
//Add some bones to the body
//bone 1
let parent_vector= Vector2::new(bones_flat[2],bones_flat[3]) - Vector2::new(bones_flat[0],bones_flat[1]);
let child_vector= Vector2::new(bones_flat[6],bones_flat[7]) - Vector2::new(bones_flat[4],bones_flat[5]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x); //get the angle between parent and child
let body_1 = body.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 2
let parent_vector = Vector2::new(bones_flat[6],bones_flat[7]) - Vector2::new(bones_flat[4],bones_flat[5]);
let child_vector= Vector2::new(bones_flat[10],bones_flat[11]) - Vector2::new(bones_flat[8],bones_flat[9]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_2 = body_1.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 3
let child_vector= Vector2::new(bones_flat[14],bones_flat[15]) - Vector2::new(bones_flat[12],bones_flat[13]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_3 = body_1.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 4
let parent_vector= Vector2::new(bones_flat[2],bones_flat[3]) - Vector2::new(bones_flat[0],bones_flat[1]);
let child_vector= Vector2::new(bones_flat[18],bones_flat[19]) - Vector2::new(bones_flat[16],bones_flat[17]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_4 = body.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 5
let parent_vector= Vector2::new(bones_flat[18],bones_flat[19]) - Vector2::new(bones_flat[16],bones_flat[17]);
let child_vector= Vector2::new(bones_flat[22],bones_flat[23]) - Vector2::new(bones_flat[20],bones_flat[21]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_5 = body_4.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 6
let parent_vector= Vector2::new(bones_flat[18],bones_flat[19]) - Vector2::new(bones_flat[16],bones_flat[17]);
let child_vector= Vector2::new(bones_flat[26],bones_flat[27]) - Vector2::new(bones_flat[24],bones_flat[25]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_6 = body_4.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 7
let parent_vector= Vector2::new(bones_flat[26],bones_flat[27]) - Vector2::new(bones_flat[24],bones_flat[25]);
let child_vector = Vector2::new(bones_flat[30],bones_flat[31]) - Vector2::new(bones_flat[28],bones_flat[29]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_7 = body_6.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 8
let parent_vector= Vector2::new(bones_flat[26],bones_flat[27]) - Vector2::new(bones_flat[24],bones_flat[25]);
let child_vector= Vector2::new(bones_flat[34],bones_flat[35]) - Vector2::new(bones_flat[32],bones_flat[33]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_8 = body_6.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 9
let parent_vector= Vector2::new(bones_flat[34],bones_flat[35]) - Vector2::new(bones_flat[32],bones_flat[33]);
let child_vector= Vector2::new(bones_flat[38],bones_flat[39]) - Vector2::new(bones_flat[36],bones_flat[37]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_9 = body_8.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 10
let parent_vector= Vector2::new(bones_flat[38],bones_flat[39]) - Vector2::new(bones_flat[36],bones_flat[37]);
let child_vector= Vector2::new(bones_flat[42],bones_flat[43]) - Vector2::new(bones_flat[40],bones_flat[41]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_10 = body_9.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 11
let parent_vector= Vector2::new(bones_flat[38],bones_flat[39]) - Vector2::new(bones_flat[36],bones_flat[37]);
let child_vector= Vector2::new(bones_flat[46],bones_flat[47]) - Vector2::new(bones_flat[45],bones_flat[46]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_11 = body_9.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
//bone 12
let parent_vector= Vector2::new(bones_flat[34],bones_flat[35]) - Vector2::new(bones_flat[32],bones_flat[33]);
let child_vector= Vector2::new(bones_flat[50],bones_flat[51]) - Vector2::new(bones_flat[48],bones_flat[49]);
let angle = (child_vector.y).atan2(child_vector.x) - (parent_vector.y).atan2(parent_vector.x);
let body_12 = body_8.add_child(get_revo(angle)).set_parent_shift(Vector2::new(0.0,parent_vector.norm()/2.0)).set_body_shift(Vector2::new(0.0,-child_vector.norm()/2.0));
return body;
}
//get a revolutejoint with consistent parameters except the starting angle
fn get_revo(angle: f32) -> RevoluteJoint<f32>{
let mut revo = RevoluteJoint::new(angle);
revo.enable_angular_motor();
revo.enable_max_angle(angle+0.6);
revo.set_desired_angular_motor_velocity(1.0);
revo.enable_min_angle(angle-0.6);
revo.set_max_angular_motor_torque(100.0/16.0);
return revo;
}
fn create_collision_group_creature(collides_with_self: bool) -> CollisionGroups {
let mut creature_group = CollisionGroups::new();
creature_group.set_membership(&[CREATURE_GROUP_ID]);
if !collides_with_self {
creature_group.set_blacklist(&[CREATURE_GROUP_ID]);
}
creature_group.set_whitelist(&[GROUND_GROUP_ID]);
return creature_group;
}
fn create_collision_group_ground() -> CollisionGroups {
let mut ground_group = CollisionGroups::new();
ground_group.set_membership(&[GROUND_GROUP_ID]);
ground_group.set_whitelist(&[CREATURE_GROUP_ID]);
return ground_group;
}
I hope this works, it does for me. (Or not regarding the non consistency this thread is about). If any questions or problems occur please ask. It took for me around 10 seconds before the bodies behave differently.
Thank you for the code to reproduce your issue, I now have a much better understanding of your problem!
What exactly do you mean by initial state?
At first I though you were running the same simulation twice (i.e. you run your binary twice) and get different results. So by “initial state” I meant that the setup of the physics world was exactly the same for both runs.
It turns out you are doing something different : you simulate two exactly identical multibodies inside of the same mechanical and geometrical world. Unfortunately, in this case it is expected that the two multibodies diverge. This is mostly due to numerical errors because of the iterative nature of the constraints solver responsible for computing contact forces. Because we only run a limited number of iterations (by default, 8), the final force computed may be different for both multibodies if convergence was not achieved. The difference originates from the fact that contacts can be solved in a different order for both multibodies, resulting in a different force/friction distribution among the contacts.
One way to reduce significantly the divergence between the two multibodies is to significantly increase the number of velocity solver iterations in order to reach convergence:
mechanical_world.integration_parameters.max_velocity_iterations = 50;
Thank you once again.
If i would instantiate multiple worlds, all using the same parameters and only add one creature to each world, would (or at least should) they all behave the same in their own worlds like the other do in theirs?
Or in another wait, if i loop not over the creatures and adding them, but rather over the code, inititiating each world only with one multibody would it always lead to the same results?
Also just out of interest if you would care to elaborate why exactly these creatures diverge i would be grateful. It is not necessary just something i would like to know. I did not quite get it from your answer above.
Thank you
Yes. With two independent worlds initialized with the same code, you will get the same behavior in both.
Let’s say you are simulating 2D cube lying on the floor. There are two contacts (each on a corner of the cube). Name them { C1, C2 }. Now we have to compute the forces acting at each of those contacts. Clearly the forces applied at each contact is dependent from the force applied at the other contact. To resolve this inter-dependency, we use an iterative approach. Here is an example of iteration:
Here is another example of iteration. Here we choose to deal with the contact C2 first:
This shows that the order chosen to compute forces affect the end-result. In our first scenario computing F2 takes into account F1 while in our second scenario computing F2 does not take the value of F1 into account.
So the order of force computation affects the values of the computed forces (and thus the trajectory of the simulated objects). To prevent this issue, the force solver will repeat this two-steps iteration multiple times, always taking previously computed values into account.
With your robots, the problem is that for two robots in the same physical world, their contacts are not guaranteed to be solved in the same order, hence you are not guaranteed to obtain the same resulting contact forces for both. Increasing the number of solver iteration will reduce this non-determinism.