Iteration time

[retter-berkeley]

I am using py-pde in a non-standard way to try and use it to explore (through RL) controllers for PDE/distributed parameter systems. To allow the controller to act on the system, I am looping over the pde with a short time horizon, then applying the control inputs to the resulting state, and using that state as the basis for another short pde run, as in the following pseudo code:

while error > epsilon:
   use current state to determine control input
   state = state + control input
   solve pd (eq.solve) on updated state, time range = 0.2, adaptive time step, no tracker
   update state with new pde solution
   compare current state with desired state

The problem is the amount of time it takes py-pde to run eq.solve (currently using Diffusion) grows (like 5-10x) as I go through thousands of these iterations. I cannot figure out why. I start with taking about 30 sec to go through 50 iterations and after a couple hundred of these I am up to around 100-120 sec to complete 50 iterations. I need such high iterations for training reinforcement learning. I have tried garbage collect, and deleting variables like grid, state, etc. I know the time . The only way I have found to get back down to 30s iteration time is to restart the kernel in my notebook.

I know py-pde wasn't designed to be used (abused) like this, but is there something I can do to keep iteration time down? Is there some variable I can delete, or a memory leak somewhere? I can't figure out why the iteration time grows, or a more workflow-friendly way to restore it.

[david-zwicker]

Thanks for reporting – this sounds like a memory leak due to caching, potential inside numba. You might be able to confirm this if you watch the memory consumption of the process.

In any case, there might be a better work-around, which likely will also lead to better performance: You could create a stepper function directly:

grid = UnitGrid([32, 32])  # generate grid
state = ScalarField.random_uniform(grid, 0.2, 0.3)  # generate initial condition
eq = pde.DiffusionPDE()
solver = pde.ExplicitSolver(eq, scheme="euler", adaptive=True)  # initialize a solver
stepper = solver.make_stepper(state, dt=1e-3)  # initial time step

The stepper can then be used inside your loop:

while error > epsilon:
   use current state to determine control input
   state = state + control input
   t_final = stepper(state, t_current, t_current+0.2)
   compare current state with desired state

Note that the stepper function only needs to be created once. The stepper updates the state directly, so if you need the state before the stepper, you need to copy it. The stepper returns the final time, which might deviate slightly from the requested end time, particularly if fixed time steps are used.