#  ScalSALE: Scalable SALE Benchmark Framework for Supercomputers

Re’em Harel$1,2,3$, Matan Rusanovsky$2,3$, Ron Wagner$4,6$, Harel Levin$2,3$
and Gal Oren$2,5\hskip 2.84544pt$${}^{\textrm{{\char 0\relax}}}$
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
<EMAIL_ADDRESS>$1$Department of Computer Science, Ben-Gurion
University of the Negev, Israel  $2$Scientific Computing Center, Nuclear
Research Center – Negev, Israel  $3$Department of Physics, Nuclear Research
Center – Negev, Israel  $4$Israel Atomic Energy Commission  $5$Department of
Computer Science, Technion – Israel Institute of Technology, Israel  $6$School
of Physics and Astronomy, Tel Aviv University, Israel

###### Abstract

Supercomputers worldwide provide the necessary infrastructure for
groundbreaking research. As the demand for enhancing their performance
constantly increases, new systems are built worldwide to supply this demand at
an increased pace, volume, and power. However, most supercomputers are not
designed equally due to different desired figure of merit, which is derived
from the computational bounds of the targeted scientific applications’
portfolio. In turn, the design of such computers becomes a budget-constrained
optimization process that strives to achieve the best performances possible in
a multi-parameters search space, including heterogeneous hardware, middleware,
and software. Therefore, verifying and evaluating whether a supercomputer can
achieve its desired goal becomes a tedious and complex task. For this purpose,
many full, mini, proxy, and benchmark applications have been introduced in the
attempt to represent scientific applications partially. Nevertheless, as these
benchmarks are hard to expand, update technologically, and most importantly,
are over-simplified compared to scientific applications that tend to couple
multiple scientific domains, they fail to represent the true scaling
capabilities. We suggest a new physical scalable benchmark framework, namely
ScalSALE, based on the well-known SALE scheme. ScalSALE’s main goal is to
provide a simple, flexible, scalable infrastructure that can be easily
expanded to include multi-physical schemes while maintaining scalable and
efficient execution times. By expanding ScalSALE, the gap between the over-
simplified benchmarks and scientific applications can be bridged. To achieve
this goal, ScalSALE is implemented in Modern Fortran with simple object-
oriented design patterns and supported by transparent MPI-3 blocking and non-
blocking communication that allows such a scalable framework. In this manner,
ScalSALE easily achieves almost identical performances in similar settings to
common original benchmark suites of its kind while allowing for expansion and
support for hardware and software modernization. The usage of ScalSALE is
demonstrated with the multi-bounded representative Sedov-Taylor blast wave
problem and compared to the well-known Lulesh benchmark using strong and weak
scaling tests. ScalSALE is executed and evaluated with both rezoning options –
Lagrangian and Eulerian, to demonstrate how the benchmark can be expanded to
solve additional problems while maintaining scalable execution times. The
framework, as well as the results and extensions, are available at:
https://github.com/Scientific-Computing-Lab-NRCN/ScalSALE.

###### Index Terms:

LULESH, Modern Fortran, MPI-3, HPC, Parallel OOP, Scalabale Framework,
Benchmark, oneAPI

##  I Introduction

###  I-A Objectives and Designs of Supercomputers

The combination of high-performance computing (high resolution and fidelity)
and storage (large outputs to be analyzed) in supercomputers provide the
essential infrastructure for groundbreaking research, simulation modeling, and
data analysis in various scientific domains [1].

Most, if not all, scientific applications’ performance tends to be bounded by
a computer component. For example, CPU-bound applications are limited by the
processing speed of the CPUs, thus, they can be sped up by adding more cores
to the machine; I/O-bound applications are limited by the amount of time the
application spends waiting for I/O operations; memory-bound applications are
limited by the memory access speed and the amount of available memory;
latency-bound applications are limited by the time it takes to retrieve data;
bandwidth-bound are limited by the amount of data that needs to be retrieved
(and transferred).

It is common to classify scientific applications into several algorithmic
classes or kernels that can be associated with a specific bound [2]. For
example, dense linear algebra [3] and structured grids with high flop/s rate
are classified as compute-bound; sparse linear algebra and particle methods
[4] with high-performance memory systems are classified as memory-bound; FFTs
[5] with high bisection bandwidth is classified as bandwidth-bound;
unstructured grid or adaptive mesh refinement (AMR) [6] with a low latency
network is classified as latency-bound; and data-intensive computations with
storage and network access can be classified as memory and I/O-bound.
Nevertheless, most applications are associated with multiple algorithmic
classes, thus making them multi-bound.

As the demand for performance enhancement of these algorithmic classes — and
in turn, the associated scientific applications — constantly increases, new
hardware architectures and designs have emerged over the years. For example,
the introduction of Infiniband [7] — a high bandwidth and low-latency network,
improved the performance of scientific applications that use unstructured or
adaptive mesh refinement (AMR), such as in fusion [8], climate [9] and
astrophysics [10]. In the same way, general-purposed GPUs (GPGPUs) [11]
reduced the execution time of dense linear algebra and FFTs calculations
common in material science, chemistry, astrophysics, and deep learning
applications. In addition, the introduction of high bandwidth memory (HBM)
[12] and DDR5 enhanced memory-bound kernels such as particle methods.
Moreover, data-intensive calculations were boosted with the introduction of
burst buffer storage [13]. Thus, it can be claimed that the targeted
applications’ work fashion drives supercomputers’ designs and architectures.

Due to the wide variety of scientific applications and their corresponding
optimal hardware, designing a supercomputer has become much more complex and
diverse compared to previous generations [14] (as can also be seen in the
Top500 list [15]), especially if it is designed under a budget. For example,
the Fugaku [16] supercomputer, which holds the first spot in the HPCG list
[17], contains only CPUs and a low-latency network which makes it optimal for
tightly-coupled simulations and modeling but not for artificial intelligence
applications [18]. On the contrary, Frontier [19], the first exascale
supercomputer, is heterogeneous — containing both CPUs and GPUs. More so,
Facebook’s supercomputer [20] is intended for AI modeling exclusively. Thus,
it contains almost only GPUs. Furthermore, the new HPC as a service (HPCaaS)
allows the formation of large-scale cloud clusters with low network
capabilities (high latency) [21, 22], best for embarrassingly parallel and
loosely coupled computations (with minor to moderate communication operations)
[23]. Such a system is suited for deep learning applications and small-scale
traditional simulations and modeling.

###  I-B Assessment and Evaluation of Supercomputers

Specs Synthetic Benchmarks Kernel Benchmarks Mini/proxy-apps
ApplicationsScalSALEComplexityRepresentative Figure 1: The steps taken in
order to evaluate and assess supercomputers as a function of their complexity
and representation of scientific applications — based on [24]. ScalSALE
extends the mini/proxy apps, such that it will maintain the same complexity as
most mini-apps while overlapping and representing scientific applications.

Over the years, the goals of supercomputers have become diverse — supporting a
wide portfolio of scientific applications — and as such, the design and
architectures are also varied. Thus, the crucial step of assessing,
evaluating, and comparing supercomputers has become hard and complex.

To cope with this hard task, the evaluation is usually broken down into
different steps, associated with its complexity, and representation of the
entire app as a whole — from the most simple test to executing a whole
application (an illustration of these steps can be seen in Figure 1):

  1. 1.

Basic parameters and specifications — include: vendor (AMD, Intel, Cray),
efficiency, Rpeak, architecture (SMP, SIMD), operating system (RHEL, SUSE,
Ubuntu), interconnect family (Ethernet, Infiniband), number of cores,
accelerators family, application area (academic, industry), and more. These
parameters provide a general overview of the machine and its composition. As
this step is not complex at all, it is represented in Figure 1 below the
x-axis.

  2. 2.

Simple Synthetic Benchmarks — performed to test and validate a single-system
component’s basic functionality and behavior. Thus, they are less
representative yet, easy to maintain and highly scalable. For example, IOR
[25] (random I/O operations), ALCF [26] or MPIGraph [27] (basic MPI
operations) and many more. These tests are simple and less representative of
the supercomputer, hence its position in Figure 1.

  3. 3.

