Merge pull request #97 from ModiaSim/otter_next_main_2022_02_06

Otter next main 2022 02 06

Author: AndreaNeumayr

Manifest.toml

@@ -1,7 +1,7 @@
 name = "Modia3D"
 uuid = "07f2c1e0-90b0-56cf-bda7-b44b56e34eed"
 authors = ["Andrea Neumayr <[email protected]>", "Martin Otter <[email protected]>", "Gerhard Hippmann <[email protected]>"]
-version = "0.8.2"
+version = "0.9.0"
 
 [deps]
 Colors = "5ae59095-9a9b-59fe-a467-6f913c188581"
@@ -32,7 +32,7 @@ FileIO = "1"
 JSON = "0.21"
 Measurements = "2"
 MeshIO = "0.4.10"
-ModiaLang = "0.10.0"
+ModiaLang = "0.11.2"
 MonteCarloMeasurements = "1"
 OrderedCollections = "1"
 Reexport = "1.0"

Project.toml

@@ -48,6 +48,24 @@ Download and install the free DLR SimVis Community Edition, e.g. with [https://v
 
 ## Release Notes
 
+### Version 0.9.0
+
+Non-backwards compatible changes
+
+- Operator `buildModia3D(..)` is removed. Instead, the new constructor `Model3D(..)` must be used at the top level of a
+  Modia3D definition. It is now possible to define several, independent multibody systems 
+  (currently, only one of them can have animation and animation export).
+- If init/start vectors are defined (e.g. initial state of a FreeMotion joint), they must be defined as SVector{3,Float64}(..).
+  Otherwise, errors occur during compilation. 
+
+Other changes
+
+- All test models changed, due to the non-backwards compatible change.
+- Code generation changed significantly, in order that the interface to Modia3D functions is type stable.
+  As a result, simulation is more efficient and much less memory is allocated during simulation.
+- Efficiency improvements for collisions (less memory is allocated during simulation + faster simulation).
+
+
 ### Version 0.8.2
 
 - Fix Cone contact detection

docs/src/index.md

@@ -14,7 +14,7 @@ The simplest technique is to put `@time` in front of the statement that should b
 
 module BouncingSphereSimulation_time
 
-@time begin BouncingSphere = Model(
+@time begin BouncingSphere = Model3D(
     ...
     )
 end
@@ -30,7 +30,7 @@ This gives the following output:
 ```julia
   0.000116 seconds (199 allocations: 12.344 KiB)
 
-Instantiating model Main.BouncingSphereSimulation.buildModia3D(BouncingSphere)
+Instantiating model Main.BouncingSphereSimulation.BouncingSphere
   1.106054 seconds (2.40 M allocations: 150.386 MiB, 3.92% gc time, 62.47% compilation time)
   0.322050 seconds (579.29 k allocations: 32.774 MiB, 5.81% compilation time)
   0.077613 seconds (1.74 k allocations: 97.844 KiB)
@@ -91,7 +91,7 @@ function `getDerivatives!(..)` is called twice to force compilation of this func
 The `logTiming=true` flag generates the following output in this case:
 
 ```
-Instantiating model Main.MobileWithLogTiming.buildModia3D(Mobile)
+Instantiating model Main.MobileWithLogTiming.Mobile
 
 Execute getDerivatives
 First executions of getDerivatives
@@ -101,9 +101,9 @@ First executions of getDerivatives
 ... first simulation:
   8.341341 seconds (37.34 M allocations: 2.236 GiB, 3.36% gc time)
 ... second simulation:
-... Simulate model buildModia3D(Mobile)
+... Simulate model Mobile
       Initialization at time = 0.0 s
-      Termination of buildModia3D(Mobile) at time = 5.0 s
+      Termination of Mobile at time = 5.0 s
         cpuTime         = 11.7 s
         allocated       = 2290.0 MiB
         algorithm       = CVODE_BDF
@@ -123,7 +123,7 @@ First executions of getDerivatives
         nStateEvents    = 0
         nRestartEvents  = 0
 
-... Timings for simulation of buildModia3D(Mobile):
+... Timings for simulation of Mobile:
  ──────────────────────────────────────────────────────────────────────────────────────────────
                                                        Time                   Allocations
                                                ──────────────────────   ───────────────────────

docs/src/internal/Profiling.md

@@ -11,7 +11,7 @@ module BouncingSphere3D
 
 using Modia3D
 
-BouncingSphere = Model(
+BouncingSphere = Model3D(
     boxHeigth = 0.1,
     world      = Object3D(feature=Scene(enableContactDetection = true,
                     animationFile = "BouncingSphere.json")),
@@ -28,7 +28,7 @@ BouncingSphere = Model(
     free       = FreeMotion(obj1=:world, obj2=:sphere, r=Var(init=[0.0, 1.0, 0.0]))
 )
 
-bouncingSphere = @instantiateModel(buildModia3D(BouncingSphere), unitless=true)
+bouncingSphere = @instantiateModel(BouncingSphere, unitless=true)
 simulate!(bouncingSphere, stopTime=2.2, dtmax=0.1)
 
 @usingModiaPlot

docs/src/tutorial/CollisionHandling.md

@@ -15,14 +15,14 @@ module Pendulum1
 
 using Modia3D
 
-Pendulum = Model(
+Pendulum = Model3D(
     world     = Object3D(feature=Scene()),
     body      = Object3D(feature=Solid(massProperties=MassProperties(mass=1.0))),
     bodyFrame = Object3D(parent=:body, translation=[-0.5, 0.0, 0.0]),
     rev       = Revolute(obj1=:world, obj2=:bodyFrame)
 )
 
-pendulum = @instantiateModel(buildModia3D(Pendulum), unitless=true)
+pendulum = @instantiateModel(Pendulum, unitless=true)
 simulate!(pendulum, stopTime=3.0)
 
 @usingModiaPlot     # use the plot package defined by ENV["MODIA_PLOT"]
@@ -38,7 +38,7 @@ julia> include("$(Modia3D.path)/test/Tutorial/Pendulum1.jl")
 
 The `world` [Object3D](@ref) has feature [Scene](@ref) and is therefore the inertial system. The `body` Object3D is a [Solid](@ref) and defines the pendulum as a mass point with `mass = 1.0`. The `bodyFrame` Object3D defines a coordinate system on the `body` that is translated along the x-axis. A revolute joint connects `world` with `bodyFrame`.
 
-With command `buildModia3D(Pendulum)`, the model definition is inspected and a few lines of code included, in order that joint variables are communicated between the Modia equations and the Modia3D multibody program. Keyword `unitless=true` defines that code generation is performed without units (because Modia3D does not yet fully support units in all components).
+With command `Model3D(..)`, a few lines of code are included, in order that joint variables are communicated between the Modia equations and the Modia3D multibody program. Keyword `unitless=true` defines that code generation is performed without units (because Modia3D does not yet fully support units in all components).
 
 The commands above generate an instance of the model, simulate it and generate the following plot:
 
@@ -54,7 +54,7 @@ module Pendulum2
 
 using Modia3D
 
-Pendulum = Model(
+Pendulum = Model3D(
     world = Object3D(feature=Scene(animationFile="Pendulum2.json")),
     obj1  = Object3D(feature=Solid(shape=Beam(axis=1, length=1.0, width=0.2, thickness=0.2),
                 solidMaterial="Steel", visualMaterial=VisualMaterial(color="Blue"))),
@@ -63,7 +63,7 @@ Pendulum = Model(
     rev   = Revolute(obj1=:world, obj2=:obj2)
 )
 
-pendulum = @instantiateModel(buildModia3D(Pendulum), unitless=true)
+pendulum = @instantiateModel(Pendulum, unitless=true)
 simulate!(pendulum, stopTime=3.0)
 
 @usingModiaPlot
@@ -97,7 +97,7 @@ using Modia3D
 
 include("$(Modia3D.modelsPath)/AllModels.jl")
 
-Pendulum = Model(
+Pendulum = Model3D(
     world = Object3D(feature=Scene(animationFile="Pendulum3.json")),
     obj1  = Object3D(feature=Solid(shape=Beam(axis=1, length=1.0, width=0.2, thickness=0.2),
                 solidMaterial="Steel", visualMaterial=VisualMaterial(color="Blue"))),
@@ -111,7 +111,7 @@ Pendulum = Model(
                 (damper.flange_a, fixed.flange)]
 )
 
-pendulum = @instantiateModel(buildModia3D(Pendulum), unitless=true)
+pendulum = @instantiateModel(Pendulum, unitless=true)
 simulate!(pendulum, stopTime=3.0)
 
 @usingModiaPlot

docs/src/tutorial/GettingStarted.md

@@ -119,7 +119,7 @@ Bushing(; kwargs...) = Bushing{Float64}(; kwargs...)
 
 
 # Compute deformation angles from rotation matrix
-function anglesFromRotation(largeAngles::Bool, R12::SMatrix{3,3,F}, w12::SVector{3,F}) where F <: Modia3D.VarFloatType
+function anglesFromRotation(largeAngles::Bool, R12::SMatrix{3,3,F,9}, w12::SVector{3,F})::Tuple{SVector{3,F},SVector{3,F},SMatrix{2,2,F,4}} where F <: Modia3D.VarFloatType
     if largeAngles
         sbe = clamp(R12[3,1], F(-1.0), F(1.0))
         cbe = sqrt(F(1.0) - sbe*sbe)
@@ -132,10 +132,10 @@ function anglesFromRotation(largeAngles::Bool, R12::SMatrix{3,3,F}, w12::SVector
             ald = w12[1] + (sal*w12[2] - cal*w12[3])*sbe/cbe
             bed = cal*w12[2] + sal*w12[3]
             gad = (-sal*w12[2] + cal*w12[3])/cbe
-            return (SVector{3,F}(al, be, ga), SVector{3,F}(ald, bed, gad), SMatrix{2,2,F}(sal, cal, sbe, cbe))
+            return (SVector{3,F}(al, be, ga), SVector{3,F}(ald, bed, gad), SMatrix{2,2,F,4}(sal, cal, sbe, cbe))
         else
             @error("Gimbal lock of Bushing transformation.")
-            return (SVector{3,F}(0.0, 0.0, 0.0), SVector{3,F}(0.0, 0.0, 0.0), SMatrix{2,2,F}(0.0, 0.0, 0.0, 0.0))
+            return (SVector{3,F}(0.0, 0.0, 0.0), SVector{3,F}(0.0, 0.0, 0.0), SMatrix{2,2,F,4}(0.0, 0.0, 0.0, 0.0))
         end
     else
         al = R12[2,3]
@@ -147,12 +147,12 @@ function anglesFromRotation(largeAngles::Bool, R12::SMatrix{3,3,F}, w12::SVector
         if (max(abs(al), abs(be), abs(ga)) > 0.174)
             @warn("Bushing angle exceeds 10 deg.")
         end
-        return (SVector{3,F}(al, be, ga), SVector{3,F}(ald, bed, gad), SMatrix{2,2,F}(0.0, 0.0, 0.0, 0.0))
+        return (SVector{3,F}(al, be, ga), SVector{3,F}(ald, bed, gad), SMatrix{2,2,F,4}(0.0, 0.0, 0.0, 0.0))
     end
 end
 
 # Compute torque vector from force law moments
-function torqueFromMoments(largeAngles::Bool, moments::SVector{3,F}, sico::SMatrix{2,2,F}) where F <: Modia3D.VarFloatType
+function torqueFromMoments(largeAngles::Bool, moments::SVector{3,F}, sico::SMatrix{2,2,F,4})::SVector{3,F} where F <: Modia3D.VarFloatType
     if largeAngles
         tx = moments[1] + sico[1,2]*moments[3]
         ty = sico[2,1]*moments[2] - sico[1,1]*sico[2,2]*moments[3]
@@ -164,13 +164,13 @@ function torqueFromMoments(largeAngles::Bool, moments::SVector{3,F}, sico::SMatr
 end
 
 
-function initializeForceElement(force::Bushing{F}) where F <: Modia3D.VarFloatType
+function initializeForceElement(force::Bushing{F})::Nothing where F <: Modia3D.VarFloatType
     force.obj1.hasForceElement = true
     force.obj2.hasForceElement = true
     return nothing
 end
 
-function evaluateForceElement(force::Bushing{F}) where F <: Modia3D.VarFloatType
+function evaluateForceElement(force::Bushing{F})::Nothing where F <: Modia3D.VarFloatType
     R12 = measFrameRotation(force.obj2; frameOrig=force.obj1)
     r12 = measFramePosition(force.obj2; frameOrig=force.obj1, frameCoord=force.obj1)
     w12 = measFrameRotVelocity(force.obj2; frameOrig=force.obj1, frameCoord=force.obj1)
@@ -190,6 +190,6 @@ function evaluateForceElement(force::Bushing{F}) where F <: Modia3D.VarFloatType
     return nothing
 end
 
-function terminateForceElement(force::Bushing{F}) where F <: Modia3D.VarFloatType
+function terminateForceElement(force::Bushing{F})::Nothing where F <: Modia3D.VarFloatType
     return nothing
 end

src/Composition/ForceElements/Bushing.jl

@@ -20,6 +20,8 @@ Andrea Neumayr and Martin Otter, [DLR - Institute of System Dynamics and Control
 """
 module Composition
 
+export MultibodyData
+
 export initialize, initAnalysis!, performAnalysis!, closeAnalysis!, visualize!, visualizeWorld!
 export updatePosition!, update!
 export RotationVariables, RCardan123
@@ -29,6 +31,13 @@ export FixTranslation
 export setAngle!, connect
 export setDistance!
 
+export openModel3D!
+export setStatesRevolute!  , setAccelerationsRevolute!  , getGenForcesRevolute
+export setStatesPrismatic! , setAccelerationsPrismatic! , getGenForcesPrismatic
+export setStatesFreeMotion!, setAccelerationsFreeMotion!, getGenForcesFreeMotion
+export setStatesFreeMotion_isrot123!
+
+export computeGeneralizedForces!
 
 export distanceAndAngles, distance, planarRotationAngle
 export measFrameRotation, measFramePosition, measFrameDistance
@@ -81,6 +90,7 @@ import Modia3D.Shapes
 import JSON
 import Printf
 import ModiaLang
+import ModiaLang: ModiaBase
 import TimerOutputs
 import MonteCarloMeasurements