The task of representing thermal systems in generic energy models is not straightforward.
There is a considerable literature on the modeling and optimization of thermal systems themselves — what we might refer to as subsystems. For instance, the optimization of heat-exchanger networks (HEN) is a longstanding research topic. Then there is more applied work, for instance, the modeling and optimization of district heating networks and CHP plant. These models tend to be handcrafted in GAMS and don’t sit well with the design of a generic modeling framework, implemented in an object-oriented language like C++ or python. One thermal model I am a somewhat familiar with is a GAMS program written by Andreas Christidis to optimize the operation of CHP plant with heat accumulators (thermal storage) while selling into a volatile electricity market.
I must refer again to a (now abandoned) modeling framework called deeco, which I maintained for several years. Some of the questions posed in this forum were satisfactorily addressed by deeco, in terms of both overarching design and the implementation of specific technologies. Plant modules (concrete classes) covered a range of thermal processes including heat pumps, seasonal storage, and concentrating solar power (CSP). The CSP module relied on dimensional analysis and Nusselt, Prandtl, and Reynolds numbers to characterize the underlying fluid dynamics. Much of this is work is now archived on GitHub for reference (as Tom noted). Not all of the background thesis work is present but I can probably get hold of missing reports on request.
The key theoretical point is that the intensive and extensive states of the system are orthogonal and that the intensive state (including temperatures) can be set before the extensive state (including plant duties) is determined.
Clearly deeco is not the only or even the best solution on offer, but it does represent a significant attempt at dealing with this design issue. Moreover the model is well documented and mostly well coded. In many respects, the use of storage, be it thermal, hydro, or otherwise, is more challenging, conceptually and numerically, than the representation of thermal subsystems.