Kernel Benchmarks — compute specific dedicated intensive calculations of
actual scientific applications are executed to test the coupling and synergy
between various components. The kernel benchmarks are supposed to be simple
yet represent a small fragment of scientific applications. For example, SPEC
[28] (tests shared-memory parallelization scheme), BTIO [29] (MPI operations
with I/O operations), high-performance Linpack (HPL, computes basic linear
algebra operations) [30], high-performance conjugate gradients (HPCG, performs
the conjugate gradients algorithm) [17] and more. The kernel benchmarks are
slightly complex and more representative, as seen in Figure 1.

  4. 4.

Mini/proxy apps — The kernel benchmarks fail to unveil the supercomputer’s
capabilities in real scientific applications [31]. For example, contrary to
the HPL benchmark (a highly optimized compute-bound kernel), real scientific
applications are usually multi-bound and not necessarily as optimized or as
highly scalable as these benchmarks. This discrepancy only increases in time
due to the rapid advancement in computation power and a relatively slow
advancement in memory bandwidth and latency. As a result, proxy or mini-apps
were devised over the years to bridge the gap between kernel benchmarks and
the necessary evaluation of a supercomputer [32, 33, 34]. These mini-apps
attempt to mimic a portion of the computational workflow of scientific
applications. For example, AMG2013 [35] mimics a multi-grid linear system
solver for unstructured mesh in physical applications; QuickSilver [36] mimics
a portion of the Mercury code [37]; and Lulesh [38] that attempts to mimic a
portion of ALE3D [39].

Nonetheless, most mini-apps still contain intrinsic differences compared to
scientific applications. For example, mini-apps do not necessarily couple
multiple scientific domains, which is common in scientific applications,
especially in their corresponding physical simulations known as multi-physical
simulations. Coupling multiple domains under one application is crucial as it
increases memory consumption — as more data structures are needed to describe
the problem. Hence, the required amount of synchronization (MPI) operations
that depend on the data structures per cycle increases. Furthermore, in some
cases, coupling several domains might heavily change the scheme of the code.
As a result, the scaling capabilities of mini-apps hardly reflect their
counterpart application.

###  I-C Lulesh

The previously mentioned, well-known Lulesh mini-app is highly simplified and
hard-coded to simulate the Sedov-Taylor problem [40] — it solves hydrodynamic
equations (Navier-Stokes) explicitly via an unstructured grid. Hence, it is
latency and bandwidth bound [41]. In addition, Lulesh supports multiple
parallelization schemes and programming paradigms: MPI, OpenMP, CUDA, OpenACC
and OpenCL [38]. Lulesh was developed as part of one of the five challenge
problems defined by the DARPA UHPC program [42, 38] and attempts to mimic a
small portion of various multi-physics applications (specifically, ALE3D [39])
that consume up to 30% of the computing resources of the DoD and DoE [43, 44,
45].

Due to its popularity, Lulesh is commonly used for testing new hardware,
applying optimizations, parallelization schemes, APIs, and more [32, 43]. For
example, Laney et al. [46] assessed the effects of data compression on the
performance of Lulesh and several other benchmarks — as they are the most
representative of real scientific applications. Furthermore, Bercea et al.
[43] evaluated the effectiveness of the OpenMP version on a GPU.

Nevertheless, there are still many crucial distinctions between Lulesh and
other common scientific applications such as the ALE3D (which Lulesh is based
upon). For example, contrary to Lulesh, ALE3D simulates multi-physical
phenomena. Therefore, ALE3D has a higher memory consumption corresponding to
the number of physical quantities necessary to describe multi-physical
phenomena. Moreover, multi-physical simulations also contain additional
calculations for each cycle (corresponding to the additional equations), such
as multi-phase transitions, thermal, advection, chemistry-related
calculations, implicit schemes, and many more. These disparities were
expressed and even quantified in [44]. Lulesh was tuned and optimized on
modern hardware, which resulted in an average of 60% decrease in total
execution time compared to an unoptimized run. However, applying the same
optimizations to ALE3D results in only a 20% time decrease for the same
problem, suggesting that although Lulesh provides a reasonable benchmark for
multi-physics applications, it is far from representing the actual nature of
such applications in supercomputers.

Therefore, similarly to kernel benchmarks, most mini-apps also fail to unveil
the true scaling capabilities of supercomputers regarding common scientific
applications such as the ALE3D, which, as mentioned, represents 30% of
calculations in the DoD and DoE.

###  I-D Contribution — ScalSALE

We propose a scalable mini-app benchmark framework, namely ScalSALE, that will
attempt to bridge the gap between benchmarks and actual scientific
applications — or at the very least, provide a framework that will allow
bridging this gap. The essence of this benchmark is to create a framework or
infrastructure that can be easily expanded to include additional numerical
algorithms, physical models, and alternative parallelization schemes with a
special focus on scalability, i.e., the ability of this benchmark to maintain
a scalable speedup and short execution times, especially compared to currently
available benchmarks.

ScalSALE utilizes Modern Fortran [47] as the programming language. Modern
Fortran provides excellent, if not the best, execution times for numerical
calculations compared to most programming languages while also enabling modern
object-oriented design patterns [48, 49, 50]. Although object-oriented design
patterns provide valuable concepts that help the code be more maintainable,
easy to modify, and expand (which is beneficial for ScalSALE as a framework),
it does not go hand to hand with performance. Thus, ScalSALE employs a
performance object-oriented programming paradigm that aims to utilize design
patterns without harming the performance, i.e., favoring performance in the
trade-off between design patterns and high-performance. Moreover, there are
many legacy codes (that are still widely used) written in Fortran [51, 52, 53]
(hence its popularity among scientists), contrary to the number of benchmarks
in Fortran.

ScalSALE’s infrastructure is composed of two independent layers that
complement one another. The first and most important layer is the kernel. This
layer implements the basic data structures, for example, a 3D array that can
be accessed via a function. Moreover, this layer provides a transparent
interface to the distributed parallelization scheme via transparent parallel
objects implemented with the MPI-3 standard that further supports ScalSALE
scaling capabilities. The layer is independent of the other layer and can be
applied in different codes and projects. The second layer, the application,
contains physical quantities and other related calculations, such as equation-
of-state (EoS), and provides an abstraction layer to the data structures
implemented in the kernel layer. Furthermore, the application utilizes the
physical quantities to implement the physical and numerical models. The
application layer is generic and can be replaced or added by any physical
model, relying on the kernel as an infrastructure.

Currently, ScalSALE implements 1D, 2D and 3D multi-material hydrodynamic
scheme. The hydrodynamic scheme is based on the well-known simplified
arbitrary Lagrangian Euler (SALE) code [54, 55]. The SALE’s simple structure
is easy to follow, maintain, modify and expand. In order to support ScalSALE
as a framework, the numerical calculations and functions in SALE were adapted
to object-oriented design patterns and optimized to modern architectures and
compilers.

(a) 2D Lagrangian Sod simulation.

(b) 3D Eulerian Sedov-Taylor simulation.

Figure 2: ScalSALE simulation of Sod and Sedov-Taylor problems.

In order to test the scaling capabilities, it was compared to the well-known
hydrodynamic mini-app — Lulesh. As mentioned, Lulesh is extensively used,
especially for evaluating supercomputers, which makes it highly suitable for
evaluating ScalSALE. Moreover, both codes perform similar calculations: both
solve the Navier-Stokes equations in 3D; can simulate the Sedov-Taylor problem
(presented in 2(b)); apply the ideal gas equation-of-state and more.

Nevertheless, ScalSALE is compared to Lulesh only in terms of scaling
capabilities, as the main goal of ScalSALE is to provide a scalable framework
that can be expanded easily, especially for next-generation supercomputers in
the exascale era. Contrary to Lulesh, ScalSALE can also solve additional
physical problems such as Sod shock-tube [56] (as seen in 2(a) for Lagrangian
rezoning), spherical divergent shocks (Noh) [57] and more. In addition,
ScalSALE can simulate problems in 1D, 2D, and 3D, with more than one material
(hence, multi-material) and its corresponding equation-of-state. Furthermore,
ScalSALE solves the equations via a structured grid, with Lagrangian or
Eulerian (with an advection phase) rezoning. Thus, providing a wider range of
problems that it can solve. The code-features of Lulesh and ScalSALE are
compared and presented in Table I.

