3D rendering and animations

Multibody.jl has an automatic 3D-rendering feature that draws a mechanism in 3D. This can be used to create animations of the mechanism's motion from a solution trajectory, as well as to create interactive applications where the evolution of time can be controlled by the user.

The functionality requires the user to install and load one of the Makie backend packages, e.g.,

using GLMakie # Preferred

or

using CairoMakie
Backend choice

GLMakie and WGLMakie produce much nicer-looking animations and are also significantly faster than CairoMakie. CairoMakie may be used to produce the graphics in some web environments if constraints imposed by the web environment do not allow any of the GL alternatives. CairoMakie struggles with the Z-order of drawn objects, sometimes making bodies that should have been visible hidden behind bodies that are further back in the scene.

After that, the render function is the main entry point to create 3D renderings. This function has the following methods:

  • render(model, prob::ODEProblem): this method creates an interactive figure corresponding to the mechanisms configuration at the specified initial condition.
  • render(model, solution): this method creates an animation corresponding to the mechanisms evolution in a simulation trajectory.
  • scene, time = render(model, solution, t::Real): this method opens an interactive window with the mechanism in the configuration corresponding to the time t. Display scene to display the interactive window, and change the time by either dragging the slider in the window, or write to the observable time[] = new_time.

Colors

Many components allows the user to select with which color it is rendered. This choice is made by providing a 4-element array with color values in the order (RGBA), where each value is between 0 and 1. The last value is the alpha channel which determines the opacity, i.e., 1 is opaque and 0 is invisible.

Rendering the world frame

The display of the world frame can be turned off by setting world.render => false in the variable map.

Tracing the path of a frame in 3D visualizations

The path that a frame traces out during simulation can be visualized by passing a vector of frames to the render function using the traces keyword, e.g., render(..., traces=[frame1, frame2]). See the Furuta-pendulum demonstration Going 3D for an example of this.

Camera controls

The camera controls are inherited from Makie, see their documentation for more information. Of particular interest may be the keyboard shortcuts x, y, z, by holding one of these keys and dragging the mouse, the camera will rotate around the corresponding axis. Use keyword argument show_axis = true to function render or pass parameter world.render => true to ODEProblem to display plot axes and/or world axes in the plot.

Rendering API

Multibody.renderFunction
scene       = render(model, prob)
scene, time = render(model, sol, t::Real; framerate = 30, traces = [])
path        = render(model, sol, timevec = range(sol.t[1], sol.t[end], step = 1 / framerate); framerate = 30, timescale=1, display=false, loop=1)

Create a 3D animation of a multibody system

Arguments:

  • model: The unsimplified multibody model, i.e., this is the model before any call to structural_simplify.
  • prob: If an ODEProblem is passed, a static rendering of the system at the initial condition is generated.
  • sol: If an ODESolution produced by simulating the system using solve is passed, an animation or dynamic rendering of the system is generated.
  • t: If a single number t is provided, the mechanism at this time is rendered and a scene is returned together with the time as an Observable. Modify time[] = new_time to change the rendering.
  • timevec: If a vector of times is provided, an animation is created and the path to the file on disk is returned.
  • framerate: Number of frames per second.
  • timescale: Scaling of the time vector. This argument can be made to speed up animations (timescale < 1) or slow them down (timescale > 1). A value of timescale = 2 will be 2x slower than real time.
  • loop: The animation will be looped this many times. Please note: looping the animation using this argument is only recommended when display = true for camera manipulation purposes. When the camera is not manipulated, looping the animation by other means is recommended to avoid an increase in the file size.
  • filename controls the name and the file type of the resulting animation
  • traces: An optional array of frames to show the trace of.
  • show_axis = false: Whether or not to show the plot axes, including background grid.

Camera control

The following keyword arguments are available to control the camera pose:

  • x = 2
  • y = 0.5
  • z = 2
  • lookat = [0,0,0]: a three-vector of coordinates indicating the point at which the camera looks.
  • up = [0,1,0]: A vector indicating the direction that is up.
  • display: if true, the figure will be displayed during the recording process and time will advance in real-time. This allows the user to manipulate the camera options using the mouse during the recording.

See also loop_render

source
Multibody.render!Function
did_render::Bool = render!(scene, ::typeof(ComponentConstructor), sys, sol, t)

Each component that can be rendered must have a render! method. This method is called by render for each component in the system.

This method is responsible for drawing the component onto the scene the way it's supposed to look at time t in the solution sol. t is an Observable. It's recommended to follow the pattern

thing = @lift begin
    acces relevant coordinates from sol at time t
    create a geometric object that can be rendered
end
mesh!(scene, thing; style...)

Returns

A boolean indicating whether or not the component performed any rendering. Typically, all custom methods of this function should return true, while the default fallback method is the only one returning false.

source