Replicating off-the-shelf ODE solvers

I’ve been working through some literature on ODEs and keep coming across the below pseudo-code expression and am having trouble interpreting it for use in the Wolfram Language. I feel like I’ve solved many initial value problems over the years with NDSolve, but now think I need some help from someone who works with ODEs in a more robust capacity to describe the differences.

The above is from here which states that h(t0) is the initial state of the system, G is the function that defines the system. t0 & t1 are the start and end times and Theta is the parameters. (Pytorch also has an implementation I’ve been looking at).

As I researched ways to understand the parameters, I kept coming up with pseudo-code examples using ODESolve. Should I be using a different MMA function? Maybe there is an option I’m not understanding? What’s the most straightforward way to think about this pseudo-code example in Mathematica?