Richard G. Weiss

OpenFOAM is a CFD software widely used in both industry and academia. The exaFOAM project aims at enhancing the HPC scalability of OpenFOAM, while identifying its current bottlenecks and proposing ways to overcome them. For the assessment of the software components and the code profiling during the code development, lightweight but significant benchmarks should be used. The answer was to develop microbenchmarks, with a small memory footprint and short runtime. The name microbenchmark does not mean that they have been prepared to be the smallest possible test cases, as they have been developed to fit in a compute node, which usually has dozens of compute cores. The microbenchmarks cover a broad band of applications: incompressible and compressible flow, combustion, viscoelastic flow and adjoint optimization. All benchmarks are part of the OpenFOAM HPC Technical Committee repository and are fully accessible. The performance using HPC systems with Intel and AMD processors (x86_64 architecture) and Arm processors (aarch64 architecture) have been benchmarked. For the workloads in this study, the mean performance with the AMD CPU is 62% higher than with Arm and 42% higher than with Intel. The AMD processor seems particularly suited resulting in an overall shorter time-to-solution.

This paper presents a new mesh data layout for parallel I/O of linear unstructured polyhedral meshes. The new mesh representation infers coherence across entities of different topological dimensions, i.e., grid cells, faces, and points. The coherence due to cell-to-face and face-to-point connectivities of the mesh is formulated as a tree data structure distributed across processors. The mesh distribution across processors creates consecutive and contiguous slices that render an optimized data access pattern for parallel I/O. A file format using the coherent mesh representation, developed and tested with OpenFOAM, enables the usability of the software at unprecedented scales. Further implications of the coherent and sliceable mesh representation arise due to simplifications in partitioning and diminished pre- and post-processing overheads.

The performance of OpenFOAM in strong scaling tests on HPC systems with AMD EPYC processors exhibits a pronounced superlinear speedup. Simple test cases show superlinear speedups of over 300%, which significantly impacts the efficient use of computing resources. With the last generation of HPC architectures, a superlinear speedup of about 10% to 20% was well expected and accepted by CFD users. The measured superlinear speedup is much more pronounced and withstands the communication overhead to even larger scales. A detailed performance analysis of OpenFOAM follows, employing various High-Performance Computing (HPC) architectures, including AMD, ARM and Intel systems. The performance metric FVOPS (Finite VOlumes solved Per Second) is introduced to compare the performance of Computational Fluid Dynamics (CFD) applications when varying the grid size, as occurs in a strong scaling test. The achievable FVOPS depends on various factors, including the simulation type, boundary conditions, and especially the grid size of a use case. Analysing FVOPS on a single node level with varying grid size shows a significant difference in performance and cache utilization, which explains the large superlinear speedups seen in the strong scaling tests. Furthermore, FVOPS can be used as a simple benchmark to determine the optimal number of grid elements per rank to simulate a given use case at peak efficiency on a given platform, resulting in time, energy, and cost savings. The FVOPS metric also facilitates the direct comparison between different HPC architectures. The tests using AMD, ARM, and Intel processors show a peak in performance when employing around 10,000 grid elements per core. The presence of a large L3 cache on AMD processors is particularly advantageous, as indicated by L3 cache miss rates observed on AMD EPYC processors. Our results suggest that future HPC architectures with larger caches and higher memory bandwidth would benefit the CFD community.

OpenFOAM is a CFD software widely used in both industry and academia. The exaFOAM project aims at enhancing the HPC scalability of OpenFOAM, while identifying its current bottlenecks and proposing ways to overcome them. For the assessment of the software components and the code profiling during the code development, lightweight but significant benchmarks should be used. The answer was to develop microbenchmarks, with a small memory footprint and short runtime. The name microbenchmark does not mean that they have been prepared to be the smallest possible test cases, as they have been developed to fit in a compute node, which usually has dozens of compute cores. The microbenchmarks cover a broad band of applications: incompressible and compressible flow, combustion, viscoelastic flow and adjoint optimisation. All benchmarks are part of the OpenFOAM HPC Technical Committee repository and are fully accessible. The performance using HPC systems with Intel and AMD processors (x86_64 architecture) and Arm processors (aarch64 architecture) have been benchmarked. For the workloads in this study, the AMD processor seems particularly suited resulting in an overall shorter time-to-solution.

This paper presents a new mesh data layout for parallel I/O of linear unstructured polyhedral meshes. The new mesh representation infers coherence across entities of different topological dimensions, i.e. grid cells, faces, and points. The coherence due to cell-to-face and face-to-point connectivities of the mesh is formulated as a tree data structure distributed across processors. The mesh distribution across processors creates consecutive and contiguous slices that render an optimized data access pattern for parallel I/O. A file format using the coherent mesh representation, developed and tested with OpenFOAM, enables the usability of the software at unprecedented scales. Further implications of the coherent and sliceable mesh representation arise due to simplifications in partitioning and diminished pre- and post-processing overheads.