Feature | Lulesh | \+ / | ScalSALE
---|---|---|---
Mesh | Unstructured | / | Structured
Interoperability | – | / | C, C++
Parallelization | Blocking | + | Non-blocking
Paradigm | Variable-oriented | / | OOP
Code Structure |  | Simplified
---
hard-coded
/ | SALE
| Numerical
---
kernels
Implicit | + |  | Lagrangian rezoning
---
Eulerian rezoning
Advection
Equations | Hydrodynamics | + |  | Multi-material
---
Multiple EoS
Dimension | 3D | + | 1D and 2D
Problems | Sedov-Taylor | + | Shock-tube

TABLE I: Comparison between ScalSALE and Lulesh code features. The ”+” refers
to additional features, and ”/” refers to different features.

##  II Foundations of ScalSALE

A modern, scalable scientific benchmark, as a framework, should be readable,
maintainable, extensible, well-tested, well-documented, portable to modern
architectures, and most importantly — be constructed in a way that will allow
it to take advantage of the optimal parallelism schemes [58]. A benchmark that
copes with these requirements must, on the one hand, be written in a
programming language with high-enough abstraction capabilities while, on the
other hand, be as light-weighted as possible. This section presents the
foundations of ScalSALE, that enable such a scalable benchmark framework,
based on [59].

###  II-A High-performance Oriented Design and Code

The programming languages’ abstraction–performance trade-off consists of very
fast languages in one extreme and high abstraction level ones in the other. An
example of the former languages is Assembly, which is very fast since it is
close to the machine code and does not require many translations. However,
programming in Assembly requires a very high effort and expertise and does not
support OOP. Most importantly, it is not suitable as a benchmark since it does
not represent most applications written in higher abstraction level
programming languages. An example of the latter language is Python, which is
very expressive and easy to use and maintain. Furthermore, it supports OOP but
is much slower since it translates the programming abstractions into machine
code in run-time. Between these two extremes, there is a range of languages,
starting from Assembly, with a decreasing performance power coupled with
increasing expressive power. This range includes, among others, Assembly,
Fortran 66, Fortran 77, Fortran 90, C, Modern Fortran, C++, and Python.
Fortran is used in many benchmarks and codes [60, 61, 62, 52, 53] mainly due
to its low execution times while supporting sufficient programming
abstraction.

Fortran underwent a major revision in 2003 and subsequently was re-branded as
Modern Fortran. Modern Fortran supports object-oriented design patterns, thus
enabling system design through modern design concepts [63, 64], while still
benefiting from the fast execution time of the language. In addition, the
language features the basic concepts of object-oriented programming (OOP):
objects, inheritance, polymorphism, dynamic type allocation, type-bound
procedures, and many more. In 2008, Modern Fortran was updated to support
additional sub-modules, coarrays, _DO CONCURRENT_ construct that hints to the
compiler that a loop can be parallelized, and more. The latest version of
Modern Fortran was released in 2018 and included further interoperability with
C and additional parallelization features in Fortran.

Over the years, Fortran and C++ were compared to decide which language is best
suited for large-scale codes, both in execution and in numerical capabilities,
e.g., in multi-physics and scientific applications [49, 65, 48, 66, 67, 68].
The initial conclusion was that even with the faster execution times of
Fortran, C++ is preferred due to the poor support in OOP design patterns and
the programming effort Fortran requires. Then, Modern Fortran was upgraded to
support more design patterns, resulting in a preference for Modern Fortran
over C++. Currently, Fortran is regarded as the indisputable leader for
carrying out highly optimized scalable scientific applications in the shortest
execution times — even more so since the introduction of Modern Fortran.
Moreover, Fortran nowadays provides excellent interoperability with C and C++.
Thus, applications written in Fortran, specifically with Modern Fortran
object-oriented design patterns, can include modules written in C and C++ —
providing the best from both worlds.

Although Modern Fortran does not fully support all object-oriented features
such as generic types, pointers, and polymorphic arrays, these design patterns
are not necessarily needed. Instead, they may be replaced with other available
programming solutions — which usually require slightly higher programming
effort. On the other hand, the lack of available abstractions offered in
Modern Fortran guarantees the relatively thin implementation for OOP support.
This abstraction deficiency demands higher programming effort but allows for
an execution time reduction in general and in scientific codes specifically
[69].

###  II-B Transparent Parallel Class

Since the introduction of object-oriented programming, principles such as
modularity, abstraction, encapsulation, and hierarchy have been encouraged.
These principles make it easier to maintain and modify large software systems.
Modular programming requires that the code follow decomposability and
composability, among other concepts, which promote coupling data structures
with their related functions within classes. These concepts allow new
scientific applications’ developers and users to easily understand these
classes’ functionality and the code’s general structure, subsequently
expanding them with additional functionalities [70, 71].

Modularity is beneficial for developing and composing new scientific
applications and for decomposing conceptually independent implementations from
each other. This is significantly useful for separating the implementations of
parallel modules from the rest of the scientific application. Following this
idea, the parallel implementation is abstracted by a class that can be updated
locally instead of along the entire code to support, e.g., a new MPI standard.
Abstraction of this sort also allows to optimize the parallel scheme and tune
it tightly to given hardware in a single code region.

OOP is a powerful paradigm, but as a rule of thumb, it does not necessarily
conform with performance. For example, abstract classes are great for defining
functionality without a concrete implementation (such as EoS calculations).
However, if a class that implements an abstract class is called multiple times
in a single loop iteration, the memory and computation overhead, derived from
the inheritance, becomes significant, which harms the total execution time. In
contrast, a simple switch-case clause is thin computationally but lacks the
benefits of OOP. Thus, creating a scalable scientific code based on OOP
principles must be done carefully and with the right balance between design
patterns and scalable execution time.

ScalSALE follows Object-Oriented Parallel Programming (OOPP) ideas which
resemble other OOP-based scientific applications such as [72, 73], in which
objects, specifically the data structures, are intrinsically parallel. Thus,
developers can synchronize the data structures with ease. Nevertheless,
developing an efficient and scalable scientific code may and usually will
violate important OOP principles, such as object encapsulation, functional
independence, and modularity. This trade-off is investigated with caution, but
an OOP solution is often preferred if it does not come with substantial
running-time growth. Other works show OOPP variants that emphasize the
importance of introducing OOP into HPC applications [74, 75], and thus keeping
the benefits of OOP (without breaking too-many OOP principles) while
preserving the efficiency required from modern scientific codes.

##  III Implementation

ScalSALE is designed as a computational infrastructure for scalable multi-
material and multi-physical numerical benchmarks. ScalSALE is written in
Modern Fortran [76] (See subsection II-A) and implemented using simple yet
practical design patterns from OOPP (See subsection II-B). Further
implementation details and important concepts are presented.

###  III-A High-performance Code Written in OOP

As discussed in section II, the foundations of ScalSALE paradigm is a balance
between high programming effort and fast codes to easier to use and slower
codes. ScalSALE as a benchmark was developed following general OOP principles,
with exceptions that favor high-performance over OOP paradigms. Following are
representative exceptions for each OOP main principle.

  1. 1.

Encapsulation: Usually, only the object has the authority to access and modify
its data attributes directly. The access to the object’s data from other
objects is done via public methods such as get/set. These functions might
incur significant time and memory overheads since they require copying and
setting the entire data attribute for each access. Since data attributes in
scientific applications may be very large, this functionality is not desired.
A workaround adapted in ScalSALE is to use Point_to_data functions, which
return a pointer directly to the data attribute from outside the object. This
eliminates the need to copy the data and serves as a hand-shake between the
object that owns the data attribute with the code segment that accesses the
data.

  2. 2.

Polymorphism: Modern Fortran does not support array of polymorphic types. This
functionality is useful for defining the boundary conditions scheme for each
edge of the mesh. A fully polymorphic solution would have had an abstract base
class of cell boundary conditions and concrete implementations of cell
boundary conditions’ calculations such as slip cell and free surface classes.
Since cell boundary conditions on a given edge of the mesh do not change
during the simulation, it might be convenient to initialize the specific cell
boundary condition class at the beginning of the program and store an array
with a polymorphic pointer to the specific cell boundary condition class for
each edge of the mesh. A simple workaround would be to define a wrapper type
with a single attribute, a polymorphic pointer, to the base cell boundary
conditions class. An array of these wrapper types can implement the edge
boundary conditions array.

  3. 3.

