Generally speaking, the timescales targeted by numerical ice sheet model developers have decreased over time, from the Milankovic-forcing scales of early conceptual models, to the century-scales targeted by sophisticated modern codes simulating the evolution of ice streams. Most of the models used today to make large-scale projections of the Earth's ice sheets were built in the 2000s and 2010s using the state-of-the-art computational fluid dynamics methods of the time. Many of these, such as the BISICLES ice sheet model, were aimed at covering timescales from centuries to multiple millennia. As such, to facilitate these longer timescales, they were built around technology that guaranteed high performance, such as libraries enabling adaptive mesh refinement. To ensure accuracy on the shorter timescales, many also include sophisticated methods of model initialisation, such as solving inverse problems using gradient-based optimisation. Since then, the proliferation of satellite data has increased demand on models such as BISICLES for diagnostic modelling and projections on the sub-century timescale. However, BISICLES was not built for this task and is not so easily adapted to it. Thankfully, there have been advances in numerical methods for solving CFD problems in the time BISICLES has been around, and we can hope to take advantage of these to build something that is. In this work, we are building a finite volume code, using the numerical computation library JAX with the aim of producing something similar to BISICLES but fully equipped for the satellite era. The focus is on differentiability to facilitate accurate higher-order derivatives, the optimisation of more complex internal models, the replacement of adaptive mesh methods, and taking advantage of GPUs. In this talk, we will cover our approach, review the state of progress and discuss the suitability of JAX for this task.