The calculation of relative free-energy differences between different compounds plays an important role in drug design to identify potent binders for a given protein target. Most rigorous methods based on molecular dynamics simulations estimate the free-energy difference between pairs of ligands. Thus, the comparison of multiple ligands requires the construction of a “state graph”, in which the compounds are connected by alchemical transformations. The computational cost can be optimized by reducing the state graph to a minimal set of transformations. However, this may require individual adaptation of the sampling strategy if a transformation process does not converge in a given simulation time. In contrast, path-free methods like replica-exchange enveloping distribution sampling (RE-EDS) allow the sampling of multiple states within a single simulation without the pre-definition of alchemical transition paths. To optimize sampling and convergence, a set of RE-EDS parameters needs to be estimated in a pre-processing step. Here, we present an automated procedure for this step that determines all required parameters, improving the robustness and ease of use of the methodology. To illustrate the performance, the relative binding free energies are calculated for a series of checkpoint kinase 1 inhibitors containing challenging transformations in ring size, opening/closing, and extension, which reflect changes observed in scaffold hopping. The simulation of such transformations with RE-EDS can be conducted with conventional force fields and, in particular, without soft bond-stretching terms.

Secreted proteins of eukaryotes are decorated with branched carbohydrate oligomers called glycans. This fact is only starting to be considered for in silico investigations of protein dynamics. Using all-atom molecular dynamics (MD) simulations and Markov state modeling (MSM), we unveil the influence of glycans on the conformational flexibility of the multidomain protein disulfide isomerase (PDI), which is a ubiquitous chaperone in the endoplasmic reticulum (ER). Yeast PDI (yPDI) from Saccharomyces cerevisiae is glycosylated at asparagine side chains and the knowledge of its five modified sites enables a realistic computational modeling. We compare simulations of glycosylated and unglycosylated yPDI and find that the presence of glycan–glycan and glycan–protein interactions influences the flexibility of PDI in different ways. For example, glycosylation reduces interdomain interactions, shifting the conformational ensemble toward more open, extended structures. In addition, we compare our results on yPDI with structural information of homologous proteins such as human PDI (hPDI), which is natively unglycosylated. Interestingly, hPDI lacks a surface recess that is present in yPDI. We find that glycosylation of yPDI facilitates its catalytic site to reach close to this surface recess. Hence, this might point to a possible functional relevance of glycosylation in yeast to act on substrates, while glycosylation seems redundant for the human homologous protein. We conclude that glycosylation is fundamental for protein dynamics, making it a necessity for a truthful representation of the flexibility and function in in silico studies of glycoproteins.

A hallmark of N-linked glycosylation in the secretory compartments of eukaryotic cells is the sequential remodeling of an initially uniform oligosaccharide to a site-specific, heterogeneous ensemble of glycostructures on mature proteins. To understand site-specific processing, we used protein disulfide isomerase (PDI), a model protein with five glycosylation sites, for molecular dynamics (MD) simulations and compared the result to a biochemical in vitro analysis with four different glycan processing enzymes. As predicted by an analysis of the accessibility of the N-glycans for their processing enzymes derived from the MD simulations, N-glycans at different glycosylation sites showed different kinetic properties for the processing enzymes. In addition, altering the tertiary structure of the glycoprotein PDI affected its N-glycan remodeling in a site-specific way. We propose that the observed differential N-glycan reactivities depend on the surrounding protein tertiary structure and lead to different glycan structures in the same protein through kinetically controlled processing pathways.

The combination of Markov state modeling (MSM) and molecular dynamics (MD) simulations has been shown in recent years to be a valuable approach to unravel the slow processes of molecular systems with increasing complexity. While the algorithms for intermediate steps in the MSM workflow such as featurization and dimensionality reduction have been specifically adapted to MD datasets, conventional clustering methods are generally applied to the discretization step. This work adds to recent efforts to develop specialized density-based clustering algorithms for the Boltzmann-weighted data from MD simulations. We introduce the volume-scaled common nearest neighbor (vs-CNN) clustering that is an adapted version of the common nearest neighbor (CNN) algorithm. A major advantage of the proposed algorithm is that the introduced density-based criterion directly links to a free-energy notion via Boltzmann inversion. Such a free-energy perspective allows a straightforward hierarchical scheme to identify conformational clusters at different levels of a generally rugged free-energy landscape of complex molecular systems.

Thermally driven processes of molecular systems include transitions of energy barriers on the microsecond timescales and higher. Sufficient sampling of such processes with molecular dynamics simulations is challenging and often requires accelerating slow transitions using external biasing potentials. Different dynamic reweighting algorithms have been proposed in the past few years to recover the unbiased kinetics from biased systems. However, it remains an open question if and how these dynamic reweighting approaches are connected. In this work, we establish the link between the two main reweighting types, i.e., path-based and energy-based reweighting. We derive a path-based correction factor for the energy-based dynamic histogram analysis method, thus connecting the previously separate reweighting types. We show that the correction factor can be used to combine the advantages of path-based and energy-based reweighting algorithms: it is integrator independent, more robust, and at the same time able to reweight time-dependent biases. We can furthermore demonstrate the relationship between two independently derived path-based reweighting algorithms. Our theoretical findings are verified on a one-dimensional four-well system. By connecting different dynamic reweighting algorithms, this work helps to clarify the strengths and limitations of the different methods and enables a more robust usage of the combined types.