Merge pull request #117 from JuliaComputing/nonegw

rm some uses of neg_w=false

docs/src/examples/pendulum.md

@@ -139,7 +139,7 @@ nothing # hide
 ### Why do we need a joint?
 In the example above, we introduced a prismatic joint to model the oscillating motion of the mass-spring system. In reality, we can suspend a mass in a spring without any joint, so why do we need one here? The answer is that we do not, in fact, need the joint, but if we connect the spring directly to the world, we need to make the body (mass) the root object of the kinematic tree instead:
 ```@example pendulum
-@named root_body = Body(; m = 1, isroot = true, r_cm = [0, 1, 0], phi0 = [0, 1, 0], neg_w=false)
+@named root_body = Body(; m = 1, isroot = true, r_cm = [0, 1, 0], phi0 = [0, 1, 0])
 @named multibody_spring = Multibody.Spring(c=10)
 
 connections = [connect(world.frame_b, multibody_spring.frame_a)
@@ -148,11 +148,11 @@ connections = [connect(world.frame_b, multibody_spring.frame_a)
 @named model = ODESystem(connections, t, systems = [world, multibody_spring, root_body])
 ssys = structural_simplify(IRSystem(model))
 
-defs = Dict()
+defs = Dict(collect(root_body.r_0) .=> [0, 1e-3, 0]) # The spring has a singularity at zero length, so we start some distance away
 
 prob = ODEProblem(ssys, defs, (0, 10))
 
-sol = solve(prob, Rodas4(), u0 = prob.u0 .+ 1e-5 .* randn.())
+sol = solve(prob, Rodas4())
 plot(sol, idxs = multibody_spring.r_rel_0[2], title="Mass-spring system without joint")
 ```
 Here, we used a [`Multibody.Spring`](@ref) instead of connecting a `Translational.Spring` to a joint. The `Translational.Spring`, alongside other components from `ModelingToolkitStandardLibrary.Mechanical`, is a 1-dimensional object, whereas multibody components are 3-dimensional objects.

docs/src/examples/swing.md

@@ -102,7 +102,7 @@ Next, we create the full swing assembly
         body_left  = BodyShape(m=0.1, r = [0, 0, w])
         body_right = BodyShape(m=0.1, r = [0, 0, -w])
 
-        body  = Body(m=6, isroot=true, r_cm = [w/2, -w/2, w/2], neg_w=true, cylinder_radius=0.01)
+        body  = Body(m=6, isroot=true, r_cm = [w/2, -w/2, w/2], cylinder_radius=0.01)
 
         damper = Damper(d=0.5)
     end

test/test_quaternions.jl

@@ -1,4 +1,4 @@
-# test utils
+using Test
 import Multibody.Rotations.QuatRotation as Quat
 import Multibody.Rotations
 import Multibody.Rotations: RotXYZ
@@ -149,7 +149,7 @@ end
     t = Multibody.t
     world = Multibody.world
 
-    @named joint = Multibody.FreeMotion(isroot = true, state=true, quat=true, neg_w=false)
+    @named joint = Multibody.FreeMotion(isroot = true, state=true, quat=true, neg_w=true)
     @named body = Body(; m = 1, r_cm = [0.0, 0, 0])
 
     connections = [connect(world.frame_b, joint.frame_a)
@@ -196,15 +196,15 @@ end
 ## Spherical joint pendulum with quaternions
 # ============================================================
 
-@testset "Spherical joint with quaternion state" begin
+# @testset "Spherical joint with quaternion state" begin
     using LinearAlgebra, ModelingToolkit, Multibody, JuliaSimCompiler
     t = Multibody.t
     world = Multibody.world
 
 
     @named joint = Multibody.Spherical(isroot=false, state=false, quat=false)
     @named rod = FixedTranslation(; r = [1, 0, 0])
-    @named body = Body(; m = 1, isroot=true, quat=true, air_resistance=0.0, neg_w=false)
+    @named body = Body(; m = 1, isroot=true, quat=true, air_resistance=0.0, neg_w=true)
 
     # @named joint = Multibody.Spherical(isroot=true, state=true, quat=true)
     # @named body = Body(; m = 1, r_cm = [1.0, 0, 0], isroot=false)
@@ -217,45 +217,42 @@ end
     @named model = ODESystem(connections, t,
                             systems = [world, joint, body, rod])
     irsys = IRSystem(model)
-    @test_skip begin # https://github.com/JuliaComputing/JuliaSimCompiler.jl/issues/298
-        ssys = structural_simplify(irsys)
+    ssys = structural_simplify(irsys)
 
 
-        D = Differential(t)
-        # q0 = randn(4); q0 ./= norm(q0)
-        # q0 = [1,0,0,0]
-        prob = ODEProblem(ssys, [
-            collect(body.w_a) .=> [1,0,0];
-            # collect(body.Q) .=> q0;
-            # collect(body.Q̂) .=> q0;
-            ], (0, 30))
+    D = Differential(t)
+    # q0 = randn(4); q0 ./= norm(q0)
+    # q0 = [1,0,0,0]
+    prob = ODEProblem(ssys, [
+        # collect(body.w_a) .=> [1,0,0];
+        # collect(body.Q) .=> q0;
+        body.k => 0.1;
+        collect(body.Q̂d) .=> [0,0,0,0];
+        ], (0, 30))
 
-        using OrdinaryDiffEq, Test
-        sol = solve(prob, Rodas4(); u0 = prob.u0 .+ 0 .* randn.())
-        @test SciMLBase.successful_retcode(sol)
-        # doplot() && plot(sol, layout=21)
+    using OrdinaryDiffEq, Test
+    sol = solve(prob, Rodas5Pr(); u0 = prob.u0 .+ 0 .* randn.())
+    @test SciMLBase.successful_retcode(sol)
+    # doplot() && plot(sol, layout=21)
 
 
-        ts = 0:0.1:2pi
-        Q = Matrix(sol(ts, idxs = [body.Q...]))
-        Qh = Matrix(sol(ts, idxs = [body.Q̂...]))
-        n = Matrix(sol(ts, idxs = [body.n...]))
-        @test mapslices(norm, Qh, dims=1).^2 ≈ n
-        @test Q ≈ Qh ./ sqrt.(n) atol=1e-2
-        @test norm(mapslices(norm, Q, dims=1) .- 1) < 1e-2
+    ts = 0:0.1:2pi
+    Q = Matrix(sol(ts, idxs = [body.Q...]))
+    Qh = Matrix(sol(ts, idxs = [body.Q̂...]))
+    n = Matrix(sol(ts, idxs = [body.n...]))
+    @test mapslices(norm, Qh, dims=1).^2 ≈ n
+    @test Q ≈ Qh ./ sqrt.(n) atol=1e-2
+    @test norm(mapslices(norm, Q, dims=1) .- 1) < 1e-2
 
-        Matrix(sol(ts, idxs = [body.w_a...]))
+    Matrix(sol(ts, idxs = [body.w_a...]))
 
-        @test get_R(sol, joint.frame_b, 0pi) ≈ I
-        @test get_R(sol, joint.frame_b, sqrt(9.81/1)) ≈ diagm([1, -1, -1]) atol=1e-3
-        @test get_R(sol, joint.frame_b, 2pi) ≈ I atol=1e-3
+    @test get_R(sol, joint.frame_b, 0pi) ≈ I
 
 
-        # Matrix(sol(ts, idxs = [joint.w_rel_b...]))
+    # Matrix(sol(ts, idxs = [joint.w_rel_b...]))
 
-        # render(model, sol)
-    end
-end
+    # render(model, sol)
+# end
 
 
 # ==============================================================================