Lean 4 as a scripting language in Houdini

Stream: lean4

[Tomas Skrivan]

I got really into FFI :) This time I have managed to use Lean 4 as a scripting language in Houdini, a software similar to Blender.

Getting this to work required lots of head scratching and wall banging. I would like to share my experience for anyone who might attempt anything similar.

preview.png

---

Getting Lean as a scripting language in another application(e.g. Houdini) requires creating two dynamic libraries:

1. A plugin for the application, in my case libHouLeanCore.so, that gets loaded at the start up and is responsible for loading/unloading Lean and executing Lean code snippets.
2. A library for the Lean code snippet, in my case Main.so that has to define a function: def run : Sop Unit := ...
   Where Sop is a monad holding application's context information. If the Lean code is supposed to do anything useful it has to somehow interact with the application. Using monads is likely the best way to do it. For Houdini, Sop stands for "Shape OPeration" and it holds a reference to one mutable geometry and multiple references to immutable geometries.

Let's go over those two components in more detail.

---

Application Plugin - libHouLeanCore.so

• The job of the application plugin can be summarized as follows:

  1. Manually load library libleanshared.so from which you need to load function lean_initialize_runtime_module and maybe more if needed.
  Make sure that you load the library with RTLD_GLOBAL, i.e. dlopen(leanshared_path.c_str(), RTLD_GLOBAL | RTLD_NOW). As you need your Lean snippet library to see those symbols.

  2. Manually load library Main.so and load functions initialize_Main, l_run and mk_sop_context. The names may vary, initialize_Main depends on your lean file name and l_run is the C function name corresponding to run : Sop Unit. About mk_sop_context later.
  3. Initialize lean state. I have figured out this by trial and error so like you need to do more. But it seems to be enough to call:
  lean_initialize_runtime_module and initialize_Main

  4. Before calling run : Sop Unit we need to create the Sop context. The Sop monad is defined as

constant SopContext.nonempty : NonemptyType
def SopContext := SopContext.nonempty.type
instance : Nonempty SopContext := SopContext.nonempty.property

inductive SopError where
| missingValue
| user_error (msg : String)

abbrev Sop := EStateM SopError SopContext

where SopContext is an opaque type that holds a c++ class using lean_external_object. The trouble is that to create lean external object you need some functions from libleanshared.so and those need to be loaded manually. To reduce the headache, I have created a function mk_sop_context in the library Main.so that accepts a pointer to the c++ context class and returns a pointer lean_external_object.

5. Call run : Sop Unit, together with creating the context:

auto sopContext = new SopContext;
sopContext->time = time;
sopContext->geo = geo;

l_run(mk_sop_context((void*)sopContext));

6. Before loading new version of Main.so you have to unload Main.so and libleanshared.so

To summarize, the biggest headache and the biggest catch is: the application plugin must not link against libleanshared.so. When you are reloading the Lean snippet library Main.so you have to unload libleanshared.so otherwise you get a crash. Therefore any function from libleanshared.so that you want to use in the plugin has to be loaded manually.

---

Lean Snippet Library - Main.so

• The hardest bit in creating Main.so is to make sure that it has all the symbols it needs. It definitely needs symbols from libleanshared.so but the application plugin is responsible for making sure those are available.

For all other symbols, I took the easy way out and I link everything else statically. There are two main static libraries you want to create

1. libHouLeanCApi.a - This library contains C interface to Houdini. A simple example of a function:

extern "C" lean_object* houlean_time(lean_object* context){
  auto ctx = houlean_lean_to_external<SopContext>(context);
  return houlean_sop_mk_float(ctx->time, context);
}

Which is defined on Lean level as:

@[extern "houlean_time"]
constant time : Sop Float

2. libHouLean.a - This library is created base on lean code and is generated with lake build HouLean:staticLib.
   For example the above function used houlean_sop_mk_float that turns Float into Sop Float. You want Lean to generate such a function automatically, do this with:

@[export houlean_sop_mk_float]
def Sop.capi.float (x : Float) : Sop Float := pure x

The disadvantage of using static libraries is that for example linking against matlib4 the Main.so blows up by 6MB.

For completeness and reference, here are the exact commands I use for turning Main.lean to Main.so

LEAN_PATH=~/Documents/HouLean/build/lib: lean Main.lean -c Main.c -o Main.olean
LEAN_PATH=~/Documents/HouLean/build/lib leanc Main.c -c -o Main.o
leanc Main.o -shared -o Main.so ~/Documents/HouLean/build/lib/libHouLean.a ~/Documents/HouLean/build/cpp/libHouLeanCApi.a

Notice that LEAN_PATH is not set for the last command, I was experiencing crashes when I did provide LEAN_PATH=~/Documents/HouLean/build/lib.

---

Disclaimer: I'm a Houdini developer. If I was more familiar with Blender I would do it for Blender to make it more accessible. But Houdini has free non-commercial license, so you can play around with it anyway.

[Tomas Skrivan]

And of course the point of this exercise if to be able make cool visualizations of simulations made with SciLean
ezgif-5-b6a39d6547.gif

The whole code implementing this mass spring system is:

variable {n : Nat} [NonZero n]

def H (m k : ℝ) (l : ℝ^n) (x p : (ℝ^(3:ℕ))^n) := (1/(2*m)) * ∥p∥² + (k/2 * ∑ i, ∥∥x[i] - x[i+1]∥² - l[i]*l[i]∥²)

def solver (m k : ℝ) (l : ℝ^n) (steps : Nat)
  : Impl (ode_solve (HamiltonianSystem (H m k l))) :=
by
  -- Unfold Hamiltonian definition and compute gradients
  simp[HamiltonianSystem, H]
  autograd

  -- Apply RK4 method
  rw [ode_solve_fixed_dt runge_kutta4_step]
  lift_limit steps "Number of ODE solver steps."; admit; simp

  finish_impl