Abstraction: Following this principle, objects should be exclusively
responsible for their functionalities and data attributes. For example, it
would be preferable for each physical quantity class in ScalSALE to contain a
method called calculate. However, since each physical quantity’s calculation
should accept different quantities as parameters, and since the calculation
depends on the actual physical/numerical scheme, these calculations are kept
in a more global place. For example, Calculate_density is placed in Hydro and
not within the density class.

  4. 4.

Inheritance: This principle allows classes to inherit from base classes and
reuse or overload some basic functionalities from that base class. This
property reduces programming time and requires more thorough computation and
data analysis. However, since overloading functions come with an overhead
cost, this functionality should be used cautiously. In ScalSALE there are
scenarios in which inheritance is used and functions overloading is carried
out, and in other cases in which a switch–case clause is used to determine the
desired functionality (e.g., for determining whether to perform 2D or 3D
rezoning).

###  III-B ScalSALE as an Infrastructure

ScalSALE is implemented as two layers that are independent yet complementary
to one another: the kernel layer and the application layer. This design
provides a simple approach for expanding and modifying each of the two layers
independently. An overview of the UML of ScalSALE is presented in Figure 3.

Figure 3: ScalSALE schematic UML. The classes in blue are associated with the
application layer, while classes in red are associated with kernel.

####  III-B1 The Kernel

Handles the computational building blocks for the structured grid
calculations. It consists of the data structures that are encapsulated in the
classes _Data3d_ and _Data4d_ (as presented in Figure 3), and the
parallelization schemes that are implemented in the Communication and
CommParams classes (as presented in Figure 5, and will be explained in the
following section). The interface between the data structures and the
parallelization schemes is given by functions such as _Point_to_data_ and
_Exchange_blocking_. Furthermore, the kernel layer contains third-party
libraries, enums, linked list implementation, and more.

####  III-B2 The Application

Handles all the physical-related computations such as thermodynamics,
hydrodynamics, mesh rezoning, advection, and future models (as seen in Figure
3). Generally, each physical module should be associated with this layer and
be implemented as an independent class. Associating physical modules with
classes allows easily expanding these physical modules or overloading them
with different schemes. Furthermore, every physical module utilizes different
physical quantities represented by pointers to the relevant objects of
physical quantities. The physical quantities and their related calculations
further assist the encapsulation of the kernel layer from the physical models
in the application layer, enhancing the simplicity and readability of both
layers. Moreover, the numerical computations are written with minimal OOP
intervention for minimal execution time reduction.

###  III-C Communication Implementation

The parallelization scheme in ScalSALE distributes the mesh into cubes or
tiles and uses ghost cells [77] for stencil calculations (in which each cell
depends on its nearest neighbors). Each MPI process has 27, 8, or 2 neighbors
in 3D, 2D, and 1D, respectively. When the processes synchronize in a 3D
setting, $O(N^{2})$ while $N$ is the problem size, values are sent to the
closest face neighbors (presented in green in Figure 4); $O(N)$ values are
sent to the edge neighbors (presented in red in Figure 4) and $O(1)$ values
are sent to the corner neighbors (presented in blue in Figure 4).

03614725891215101316111417182124192225202326 Figure 4: 3D communication
topology structure of an MPI rank — 27 neighbors. The 2D communication
topology can be derived straightforwardly by limiting the ranks to a 2D
rectangle — ranks 9 until 17, for example.

The communication-related operations utilize the MPI-3 standard. The MPI-3
standard introduced important features such as an extension of collective
operations to non-blocking versions, wider support for Fortran, C and C++
binding, point-to-point communications, user-defined communication topologies,
and thread-safe operations. Each feature contributes to the scalability of a
parallel (even hybrid parallelization schemes) framework written in Modern
Fortran.

ScalSALE intent is to provide a scalable benchmark framework. Therefore, the
synchronization overhead should be minimized as much as possible. To comply
with this idea, non-blocking communication, and a communication topology were
introduced to ScalSALE:

####  III-C1 Non-blocking Communication

Naïvely, ghost cells should be synchronized across the ranks after each update
of a physical quantity. The update process consists of sending and receiving
MPI messages to the corresponding neighbors. The communication overhead
increases as the number of ranks increases due to the higher probability of a
delayed message in the system. A possible approach to diminish such a case is
by introducing non-blocking communication to the code — validating that the
received messages indeed arrived right before they are being used. The non-
blocking communication grants the ability to perform other calculations
between the two points of synchronization (send and receive). Thus, if an MPI
message is delayed, the program performs other calculations instead of being
idle, which minimizes the communication overhead. ScalSALE implements non-
blocking communications via the asynchronous collective MPI operation
_Ineighbor_alltoall_.

####  III-C2 Topology

In supercomputers, cores on the same node transfer messages via fast shared-
memory communication, while cores from different nodes are connected via
intercommunication devices such as Infiniband and Ethernet, which are slower
by factors or orders of magnitude, respectively. Therefore, ranks that
communicate frequently (for example, the green ranks in Figure 4) should be
assigned to cores that reside in the same node — thus, lowering the
communication overhead.

In ScalSALE, this idea is implemented via the communication topology function
that was presented in MPI-3, _MPI_Dist_graph_create_adjacent_. This function
creates a graph topology between the appropriated processes and initializes
the communicator. Furthermore, the function receives a weight array that hints
at the communication load between the ranks. In turn, the function may migrate
the ranks (change their original identifiers) such that the ranks that heavily
communicate will reside in the same node.

Figure 5: ScalSALE UML of the parallel object in ScalSALE.

The parallelization scheme is associated with the kernel layer and is
implemented via three classes (see Figure 5 a schematic overview of the three
classes):

  1. 1.

Optimized communication — _Communication_ class: Handles and abstracts the
communication operations. Each function implements a communication operation,
for example, _Send_ , _Recv_ and _Send_Recv_. These functions also handle the
synchronization process of the ghost cells, which, as mentioned, is
implemented via the _MPI_Ineighbor_alltoall_ function [78].

  2. 2.

Topology and ranks — _CommParams_ class: Initializes the communicator with the
_MPI_Dist_graph_create_adjacent_ function. Furthermore, the class handles the
relevant parameters that describe the MPI messages, i.e., size of messages,
the communication topology (communicator), and more.

  3. 3.

Transparent communication — _Data3d_ /_Data4d_ class: Implements the basic
data structure of ScalSALE. As such, it is also an inherently parallel data
class that provides a transparent interface for exchanging and synchronizing
data across the relevant ranks. Hence, encapsulating the parallel
implementation from the rest of the program. The _Data3d_ and _Data4d_ classes
provides three exchanging operations: a blocking operation
(_Exchange_blocking_), a non-blocking send operation (_Exchange_nonblocking)_
and a non-blocking receive operation (_Exchange_end_). In addition, it
contains the values of the data and the communication parameters via
_CommParams_ as a variable. The CommParams variable is provided by the
physical quantity, thus, hiding the physical quantity type from _Data3d_ and
_Data4d_ classes.

###  III-D A Case Study — Stress Module

⬇

type, extends(cell_quantity_t) :: stress_tensor_t

…

end type stress_tensor_t

(a) The class describing the new stress tensor quantity, it inherent from the
cell quantity base type.

⬇

module material_module

type, extends(material_base_t) :: material_t

integer, dimension(:) :: stress_module_type

! mat dependent module type (steinberg, perfect, etc)

contains

procedure :: Apply_stress_module

!calculations of the stress yield and shear

end type material_t

contains

subroutine Apply_stress_module(this, shear, yield)

implicit none

class(material_t), intent(in out) :: this

type(shear_modulus_t), pointer, intent(in out) :: shear

! the shear modulus for every material

type(stress_yield_t), pointer, intent(in out) :: yield

! the yield for every material

integer :: i, j, k, m

do m=1, nmats ! loop over all materials

do k=1, nz ! in the loops we calculate shear

! and yield for the steinberg model

do j=1, ny

do i=1, nx

! yield and shear calculations using the

! Steinberg stress hardening model

eta = density0(i,j,k)/density(i,j,k)

shear(m,i,j,k) = shear0+GP_steinberg(m)* &

pressure(i,j,k)*eta**(1d0/3)+GT_steinberg(m)* &

(temperature(i,j,k)-init_temperature(i,j,k))

