diff --git a/Scarb.lock b/Scarb.lock index b7727710..3f737994 100644 --- a/Scarb.lock +++ b/Scarb.lock @@ -1,10 +1,15 @@ # Code generated by scarb DO NOT EDIT. version = 1 +[[package]] +name = "alexandria_data_structures" +version = "0.1.0" +source = "git+https://github.com/keep-starknet-strange/alexandria.git?rev=f37d73d#f37d73d8a8248e4d8dc65de3949333e30bda022f" + [[package]] name = "cubit" version = "1.2.0" -source = "git+https://github.com/influenceth/cubit?rev=b459053#b4590530d5aeae9aabd36740cc2a3d9e6adc5fde" +source = "git+https://github.com/influenceth/cubit.git?rev=b459053#b4590530d5aeae9aabd36740cc2a3d9e6adc5fde" [[package]] name = "dojo" @@ -41,7 +46,16 @@ name = "origami" version = "0.0.0" dependencies = [ "cubit", - "dojo", + "orion", +] + +[[package]] +name = "orion" +version = "0.1.7" +source = "git+https://github.com/gizatechxyz/orion?rev=df6d6b9#df6d6b9e0ac3522cac14c28287a35c344840e3b1" +dependencies = [ + "alexandria_data_structures", + "cubit", ] [[package]] diff --git a/Scarb.toml b/Scarb.toml index 51ae95fa..2a1a3987 100644 --- a/Scarb.toml +++ b/Scarb.toml @@ -8,6 +8,7 @@ homepage = "https://github.com/dojoengine/origami" authors = ["bal7hazar@proton.me"] [workspace.dependencies] -cubit = { git = "https://github.com/influenceth/cubit", rev = "b459053" } +cubit = { git = "https://github.com/influenceth/cubit.git", rev = "b459053" } dojo = { git = "https://github.com/dojoengine/dojo", tag = "v0.3.15" } +orion = { git = "https://github.com/gizatechxyz/orion", rev = "df6d6b9" } origami = { path = "crates" } diff --git a/crates/Scarb.toml b/crates/Scarb.toml index 8b2cb66d..6b99eb77 100644 --- a/crates/Scarb.toml +++ b/crates/Scarb.toml @@ -7,4 +7,4 @@ homepage = "https://github.com/dojoengine/origami/tree/main/crates" [dependencies] cubit.workspace = true -dojo.workspace = true +orion.workspace = true diff --git a/crates/src/lib.cairo b/crates/src/lib.cairo index 8c56b3c6..8cd2bfb7 100644 --- a/crates/src/lib.cairo +++ b/crates/src/lib.cairo @@ -4,6 +4,10 @@ mod algebra { mod matrix; } +mod physics { + mod solver; +} + mod defi { mod auction { mod gda; diff --git a/crates/src/physics/force2.cairo b/crates/src/physics/force2.cairo new file mode 100644 index 00000000..40e3cd25 --- /dev/null +++ b/crates/src/physics/force2.cairo @@ -0,0 +1,166 @@ +// External imports + +use orion::numbers::FixedTrait; +// use orion::operators::tensor::{Tensor, TensorTrait, implementation::tensor_fp64x64::FP64x64Tensor}; + +mod errors { + const FORCE2_INVALID_SIZE: felt252 = 'Force2: invalid size'; +} + +#[derive(Drop, Copy)] +struct Force2 { + /// The linear force. + linear: Span, + /// The torque. + angular: N, +} + +trait Force2Trait { + fn new(linear: Span, angular: N) -> Force2; + fn zero() -> Force2; + fn from_span(data: Span) -> Force2; + fn from_spans(linear: Span, angular: Span) -> Force2; + fn torque_from_vector(torque: Span) -> Force2; + fn torque(torque: N) -> Force2; + fn linear(linear: Span) -> Force2; +} + +impl Force2Impl, +Copy, +Drop> of Force2Trait { + /// Creates a force from its linear and angular components. + #[inline(always)] + fn new(linear: Span, angular: N) -> Force2 { + Force2 { linear, angular } + } + + /// A zero force. + #[inline(always)] + fn zero() -> Force2 { + Force2Trait::new(array![].span(), FixedTrait::ZERO()) + } + + /// Create a force from a span where the entries 0 and 1 are for the linear part and 2 for the angular part. + #[inline(always)] + fn from_span(mut data: Span) -> Force2 { + assert(data.len() == 3, errors::FORCE2_INVALID_SIZE); + let torque = data.pop_back().unwrap(); + Force2Trait::new(data, *torque) + } + + /// Creates a force from its linear and angular components, both in vector form. + #[inline(always)] + fn from_spans(linear: Span, angular: Span) -> Force2 { + Force2 { linear, angular: angular.pop_front().unwrap(), } + } + + /// Create a pure torque. + #[inline(always)] + fn torque(torque: N) -> Force2 { + Force2Trait::new(array![].span(), torque) + } + + /// Create a pure torque. + #[inline(always)] + fn torque_from_vector(torque: Span) -> Force2 { + Force2Trait::new(array![].span(), angular.pop_front().unwrap()) + } + + /// Create a pure linear force. + #[inline(always)] + fn linear(linear: Span) -> Force2 { + Force2Trait::new(linear, FixedTrait::ZERO()) + } +/// Creates the resultant of a linear force applied at the given point (relative to the center of mass). +// #[inline(always)] +// fn linear_at_point(linear: Vector2, point: &Point2) -> Force2 { +// Force2Trait::new(linear, point.coords.perp(&linear)) +// } + +// /// Creates the resultant of a torque applied at the given point (relative to the center of mass). +// #[inline(always)] +// fn torque_at_point(torque: N, point: &Point2) -> Force2 { +// Force2Trait::new(point.coords * -torque, torque) +// } + +// /// Creates the resultant of a torque applied at the given point (relative to the center of mass). +// #[inline(always)] +// fn torque_from_vector_at_point(torque: Vector1, point: &Point2) -> Force2 { +// Force2::torque_at_point(torque.x, point) +// } + +// /// The angular part of the force. +// #[inline(always)] +// fn angular_vector(&self) -> Vector1 { +// Vector1::new(self.angular) +// } + +// /// Apply the given transformation to this force. +// #[inline(always)] +// fn transform_by(&self, m: &Isometry2) -> Force2 { +// Force2Trait::new(m * self.linear, self.angular) +// } + +// /// This force seen as a slice. +// /// +// /// The two first entries contain the linear part and the third entry contais the angular part. +// #[inline(always)] +// fn as_slice(&self) -> &[N] { +// self.as_vector().as_slice() +// } + +// /// This force seen as a vector. +// /// +// /// The two first entries contain the linear part and the third entry contais the angular part. +// #[inline(always)] +// fn as_vector(&self) -> &Vector3 { +// unsafe { mem::transmute(self) } +// } + +// /// This force seen as a mutable vector. +// /// +// /// The two first entries contain the linear part and the third entry contais the angular part. +// #[inline(always)] +// fn as_vector_mut(&mut self) -> &mut Vector3 { +// unsafe { mem::transmute(self) } +// } +} + +#[cfg(test)] +mod tests { + // Extenral imports + + use orion::numbers::{FP64x64, FP64x64Impl, FixedTrait}; + + // Local imports + + use super::{Force2, Force2Trait}; + + #[test] + #[available_gas(100_000)] + fn test_force2_new() { + let one: FP64x64 = FixedTrait::new(1, false); + let two: FP64x64 = FixedTrait::new(2, false); + let three: FP64x64 = FixedTrait::new(3, false); + let force: Force2 = Force2Trait::new(array![one, two].span(), three); + assert(force.linear == array![one, two].span(), 'Force2: wrong linear'); + assert(force.angular == three, 'Force2: wrong angular'); + } + + #[test] + #[available_gas(100_000)] + fn test_force2_zero() { + let force: Force2 = Force2Trait::zero(); + assert(force.linear == array![].span(), 'Force2: wrong linear'); + assert(force.angular == FixedTrait::ZERO(), 'Force2: wrong angular'); + } + + #[test] + #[available_gas(100_000)] + fn test_force2_from_slice() { + let one: FP64x64 = FixedTrait::new(1, false); + let two: FP64x64 = FixedTrait::new(2, false); + let three: FP64x64 = FixedTrait::new(3, false); + let force: Force2 = Force2Trait::from_slice(array![one, two, three].span()); + assert(force.linear == array![one, two].span(), 'Force2: wrong linear'); + assert(force.angular == three, 'Force2: wrong angular'); + } +} diff --git a/crates/src/physics/main.cairo b/crates/src/physics/main.cairo new file mode 100644 index 00000000..8b137891 --- /dev/null +++ b/crates/src/physics/main.cairo @@ -0,0 +1 @@ + diff --git a/crates/src/physics/point2.cairo b/crates/src/physics/point2.cairo new file mode 100644 index 00000000..4834bf6e --- /dev/null +++ b/crates/src/physics/point2.cairo @@ -0,0 +1,9 @@ +struct Point2 { + x: N, + y: N, +} + +trait Point2Trait { + #[inline(always)] + fn new(x: N, y: N) -> Point2; +} diff --git a/crates/src/physics/solver.cairo b/crates/src/physics/solver.cairo new file mode 100644 index 00000000..8e3d78a0 --- /dev/null +++ b/crates/src/physics/solver.cairo @@ -0,0 +1,402 @@ +// Source: https://github.com/johnBuffer/VerletSFML/blob/main/solver.hpp + +use dict::{Felt252Dict, Felt252DictTrait}; +use nullable::{NullableTrait, nullable_from_box, match_nullable, FromNullableResult}; + +use orion::numbers::fixed_point::core::{FixedTrait}; +use orion::operators::tensor::core::{Tensor, TensorTrait}; +use orion::numbers::{FP64x64, FP64x64Impl}; +use orion::operators::tensor::implementations::tensor_fp64x64::{ + FP64x64Tensor, FP64x64TensorAdd, FP64x64TensorSub, FP64x64TensorMul, FP64x64TensorDiv +}; + +#[derive(Drop, Copy)] +struct Object { + position: Tensor, + position_last: Tensor, + acceleration: Tensor, + radius: FP64x64, +} + +impl ObjectDefault of Default { + fn default() -> Object { + let zero: FP64x64 = FixedTrait::ZERO(); + let ten: FP64x64 = FixedTrait::new_unscaled(10, false); + Object { + position: TensorTrait::new(array![2].span(), array![zero, zero].span()), + position_last: TensorTrait::new(array![2].span(), array![zero, zero].span()), + acceleration: TensorTrait::new(array![2].span(), array![zero, zero].span()), + radius: ten, + } + } +} + +#[generate_trait] +impl ObjectImpl of ObjectTrait { + fn new(x: FP64x64, y: FP64x64, radius: FP64x64) -> Object { + let zero: FP64x64 = FixedTrait::ZERO(); + let position: Tensor = TensorTrait::new(array![2].span(), array![x, y].span()); + Object { + position: position, + position_last: position, + acceleration: TensorTrait::new(array![2].span(), array![zero, zero].span()), + radius: radius, + } + } + fn update(ref self: Object, dt: FP64x64) { + let dt_squarred: Tensor = TensorTrait::constant_of_shape( + self.acceleration.shape, dt * dt + ); + let displacement = self.position - self.position_last; + self.position_last = self.position; + self.position = self.position + displacement + self.acceleration * dt_squarred; + self + .acceleration = + TensorTrait::constant_of_shape(self.acceleration.shape, FixedTrait::ZERO()); + } + + fn accelerate(ref self: Object, acceleration: Tensor) { + self.acceleration = self.acceleration + acceleration; + } + + fn set_velocity(ref self: Object, velocity: Tensor, dt: FP64x64) { + let dt: Tensor = TensorTrait::constant_of_shape(velocity.shape, dt); + self.position_last = self.position - velocity * dt; + } + + fn add_velocity(ref self: Object, velocity: Tensor, dt: FP64x64) { + let dt: Tensor = TensorTrait::constant_of_shape(velocity.shape, dt); + self.position_last = self.position_last - velocity * dt; + } + + fn get_velocity(ref self: Object, dt: FP64x64) -> Tensor { + let dt: Tensor = TensorTrait::constant_of_shape(self.position.shape, dt); + let velocity = (self.position - self.position_last) / dt; + velocity + } +} + +#[derive(Destruct)] +struct Solver { + sub_steps: u32, + gravity: Tensor, + constraint_center: Tensor, + constraint_radius: FP64x64, + count: u128, + objects: Felt252Dict>, + time: FP64x64, + frame_dt: FP64x64, +} + +impl SolverDefault of Default { + fn default() -> Solver { + let hundred: FP64x64 = FixedTrait::new_unscaled(100, false); + let earth_gravity: FP64x64 = FixedTrait::new_unscaled(981, true) / hundred; + Solver { + sub_steps: 1, + gravity: TensorTrait::new( + array![2].span(), array![FixedTrait::ZERO(), earth_gravity].span() + ), + constraint_center: TensorTrait::new( + array![2].span(), array![FixedTrait::ZERO(), FixedTrait::ZERO()].span() + ), + constraint_radius: hundred, + count: 0, + objects: Default::default(), + time: FixedTrait::ZERO(), + frame_dt: FixedTrait::ZERO(), + } + } +} + +#[generate_trait] +impl PublicImpl of PublicTrait { + fn add_object(ref self: Solver, mut object: Object) { + self.objects.insert(self.count.into(), nullable_from_box(BoxTrait::new(object))); + self.count = self.count + 1; + } + + fn update(ref self: Solver) { + self.time += self.frame_dt; + let dt = self.get_step_dt(); + let mut index = self.sub_steps; + loop { + if index == 0 { + break; + } + self.apply_gravity(); + self.check_collision(dt); + self.apply_constraint(); + self.update_objects(dt); + index -= 1; + }; + } + + fn set_simulation_update_rate(ref self: Solver, rate: u32) { + let one = FixedTrait::new_unscaled(1, false); + let rate = FixedTrait::new_unscaled(rate.into(), false); + self.frame_dt = one / rate; + } + + fn set_constraint(ref self: Solver, position: Tensor, radius: FP64x64) { + self.constraint_center = position; + self.constraint_radius = radius; + } + + fn set_sub_steps_count(ref self: Solver, sub_steps: u32) { + self.sub_steps = sub_steps; + } + + fn set_object_velocity(ref self: Solver, ref object: Object, velocity: Tensor) { + object.set_velocity(velocity, self.get_step_dt()); + } + + fn get_objects(ref self: Solver) -> Span { + let mut objects = array![]; + let mut index: felt252 = 0; + loop { + if index == self.count.into() { + break; + } + let object = self.get_object(index); + objects.append(object); + index += 1; + }; + objects.span() + } + + fn get_constraint(ref self: Solver) -> Span { + let zero_index = array![0].span(); + let one_index = array![1].span(); + array![ + self.constraint_center.at(zero_index), + self.constraint_center.at(one_index), + self.constraint_radius + ] + .span() + } + + fn get_objects_count(ref self: Solver) -> u128 { + self.count.into() + } + + fn get_time(ref self: Solver) -> FP64x64 { + self.time + } + + fn get_step_dt(ref self: Solver) -> FP64x64 { + let sub_steps: FP64x64 = FixedTrait::new_unscaled(self.sub_steps.into(), false); + self.frame_dt / sub_steps + } +} + +#[generate_trait] +impl PrivateImpl of PrivateTrait { + fn get_object(ref self: Solver, index: felt252) -> Object { + let object = match match_nullable(self.objects.get(index)) { + FromNullableResult::Null => Default::default(), + FromNullableResult::NotNull(item) => item.unbox(), + }; + assert(object.radius != FixedTrait::ZERO(), 'Solver: object not found'); + object + } + + fn apply_gravity(ref self: Solver) { + let mut index: felt252 = 0; + loop { + if index == self.count.into() { + break; + } + let mut object = self.get_object(index); + object.accelerate(self.gravity); + self.objects.insert(index, nullable_from_box(BoxTrait::new(object))); + index += 1; + }; + } + + fn check_collision(ref self: Solver, dt: FP64x64) { + // Check at least 2 objects + if self.count < 2 { + return; + } + // Constants + let hundred: FP64x64 = FixedTrait::new_unscaled(100, false); + let half: FP64x64 = FixedTrait::new_unscaled(50, false) / hundred; + let response_coef: FP64x64 = FixedTrait::new_unscaled(75, false) / hundred; + let zero_index = array![0].span(); + // Iterate on all objects + let mut outer_index: felt252 = self.count.into() - 1; + loop { + if outer_index == 0 { + break; + } + let mut object_1 = self.get_object(outer_index); + // Iterate on object involved in new collision pairs + let mut inner_index: felt252 = outer_index - 1; + loop { + let mut object_2 = self.get_object(inner_index); + let vector = (object_1.position - object_2.position); + let distance = vector.matmul(@vector).at(zero_index).sqrt(); + let min_distance = object_1.radius + object_2.radius; + // Check overlap + if distance < min_distance { + let overlap = min_distance - distance; + let distance_vectorized: Tensor = TensorTrait::constant_of_shape( + vector.shape, distance + ); + let normal = vector / distance_vectorized; + let response_1 = object_2.radius / min_distance; + let response_2 = object_1.radius / min_distance; + let delta = half * response_coef * overlap; + let delta_1: Tensor = TensorTrait::constant_of_shape( + vector.shape, delta * response_1 + ); + let delta_2: Tensor = TensorTrait::constant_of_shape( + vector.shape, delta * response_2 + ); + // Update positions + object_1.position = object_1.position + normal * delta_1; + object_2.position = object_2.position - normal * delta_2; + // Update objects + self.objects.insert(outer_index, nullable_from_box(BoxTrait::new(object_1))); + self.objects.insert(inner_index, nullable_from_box(BoxTrait::new(object_2))); + }; + if inner_index == 0 { + break; + } + inner_index -= 1; + }; + outer_index -= 1; + }; + } + + fn apply_constraint(ref self: Solver) { + let zero_index = array![0].span(); + let one_index = array![0].span(); + let mut index: felt252 = 0; + loop { + if index == self.count.into() { + break; + } + let mut object = self.get_object(index); + let vector = self.constraint_center - object.position; + let x = vector.at(zero_index); + let y = vector.at(one_index); + let distance = vector.matmul(@vector).at(zero_index).sqrt(); + let min_distance = self.constraint_radius - object.radius; + // Check overlap + if distance > min_distance { + let distance_vectorized: Tensor = TensorTrait::constant_of_shape( + vector.shape, distance + ); + let normal = vector / distance_vectorized; + let delta_vectorized: Tensor = TensorTrait::constant_of_shape( + vector.shape, min_distance + ); + // Update position + object.position = self.constraint_center - normal * delta_vectorized; + // Update object + self.objects.insert(index, nullable_from_box(BoxTrait::new(object))); + }; + index += 1; + }; + } + + fn update_objects(ref self: Solver, dt: FP64x64) { + let mut index: felt252 = 0; + loop { + if index == self.count.into() { + break; + } + let mut object = self.get_object(index); + object.update(dt); + self.objects.insert(index, nullable_from_box(BoxTrait::new(object))); + index += 1; + }; + } +} + +#[cfg(test)] +mod tests { + // Core imports + + use debug::PrintTrait; + + // External imports + + use orion::numbers::fixed_point::core::{FixedTrait}; + use orion::numbers::{NumberTrait, FP64x64, FP64x64Impl}; + use orion::operators::tensor::core::{Tensor, TensorTrait}; + use orion::operators::tensor::implementations::tensor_fp64x64::{FP64x64Tensor}; + + // Local imports + + use super::{Object, ObjectTrait, Solver, PublicTrait, PrivateTrait}; + + #[test] + #[available_gas(1_000_000)] + fn test_create_default_object() { + let object: Object = Default::default(); + } + + #[test] + #[available_gas(1_000_000)] + fn test_create_new_object() { + let x: FP64x64 = FixedTrait::new_unscaled(1, false); + let y: FP64x64 = FixedTrait::new_unscaled(2, false); + let r: FP64x64 = FixedTrait::new_unscaled(3, false); + let object: Object = ObjectTrait::new(x, y, r); + } + + #[test] + #[available_gas(1_000_000)] + fn test_solver_add_object() { + let x: FP64x64 = FixedTrait::new_unscaled(1, false); + let y: FP64x64 = FixedTrait::new_unscaled(2, false); + let r: FP64x64 = FixedTrait::new_unscaled(3, false); + let object: Object = ObjectTrait::new(x, y, r); + let mut solver: Solver = Default::default(); + solver.add_object(object); + } + + #[test] + #[available_gas(1_000_000_000)] + fn test_solver_one_object_update() { + let x: FP64x64 = FixedTrait::new_unscaled(1, false); + let y: FP64x64 = FixedTrait::new_unscaled(2, false); + let r: FP64x64 = FixedTrait::new_unscaled(3, false); + let object: Object = ObjectTrait::new(x, y, r); + let mut solver: Solver = Default::default(); + solver.add_object(object); + solver.update(); + } + + #[test] + #[available_gas(1_000_000_000)] + fn test_solver_two_objects_update() { + let mut solver: Solver = Default::default(); + let x: FP64x64 = FixedTrait::new_unscaled(1, false); + let y: FP64x64 = FixedTrait::new_unscaled(2, false); + let r1: FP64x64 = FixedTrait::new_unscaled(3, false); + let r2: FP64x64 = FixedTrait::new_unscaled(4, false); + let object_1: Object = ObjectTrait::new(x, y, r1); + solver.add_object(object_1); + let object_2: Object = ObjectTrait::new(y, x, r2); + solver.add_object(object_2); + solver.update(); + let objects = solver.get_objects(); + let zero: FP64x64 = FixedTrait::ZERO(); + let zero_index = array![0].span(); + let one_index = array![1].span(); + let updated_object_1 = *objects.at(0); + assert(updated_object_1.radius == object_1.radius, 'Solver: wrong object 1 radius'); + assert( + updated_object_1.position.at(zero_index) != object_1.position.at(zero_index), + 'Solver: wrong object Px' + ); + assert( + updated_object_1.position.at(one_index) != object_1.position.at(one_index), + 'Solver: wrong object Py' + ); + } +} diff --git a/crates/src/physics/vector2.cairo b/crates/src/physics/vector2.cairo new file mode 100644 index 00000000..9f4e48c3 --- /dev/null +++ b/crates/src/physics/vector2.cairo @@ -0,0 +1,14 @@ +// External imports + +use orion::numbers::FixedTrait; +use orion::operators::tensor::{TensorTrait, Tensor}; + +struct Vector2 { + x: N, + y: N, +} + +trait Vector2Trait { + #[inline(always)] + fn new(x: N, y: N) -> Vector2; +}