[Tomas Skrivan]

This is fun to play with :smile: same system different setup ezgif-5-dcc0edfbdf.gif

[Kevin Buzzard]

@Patrick Massot I'm assuming that your sphere eversion proof can't be turned into a concrete construction like this!

[Tomas Skrivan]

I have no clue but I doubt that the proof would be constructive to this level btw. few years back there was a paper presenting an algorithm that can turn bunny(or sphere) inside out https://youtu.be/FIVqa794w3U?t=162 . And the paper's implementation is in Houdini, see project page.

[Tomas Skrivan]

A technical question: If I have a Lean snippets library Main.so and its function run : Sop Unit. Can the main application execute multiple run : Sop Unit functions in parallel? I guess the questions is, are there any global objects that can be mutated by running run : Sop Unit? If I remember correctly, object's reference counter is set to -1 if it is shared across multiple threads. Since the main application would be doing the parallelism, those reference counters would not be correctly set to -1.

My expectation is that the main application can not execute multiple run : Sop Unit in parallel. Any ideas how do it? Maybe wrap run inside of another function that would wrap run inside of IO.Process?

This is a super technical question and probably very few people are familiar with lean internals like this, so maybe @Sebastian Ullrich knows the answer?

[Sebastian Ullrich]

Global objects are always threadsafe. If Sop carries any ("local") objects, you need to make sure to first call lean_mark_mt on them.

[Tomas Skrivan]

Right now, Sop contains only a c++ struct wrapped into lean_external_object and each execution of run gets its own unique version.

Different scenario, if I have two snippet libraries Main1.so and Main2.so then I should be able to safely execute their run function independently without any kind of synchronization, right?

[Tomas Skrivan]

I guess my worry is that they share some variables from libleanshared.so but as you said those should be thread-safe.

[Sebastian Ullrich]

Yes, it shouldn't matter if the functions are from the same shared object or not

[Patrick Massot]

Kevin, the techniques we use in the sphere eversion project can definitely be turned into a fully effective proof. Convex integration is used to draw all pictures from http://hevea-project.fr/ENIndexHevea.html. However we prove such a general theorem that part of the proof isn't fully constructive. Still, we could specialize the proof and make it constructive.

[Tomas Skrivan]

@Patrick Massot Interesting! If at some point you have something computable I'm happy to make it look pretty :smile: (assuming it is in Lean 4 ... I have no clue how the porting from Lean 3 is going along)

[Mario Carneiro]

My understanding is that a fully computable version was considered "too trivial" / uninteresting and was rejected early on. Surely it is possible to just compute a concrete, say, piecewise linear eversion of an icosahedron or perhaps an even smaller triangulated mesh, and prove that no edges are inverted in the process. I would be interested to know how few vertices you can get away with here

[Patrick Massot]

It's not uninteresting, it's something different that wouldn't achieve the goal of proving that you can formalize abstract differential topology. And what I called a constructive version wasn't a PL version it would be a smooth version involving elementary functions like trigonometric functions and integrals that may not have a nice expression. Then plotting the result would involve numerical analysis and meshing as a post-processing step.

[Tomas Skrivan]

Some more fun with Hamiltonian systems:
https://www.youtube.com/watch?v=qcE9hFPgYkg&ab_channel=Lecopivo

Macros in Lean are really cool, I can now annotate function arguments and automatically generate functions derivatives and proofs of smoothness. The Hamiltonian definition for the above system is defined as:

def LennardJones (ε minEnergy : ℝ) (radius : ℝ) (x : ℝ^(3:ℕ)) : ℝ :=
  let x' := ∥1/radius * x∥^{-6, ε}
  4 * minEnergy * x' * (x' - 1)
argument x [Fact (ε≠0)]
  isSmooth, diff, hasAdjDiff, adjDiff

def Coloumb (ε strength mass : ℝ) (x : ℝ^(3:ℕ)) : ℝ := - strength * mass * ∥x∥^{-1,ε}
argument x [Fact (ε≠0)]
  isSmooth, diff, hasAdjDiff, adjDiff

def H (n : ℕ) (ε : ℝ) (C LJ : ℝ) (r m : Idx n → ℝ) (x p : (ℝ^(3:ℕ))^n) : ℝ :=
  (∑ i, (1/(2*m i)) * ∥p[i]∥²)
  +
  ∑ i j,   Coloumb ε C (m i * m j) (x[i] - x[j])
         + LennardJones ε LJ (r i + r j) (x[i] - x[j])
argument p [Fact (n≠0)] [Fact (ε≠0)]
  isSmooth, diff, hasAdjDiff, adjDiff
argument x [Fact (n≠0)] [Fact (ε≠0)]
  isSmooth, diff, hasAdjDiff,
  adjDiff by
    simp[H]; unfold hold; simp
    simp [sum_into_lambda]
    simp [← sum_of_add]

The annotation

argument x [Fact (n≠0)] [Fact (ε≠0)]
  isSmooth, diff, hasAdjDiff,
  adjDiff by
    simp[H]; unfold hold; simp
    simp [sum_into_lambda]
    simp [← sum_of_add]

says under assumption that n≠0 and ε≠0 do:

• isSmooth - prove that function is smooth in x
• diff- generate a new function that is differential of H w.r.t. x and generate new simp theorem
• hasAdjDiff - prove that the differential has adjoint
• adjDif by ... - auto differentiation + some custom transformation to fuse sums together to produce more efficient code

I should play with elaborator that adds these annotations automatically but only if they succeed. However, I'm not sure how to deal with the additional assumptions n≠0 and ε≠0.

[Locria Cyber]

this is epic