yield(m,i,j,k) = min(Y0(m)*(1+BETA_steinberg(m)* &

plastic_strain(i,j,k)) ** N_steinberg(m) &

,YIELD_max(m))*(shear(m,i,j,k)/shear0(m))

end do

end do

end do

end do

!the MPI synchronization

call shearcall yieldend subroutine Apply_stress_module

…

end module material_module

(b) This code describes adding a new material-dependent calculation in the
material class. The numerical stress model calculations of shear and yield use
the Steinberg stress hardening model [79].

⬇

type :: stress_model

type(shear_modulus_t) , pointer :: shear

type(stress_yield_t) , pointer :: yield

type(stress_tensor_t) , pointer :: stress_tensor

contains

procedure :: Calculate_stress

procedure :: Update_acceleration

end type stress_model

(c) The class described the new physical stress module. It contains two main
functions, a function that calculates the stress in the current time step and
a function that update the acceleration values using the updated stress
tensor.

Figure 6: A schematic code of stress physical model addition to ScalSALE.

Physical attributes such as stress play an important role in describing the
motion of solids in hydrodynamic codes. In non-compressible, visco-elastic, or
visco-inelastic fluids, forces between the different layers of the fluid
result in stress forces that oppose the fluid flow. The stress is strain and
velocity dependent, influencing the fluid’s momentum distribution and other
hydrodynamic properties [80]. Hence, the stress-strain-velocity relation,
which is material-dependent, is crucial for describing the motion of visco-
fluids correctly.

As a case study, the following workflow describes how the stress module can be
easily added to ScalSALE relatively. Firstly, adding new physical quantities,
such as the stress tensor that inherits from the cell quantity base class, to
the application layer (as seen in 6(a) in lines 1–3). Afterward, implementing
the yield and shear modulus calculations [81] to the application layer. The
physical calculations are added as new procedures in the materials class as
seen in 6(b) from lines 4–44. In lines 40–41, the parallelization
synchronization is done using a blocking manner. Finally, adding the physical
stress calculation [82], and the correction to the acceleration due to stress,
in a new stress model class in the application layer, as seen in the code in
6(c) in lines 40–47.

##  IV Results

To evaluate the potential of ScalSALE as a scalable benchmark, it is compared
to Lulesh [83] — the well-known mini-app benchmark. Lulesh is a hydrodynamic
code, highly optimized, that demonstrates excellent scaling capabilities.
Moreover, both codes can simulate the same physical problem — 3D Sedov-Taylor,
providing comparable numerical workloads. Thus, it is the most suitable for
comparing and evaluating ScalSALE.

Most (if not all) applications that aim to be executed in supercomputers are
parallelized with MPI due to the need to scale enormously beyond a single
node. However, the introduction of heterogeneous systems poses a challenge in
obtaining maximum performance via MPI exclusively. As a result, the newer
standards of MPI (especially MPI-3, the standard that is used in ScalSALE)
supply comprehensive support and even inherently supports shared-memory
environments, hence its popularity among hybrid parallelization paradigms
(MPI+X) [84]. Nevertheless, evaluating the distributed parallelization scheme
scaling capabilities exclusively (MPI) is a necessary step that must be taken
independently to MPI+X. Therefore, although Lulesh implements several
parallelization paradigms, the MPI version was exclusively enabled and
compared to ScalSALE. At the same time, ScalSALE could be easily extended to
other paradigms without needing the transparent parallel class (for example,
inserting an OpenMP directive in line 21 in 6(b)).

ScalSALE can simulate physical problems with Lagrangian or Eulerian rezoning
(with an advection phase), contrary to Lulesh which can simulate with
Lagrangian rezoning exclusively. Thus, both rezoning must be considered to
evaluate the full scaling capabilities of ScalSALE. Executing ScalSALE with
Lagrangian rezoning provides a more comparable execution time to Lulesh while
Eulerian rezoning provides a more stressful and realistic calculation as the
advection phase results in additional numerical calculations and
synchronization points.

Another point of reference between the two physical schemes is that in
Lagrangian rezoning, there are approximately 45 calls to a function in
_Communication_ i.e., 45 MPI collective operations per cycle, while Eulerian
rezoning performs 60 calls per cycle. Each call to a function performs a
communication operation of a single physical quantity. In contrary, Lulesh
performs approximately 5 calls to a similar function per cycle (basic
_MPI_IRecv_ and _MPI_ISend_ functions), while each call transfers multiple
physical quantity at once. Moreover, as ScalSALE is based on the multi-
material SALE scheme, it contains additional material-related quantities. As a
result, ScalSALE not only performs more communication operations but sends
more data per cycle. Nevertheless, due to the utilization of non-blocking
collective communication and topology distribution via the MPI-3 standard,
ScalSALE achieves comparable results to Lulesh.

The weak and strong scaling tests are used to compare ScalSALE and Lulesh. The
strong scaling test is performed by creating a fixed size problem, considering
it has a sufficient workload and increasing the number of cores. The weak
scaling test is performed by creating a fixed size problem per core and
increasing the number of cores. Both scaling tests are important to understand
the scaling capabilities of the parallelization scheme and, in turn, the
communication overhead of the code. The efficiency is calculated as follows:

$Efficiency=\begin{cases}\frac{Time_{1}}{Time_{n}}&\text{Weak scaling,}\\\
\frac{Time_{1}}{n\cdot Time_{n}}&\text{Strong scaling}\end{cases}$

whereas $Time_{1}$ is the execution time of a single core, $Time_{n}$ is the
execution time where $n$ is equal to the number of cores.

In Lagrangian rezoning, the amount of numerical operations per cycle is
constant. Thus, the evaluation was carried out over only 20 cycles, without
the initialization time (similar to [85]). In contrast, in Eulerian rezoning,
the number of operations per cycle fluctuates due to the advection phase.
Hence, 200 cycles were performed for more precise estimation.

The two codes were compiled with Intel’s oneAPI [86] MPI Fortran (Modern
Fortran 2018), and C++ compilers [87]. oneAPI provides wide support to
multiple programming languages [88] and parallelization paradigms, especially
for heterogeneous nodes, thus allowing further local implementation of MPI+X
strategy using Intel’s oneAPI OpenMP* [89]. The codes were compiled with no
compiler optimization (O0) to allow a fair comparison of the raw scaling
capabilities. The codes were executed on NegevHPC [90] on Intel Xeon Gold 6130
compute nodes. Each node has 32 cores and 128GB of RAM.

###  IV-A Strong Scale

$10^{0}$$10^{1}$$10^{2}$$10^{3}$$10^{1}$$10^{2}$$10^{3}$$10^{4}$#CoresTime
(sec)Runtime, strong scaleLULESHScalSALE \- LagrangeScalSALE \- Euler (a) The
execution times of ScalSALE and Lulesh on the strong scaling test.

$0$$200$$400$$600$$800$$1{,}000$$0.6$$0.7$$0.8$$0.9$$1$#CoresEfficiencyEfficiency,
strong scaleLULESHScalSALE \- LagrangeScalSALE \- Euler (b) The efficiencies
of ScalSALE and Lulesh on the strong scaling test.

$0$$200$$400$$600$$800$$1{,}000$$200$$400$$600$#CoresTime (sec)Runtime, weak
scaleLULESHScalSALE \- LagrangeScalSALE \- Euler (c) The execution times of
ScalSALE and Lulesh on the weak scaling test.

$0$$200$$400$$600$$800$$1{,}000$$0.85$$0.9$$0.95$$1$#CoresEfficiencyEfficiency,
weak scaleLULESHScalSALE \- LagrangeScalSALE \- Euler (d) The efficiencies of
ScalSALE and Lulesh on the weak scaling test.

Figure 7: ScalSALE and Lulesh comparison on the strong and weak scaling test.

For the strong scale testing, a 3D Sedov-Taylor problem with a mesh size of
$480^{3}=110,592,000$ (480 cells in each axis) was defined. The total memory
consumption of the problem is 80GB in ScalSALE and 71GB in Lulesh. 7(a)
presents the total execution time of the strong scaling problem and 7(b)
presents the strong scaling efficiency.

As expected, the execution time of Lulesh is the shortest as it is a highly
optimized dedicated Sedov-Taylor hydrodynamic computation software.
Nevertheless, the execution time of ScalSALE with Lagrangian rezoning is
comparable to Lulesh.

Although Eulerian rezoning achieves the slowest execution time, it still
presents excellent scaling capabilities even when normalized per cycle. Most
importantly, Eulerian rezoning provides a broader range of hydrodynamic
problems to be solved while maintaining the same speedup presented for the
Sedov-Taylor problem.

A decrease in efficiency can be spotted at 1–8 cores due to the transition
from a serial execution with no communication overhead to a distributed
execution. Due to the nature of the parallelization scheme in 3D simulations,
there is another decrease in efficiency at 8–27 cores, as 27 cores are the
minimal number of processes to achieve the maximum amount of communication
(exactly one process has 27 neighbors). Next, a slight decrease in efficiency
can also be spotted at 27–64 cores for all three executions due to the
introduction of a slightly slower medium that connects processes from
different compute nodes. The final decrease in efficiency occurs at 216 cores
(dashed line). For a fixed problem size, as the number of cores increases, the
numerical workload of each MPI process decreases — i.e., the number of cells
per MPI process decreases. Although the message size also decreases (message
size is of the size $O(N^{2})$), it is less dramatic than the decrease in the
numerical workload (of size $O(N^{3})$). Hence, the communication overhead
outweighs the numerical calculations, and the efficiencies of both benchmarks
decrease.

In addition, at 216–1,000 cores, due to the additional synchronization
operations in ScalSALE, the margin between the two codes increase as the
workload decreases faster than the communication overhead.

Nonetheless, ScalSALE presents excellent strong scaling capabilities in both
Lagrangian and Eulerian rezoning compared to Lulesh.

###  IV-B Weak Scale

For the weak scale comparison, the mesh of each MPI process was set to
$120^{3}=1,728,000$ (120 cells in each axis), similarly to [91]. A mesh size
of $120^{3}$ is equivalent to approximately 1.7GB in ScalSALE and 1.1GB in
Lulesh. As each computational node contains 32 cores, at full capacity, each
node would consume up to 55GB in ScalSALE and 35GB in Lulesh. 7(c) presents
the total execution time of the weak scaling problem, and 7(d) presents the
weak scaling efficiency.

Similar to the strong scaling, Lulesh achieves the best execution time, while
ScalSALE with Lagrangian rezoning achieves a comparable execution time.
However, contrary to strong scaling, each MPI process computes the same amount
of numerical operations in weak scaling. Thus, the efficiency does not exhibit
the a decrease in the efficiency. Nevertheless, the efficiency does tend to
decrease slightly due to the additional cores that need to communicate with
one another.

Moreover, as the MPI implementation is generic and associated with the data
structures and not the physical calculation, it has the potential to exhibit
the same speedups on other problems in different domains. Thus achieving its
targeted goal as a scalable benchmark framework.

##  V Conclusion

In this work, the problem with existing benchmarks was tackled by introducing
a new scalable benchmark framework named ScalSALE. The main objective of such
a benchmark is to provide a more representative benchmark for scientific
applications while maintaining simplicity. Contrary to existing benchmarks,
ScalSALE is flexible and designed such that it is easy to introduce new
physical models, all while maintaining scalable speedups and low execution
time.

ScalSALE is based on the well-known SALE scheme and is currently capable of
solving the 1D, 2D, and 3D hydrodynamic equations in a structured grid with
multi-material ALE capabilities (with an advection phase). A code with multi-
material and ALE capabilities has the basic foundations for a code similar to
hydrodynamic scientific applications.

ScalSALE is written in Modern Fortran, utilizing object-oriented design
patterns that provide the infrastructure to further develop the benchmark as a
framework and the fast execution times of Fortran to support the code as
optimized and scalable. Furthermore, Modern Fortran provides a well-balanced
combination of fast execution times, especially for numerical computations,
and sufficient support for object-oriented design patterns that allow easy
modification and expansion of code.

In order to further support the code as a framework, ScalSALE is implemented
with a two-layer design, the kernel and the application. The former handles
the inherently parallel data structures (consequently, the physical
quantities) that provide a transparent interface to the implementation of the
distributed parallelization scheme. While the latter handles the numerical
implementation and utilizes the parallel data structures to synchronize and
keep up-to-date values.

ScalSALE parallelization scheme is based on the MPI-3 standard and utilizes
two important features necessary for large-scale executions — non-blocking
collective communication and the processes/cores topology. Moreover, the MPI-3
standard provides wide support for hybrid parallelization schemes. Along with
the separation of the numerical computations and the MPI implementation, MPI+X
schemes can be added easily.

To evaluate ScalSALE scaling capabilities, it was tested via the weak and
strong scaling tests and compared to the well-known benchmark — Lulesh.
Despite performing more communication operations, ScalSALE achieves excellent
and comparable scaling results compared to Lulesh in both rezoning options.
Moreover, as the parallelization scheme’s implementation is independent of the
numerical simulation, ScalSALE can exhibit the same scaling capabilities for
other physical domains. Hence, ScalSALE achieves its targeted goal as a multi-
material, multi-physical scalable benchmark framework.

## Acknowledgment

This work was supported by the Pazy foundation and the Lynn and William
Frankel Center for Computer Science. Computational support was provided by the
NegevHPC project [90].

## References

  * [1] K. G. Wilson, “Grand challenges to computational science,” _Future Generation Computer Systems_ , vol. 5, no. 2-3, pp. 171–189, 1989.
  * [2] A. L. C. Facility, “Preparing to program aurora at exascale.” [Online]. Available: https://www.iwocl.org/wp-content/uploads/iwocl-syclcon-2020-finkel-keynote-slides.pdf
  * [3] P. Bientinesi, J. A. Gunnels, M. E. Myers, E. S. Quintana-Ortí, and R. A. v. d. Geijn, “The science of deriving dense linear algebra algorithms,” _ACM Transactions on Mathematical Software (TOMS)_ , vol. 31, no. 1, pp. 1–26, 2005.
  * [4] W. K. Liu, S. Jun, and Y. F. Zhang, “Reproducing kernel particle methods,” _International journal for numerical methods in fluids_ , vol. 20, no. 8-9, pp. 1081–1106, 1995.
  * [5] E. O. Brigham, _The fast Fourier transform and its applications_. Prentice-Hall, Inc., 1988.
  * [6] M. J. Berger and P. Colella, “Local adaptive mesh refinement for shock hydrodynamics,” _Journal of computational Physics_ , vol. 82, no. 1, pp. 64–84, 1989.
  * [7] G. F. Pfister, “An introduction to the infiniband architecture,” _High performance mass storage and parallel I/O_ , vol. 42, no. 617-632, p. 102, 2001\.
  * [8] J. Marian, C. S. Becquart, C. Domain, S. L. Dudarev, M. R. Gilbert, R. J. Kurtz, D. R. Mason, K. Nordlund, A. E. Sand, L. L. Snead _et al._ , “Recent advances in modeling and simulation of the exposure and response of tungsten to fusion energy conditions,” _Nuclear Fusion_ , vol. 57, no. 9, p. 092008, 2017.
  * [9] W. M. Washington and C. Parkinson, _Introduction to three-dimensional climate modeling_. University science books, 2005.
  * [10] C. R. Kitchin, _Astrophysical techniques_. CRC press, 2008.
  * [11] T. M. Aamodt, W. W. L. Fung, and T. G. Rogers, “General-purpose graphics processor architectures,” _Synthesis Lectures on Computer Architecture_ , vol. 13, no. 2, pp. 1–140, 2018.
  * [12] H. Jun, J. Cho, K. Lee, H.-Y. Son, K. Kim, H. Jin, and K. Kim, “Hbm (high bandwidth memory) dram technology and architecture,” in _2017 IEEE International Memory Workshop (IMW)_. IEEE, 2017, pp. 1–4.
  * [13] N. Liu, J. Cope, P. Carns, C. Carothers, R. Ross, G. Grider, A. Crume, and C. Maltzahn, “On the role of burst buffers in leadership-class storage systems,” in _2012 IEEE 28th Symposium on Mass Storage Systems and Technologies (MSST)_. IEEE, 2012, pp. 1–11.
  * [14] A. Khan, H. Sim, S. S. Vazhkudai, A. R. Butt, and Y. Kim, “An analysis of system balance and architectural trends based on top500 supercomputers,” in _The International Conference on High Performance Computing in Asia-Pacific Region_ , 2021, pp. 11–22.
  * [15] Top500, “Top#500 website,” https://www.top500.org/., 2019, [Online].
  * [16] D. Monroe, “Fugaku takes the lead,” _Communications of the ACM_ , vol. 64, no. 1, pp. 16–18, 2020.
  * [17] J. Dongarra, M. A. Heroux, and P. Luszczek, “High-performance conjugate-gradient benchmark: A new metric for ranking high-performance computing systems,” _The International Journal of High Performance Computing Applications_ , vol. 30, no. 1, pp. 3–10, 2016.
  * [18] R. Stevens, J. Ramprakash, P. Messina, M. Papka, and K. Riley, “Aurora: Argonne’s next-generation exascale supercomputer,” Argonne National Lab.(ANL), Argonne, IL (United States), Tech. Rep., 2019.
  * [19] D. Schneider, “The exascale era is upon us: The frontier supercomputer may be the first to reach 1,000,000,000,000,000,000 operations per second,” _IEEE Spectrum_ , vol. 59, no. 1, pp. 34–35, 2022.
  * [20] Meta, “Meta’s supercomputer.” https://ai.facebook.com/blog/ai-rsc/, 2022, [Online].
  * [21] R. Walkup, S. R. Seelam, and S. Wen, “Best practices for hpc workloads on public cloud platforms: A guide for computational scientists to use public cloud for hpc workloads,” in _Proceedings of the 2022 ACM/SPEC on International Conference on Performance Engineering_ , 2022, pp. 29–35.
  * [22] A. Gupta and D. Milojicic, “Evaluation of hpc applications on cloud,” in _2011 Sixth Open Cirrus Summit_. IEEE, 2011, pp. 22–26.
  * [23] A. Ali and K. S. Syed, “An outlook of high performance computing infrastructures for scientific computing,” in _Advances in Computers_. Elsevier, 2013, vol. 91, pp. 87–118.
  * [24] O. Lascu, K. Bosworth, H. Kitzhšfer, J. K. N. Mroz, G. Noble, and F. Pizzano, _HPC Clusters using Infiniband on IBM Power Systems Servers_. IBM, 2009.
  * [25] H. Shan and J. Shalf, “Using ior to analyze the i/o performance for hpc platforms,” Lawrence Berkeley National Lab.(LBNL), Berkeley, CA (United States), Tech. Rep., 2007.
  * [26] V. Morozov, J. Meng, V. Vishwanath, J. R. Hammond, K. Kumaran, and M. E. Papka, “Alcf mpi benchmarks: Understanding machine-specific communication behavior,” in _2012 41st International Conference on Parallel Processing Workshops_. IEEE, 2012, pp. 19–28.
  * [27] A. Moody, “Mpigraph,” Lawrence Livermore National Lab.(LLNL), Livermore, CA (United States), Tech. Rep., 2007.
  * [28] G. Juckeland, W. Brantley, S. Chandrasekaran, B. Chapman, S. Che, M. Colgrove, H. Feng, A. Grund, R. Henschel, W.-M. W. Hwu _et al._ , “Spec accel: A standard application suite for measuring hardware accelerator performance,” in _International Workshop on Performance Modeling, Benchmarking and Simulation of High Performance Computer Systems_. Springer, 2014, pp. 46–67.
  * [29] P. Wong and R. Der Wijngaart, “Nas parallel benchmarks i/o version 2.4,” _NASA Ames Research Center, Moffet Field, CA, Tech. Rep. NAS-03-002_ , p. 91, 2003.
  * [30] HPL, “HPL website,” https://icl.utk.edu/hpl/index.html., [Online].
  * [31] J. Dongarra and M. A. Heroux, “Toward a new metric for ranking high performance computing systems,” _Sandia Report, SAND2013-4744_ , vol. 312, p. 150, 2013.
  * [32] V. G. Vergara Larrea, W. Joubert, M. J. Brim, R. D. Budiardja, D. Maxwell, M. Ezell, C. Zimmer, S. Boehm, W. Elwasif, S. Oral, C. Fuson, D. Pelfrey, O. Hernandez, D. Leverman, J. Hanley, M. Berrill, and A. Tharrington, “Scaling the summit: Deploying the world’s fastest supercomputer,” in _High Performance Computing: ISC High Performance 2019 International Workshops, Frankfurt, Germany, June 16-20, 2019, Revised Selected Papers_. Berlin, Heidelberg: Springer-Verlag, 2019, p. 330–351. [Online]. Available: https://doi.org/10.1007/978-3-030-34356-9_26
  * [33] LLNL, “CORAL benchmark website,” https://asc.llnl.gov/coral-2-benchmarks, [Online].
  * [34] B. Bland, “Titan-early experience with the titan system at oak ridge national laboratory,” in _2012 SC Companion: High Performance Computing, Networking Storage and Analysis_. IEEE, 2012, pp. 2189–2211.
  * [35] U. Yang, R. Falgout, and J. Park, “Algebraic multigrid benchmark,” Lawrence Livermore National Lab.(LLNL), Livermore, CA (United States), Tech. Rep., 2017\.
  * [36] D. F. Richards, R. C. Bleile, P. S. Brantley, S. A. Dawson, M. S. McKinley, and M. J. O’Brien, “Quicksilver: a proxy app for the monte carlo transport code mercury,” in _2017 IEEE International Conference on Cluster Computing (CLUSTER)_. IEEE, 2017, pp. 866–873.
  * [37] LLNL, “Mercury Code webpage.” https://wci.llnl.gov/simulation/computer-codes/mercury, 2017, [Online].
  * [38] ——, “LULESH webpage.” https://asc.llnl.gov/codes/proxy-apps/lulesh, 2012, [Online].
  * [39] C. R. Noble, A. T. Anderson, N. R. Barton, J. A. Bramwell, A. Capps, M. H. Chang, J. J. Chou, D. M. Dawson, E. R. Diana, T. A. Dunn _et al._ , “Ale3d: An arbitrary lagrangian-eulerian multi-physics code,” Lawrence Livermore National Lab.(LLNL), Livermore, CA (United States), Tech. Rep., 2017\.
  * [40] J. R. Kamm, “Evaluation of the sedov-von neumann-taylor blast wave solution,” _Astrophys. J. Suppl.(submitted)_ , p. 46, 2000.
  * [41] S. Wen, L. Cherkasova, F. X. Lin, and X. Liu, “Profdp: A lightweight profiler to guide data placement in heterogeneous memory systems,” in _Proceedings of the 2018 International Conference on Supercomputing_ , 2018, pp. 263–273.
  * [42] M. Feldman, “DARPA Sets Ubiquitous HPC Program in Motion.” https://www.hpcwire.com/2010/08/10/darpa_sets_ubiquitous_hpc_program_in_motion/, 2010, [Online].
  * [43] G.-T. Bercea, C. Bertolli, S. F. Antao, A. C. Jacob, A. E. Eichenberger, T. Chen, Z. Sura, H. Sung, G. Rokos, D. Appelhans _et al._ , “Performance analysis of openmp on a gpu using a coral proxy application,” in _Proceedings of the 6th International Workshop on Performance Modeling, Benchmarking, and Simulation of High Performance Computing Systems_ , 2015, pp. 1–11.
  * [44] I. Karlin, J. McGraw, J. Keasler, and B. Still, “Tuning the lulesh mini-app for current and future hardware,” Lawrence Livermore National Lab.(LLNL), Livermore, CA (United States), Tech. Rep., 2013.
  * [45] R. Hornung, J. Keasler, and M. Gokhale, “Hydrodynamics challenge problem,” Lawrence Livermore National Lab.(LLNL), Livermore, CA (United States), Tech. Rep., 2011.
  * [46] D. Laney, S. Langer, C. Weber, P. Lindstrom, and A. Wegener, “Assessing the effects of data compression in simulations using physically motivated metrics,” in _SC’13: Proceedings of the International Conference on High Performance Computing, Networking, Storage and Analysis_. IEEE, 2013, pp. 1–12.
  * [47] M. Metcalf, J. K. Reid, and M. Cohen, _Fortran 95/2003 Explained_. Oxford university press Oxford, 2004, vol. 3.
  * [48] J. R. Cary, S. G. Shasharina, J. C. Cummings, J. V. Reynders, and P. J. Hinker, “Comparison of c++ and fortran 90 for object-oriented scientific programming,” _Computer Physics Communications_ , vol. 105, no. 1, pp. 20–36, 1997.
  * [49] S. Alam, “Is fortran still relevant? comparing fortran with java and c++,” _arXiv preprint arXiv:1407.2190_ , 2014.
  * [50] J. L. Gustafson, “Reevaluating amdahl’s law,” _Communications of the ACM_ , vol. 31, no. 5, pp. 532–533, 1988.
  * [51] E. Loh, “The ideal hpc programming language,” _Queue_ , vol. 8, no. 6, p. 30, 2010.
  * [52] W. Vanderbauwhede, “Making legacy fortran code type safe through automated program transformation,” _The Journal of Supercomputing_ , vol. 78, no. 2, pp. 2988–3028, 2022.
  * [53] J. M. Perkel, “Ten computer codes that transformed science.” _Nature_ , vol. 589, no. 7842, pp. 344–349, 2021.
  * [54] A. Amsden, H. Ruppel, and C. Hirt, “Sale: A simplified ale computer program for fluid flow at all speeds,” _LANL Rep. LA-8095_ , 06 1980.
  * [55] A. Amsden and H. Ruppel, “Sale-3d: A simplified ale computer program for calculating three-dimensional fluid flow,” Los Alamos National Lab., NM (USA), Tech. Rep., 1981.
  * [56] G. A. Sod, “A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws,” _Journal of computational physics_ , vol. 27, no. 1, pp. 1–31, 1978.
  * [57] W. Noh, “Artificial viscosity (q) and artificial heat flux (h) errors for spherically divergent shocks,” Lawrence Livermore National Lab., CA (USA), Tech. Rep., 1983.
  * [58] M. Feathers, _Working Effectively with Legacy Code_. Prentice Hall Professional, 2004.
  * [59] M. Rusanovsky, R. Harel, L.-o. Alon, I. Mosseri, H. Levin, and G. Oren, “Backus: Comprehensive high-performance research software engineering approach for simulations in supercomputing systems,” _arXiv preprint arXiv:1910.06415_ , 2019.
  * [60] L. H. Howell, B. T. Gunney, and A. Bhatele, “Characterization of proxy application performance on advanced architectures. umt2013, mcb, amg2013,” Lawrence Livermore National Lab.(LLNL), Livermore, CA (United States), Tech. Rep., 2015.
  * [61] A. R. Conn, G. Gould, and P. L. Toint, _LANCELOT: a Fortran package for large-scale nonlinear optimization (Release A)_. Springer Science & Business Media, 2013, vol. 17.
  * [62] S. E. Koonin, _Computational physics: Fortran version_. CRC Press, 2018.
  * [63] R. J. Hanson and T. Hopkins, _Numerical Computing with Modern Fortran_. SIAM, 2013, vol. 134.
  * [64] A. M. Gorelik, “Object-oriented programming in modern fortran,” _Programming and Computer Software_ , vol. 30, no. 3, pp. 173–179, 2004.
  * [65] J. E. Moreira, S. P. Midkiff, and M. Gupta, “A comparison of java, c/c++, and fortran for numerical computing,” _IEEE Antennas and Propagation Magazine_ , vol. 40, no. 5, pp. 102–105, 1998.
  * [66] T. L. Veldhuizen and M. Jernigan, “Will c++ be faster than fortran?” in _International Conference on Computing in Object-Oriented Parallel Environments_. Springer, 1997, pp. 49–56.
  * [67] J. E. Moreira, S. P. Midkiff, and M. Gupta, “A comparison of java, c/c++, and fortran for numerical computing,” _IEEE Antennas and Propagation Magazine_ , vol. 40, no. 5, pp. 102–105, 1998.
  * [68] B. Mösli, “A comparison of c++, fortran 90 and oberon-2 for scientific programming,” in _GISI 95_. Springer, 1995, pp. 740–748.
  * [69] S. Arabas, D. Jarecka, A. Jaruga, and M. Fijałkowski, “Formula translation in blitz++, numpy and modern fortran: A case study of the language choice tradeoffs,” _Scientific Programming_ , vol. 22, no. 3, pp. 201–222, 2014\.
  * [70] B. Stroustrup, _The C++ programming language_. Pearson Education India, 2000.
  * [71] D. Rouson, J. Xia, and X. Xu, _Scientific software design: the object-oriented way_. Cambridge University Press, 2011.
  * [72] E. G. Pinho and F. H. de Carvalho Junior, “An object-oriented parallel programming language for distributed-memory parallel computing platforms,” _Science of Computer Programming_ , vol. 80, pp. 65–90, 2014.
  * [73] C. Calvin, O. Cueto, and P. Emonot, “An object-oriented approach to the design of fluid mechanics software,” _ESAIM: Mathematical Modelling and Numerical Analysis_ , vol. 36, no. 5, pp. 907–921, 2002.
  * [74] X. Fan, “Towards object-oriented openmp parallel programming,” Ph.D. dissertation, ResearchSpace@ Auckland, 2017.
  * [75] A. Radenski, “Object-oriented programming and parallelism,” 1998.
  * [76] M. Metcalf, J. Reid, and M. Cohen, _Modern fortran explained_. Oxford University Press, 2011.
  * [77] F. B. Kjolstad and M. Snir, “Ghost cell pattern,” in _Proceedings of the 2010 Workshop on Parallel Programming Patterns_ , 2010, pp. 1–9.
  * [78] MPICH, “MPI AlltoAll.” https://www.mpich.org/static/docs/v3.3/www3/MPI_Neighbor_alltoall.html, 2018, [Online].
  * [79] D. Steinberg and C. Lund, “A constitutive model for strain rates from 10- 4 to 106 s-,” _Journal of applied physics_ , vol. 65, no. 4, pp. 1528–1533, 1989\.
  * [80] R. S. Rivlin, “The hydrodynamics of non-newtonian fluids. i,” _Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences_ , vol. 193, no. 1033, pp. 260–281, 1948.
  * [81] D. Steinberg, S. Cochran, and M. Guinan, “A constitutive model for metals applicable at high-strain rate,” _Journal of applied physics_ , vol. 51, no. 3, pp. 1498–1504, 1980.
  * [82] M. L. Wilkins, “Calculation of elastic-plastic flow,” California Univ Livermore Radiation Lab, Tech. Rep., 1963.
  * [83] I. Karlin, A. Bhatele, J. Keasler, B. L. Chamberlain, J. Cohen, Z. DeVito, R. Haque, D. Laney, E. Luke, F. Wang _et al._ , “Exploring traditional and emerging parallel programming models using a proxy application,” in _2013 IEEE 27th International Symposium on Parallel and Distributed Processing_. IEEE, 2013, pp. 919–932.
  * [84] Y. He and H. Ding, “Mpi and openmp paradigms on cluster of smp architectures: the vacancy tracking algorithm for multi-dimensional array transposition,” in _SC’02: Proceedings of the 2002 ACM/IEEE Conference on Supercomputing_. IEEE, 2002, pp. 6–6.
  * [85] K. Dichev and D. S. Nikolopoulos, “Implementing efficient message logging protocols as mpi application extensions,” in _Proceedings of the 26th European MPI Users’ Group Meeting_ , 2019, pp. 1–11.
  * [86] Intel, “A New Era of Accelerated Computing.” https://www.intel.com/content/www/us/en/developer/tools/oneapi/overview.html, 2020, [Online].
  * [87] ——, “Intel® Components,” https://www.intel.com/content/www/us/en/developer/tools/oneapi/fortran-compiler.html#gs.9yqajb, 2020, [Online].
  * [88] ——, “Intel® Components,” https://www.intel.com/content/www/us/en/developer/tools/oneapi/components.html#gs.9yq7fl, 2020, [Online].
  * [89] intel, “oneAPI OpenMP Support.” https://www.intel.com/content/www/us/en/develop/documentation/oneapi-dpcpp-cpp-compiler-dev-guide-and-reference/top/optimization-and-programming/openmp-support.html, 2021, [Online].
  * [90] “NegevHPC Project,” https://www.negevhpc.com, [Online].
  * [91] M. G. Dosanjh, R. E. Grant, P. G. Bridges, and R. Brightwell, “Re-evaluating network onload vs. offload for the many-core era,” in _2015 IEEE International Conference on Cluster Computing_. IEEE, 2015, pp. 342–350.