Polymer electrolyte fuel cells (PEFCs) have attracted special attention as a source of power for vehicles and stationary power generation systems, owing to features such as high energy density and fast refueling times. They are, however, still expensive and therefore their successful commercialization is contingent on continued reduction in their price. During the last decade according to U.S. Department of Energy, the cost power generated by fuel cell (FC) systems has been reduced from $124/kW in 2006 to $53/kW in 2015, with the ultimate goal set at $30/kW. Apart from reducing the manufacturing and material processing costs, maximizing their catalytic performance directly translates into cheaper FCs. This can be done in two ways: (1) coming up with new catalysts that are either cheaper or more durable (or both) than platinum, and (2) reducing the platinum loading without affecting performance. Here, we consider the latter by optimizing – via computer simulation – the internal structure of the PEFC catalyst layer (CL), which is host to platinum nanoparticles deposited on a carbon backbone for electron transport and a thin coating of an ionomer (e.g. Nafion) for ion transport. Optimizing the internal structure, by definition, requires that a pore-scale modeling strategy to be employed rather than volume-averaged continuum methods. This challenge has two parts: (1) digital reconstruction of the CL with high fidelity such that it is representative of an actual CL, and (2) multiphysics simulation of coupled nonlinear equations associated with different transport mechanisms in different locations, i.e. electron transport in the solid phase, ion transport in the ionomer phase, and gas (oxygen) transport in the void phase together coupled with interfacial phenomena such as gas partitioning and electrochemical reaction at void-ionomer and ionomer-platinum interfaces, respectively. The length scales at which different transport mechanisms occur vary from only a few nanometers through the ionomer film up to a few tenths of a micron in inter-agglomerate pores. Such disparate length scales make it difficult to perform – in reasonable time – direct numerical simulation (DNS) on a large-enough section of the entire CL such that it is statistically representative. We propose pore-scale strategy based on a hybrid of DNS and pore network modeling (PNM). PNM is a pore-scale modeling technique in which a porous domain is mapped onto an equivalent graph such that each node represents a single pore. Such simplification reduces the computational complexity of PNMs up to 4-5 orders of magnitude, compared to an equivalent DNS. In this work, the void and the carbon backbone are each modeled by PNM, and the ionomer thin film and platinum nanoparticles are modeled by DNS. Briefly, the justification behind this decision is that transport of oxygen and protons, for the most part, are restricted by the diffusivity and ionic conductivity in the ionomer phase rather than in the void and carbon, respectively. For the digital reconstruction of the CL, we use “process-based” reconstruction in which the actual manufacturing process of the CL is mimicked. The reconstructed geometry is in the form of a three-dimensional (3D) image and consists of 4 distinct phases: carbon support, platinum, Nafion (as the ionomer), and void. Using a network extraction algorithm, the equivalent networks corresponding to the void and carbon phases are extracted. The entire digital reconstruction was performed using the open-source software OpenPNM and PoreSpy. We propose a simple algorithm to couple the two extracted networks with the remainder of the original 3D image, which now only consists of Nafion and platinum. Note that the computational bottleneck of the entire simulation is the DNS, which is only applied to the Nafion and platinum phases. Since the carbon backbone and the void usually contribute up to 80% of the total volume of the CL, the proposed hybrid scheme significantly reduces the computational time compared with the case where the entire domain is modeled using DNS. While the proposed hybrid model allows for simulating much larger sections of the CL – therefore obtaining more representative results – the DNS step still imposes a serious limitation. Ideally, one could apply PNM to model the Nafion film, eliminating the DNS step altogether. However, since the thickness of the Nafion film is very small compared to the largest length scale of the system, network extraction algorithms require ultrahigh-resolution 3D images to generate a reasonably-accurate network that represents the thin Nafion film. Therefore, this study can be regarded as the first step toward developing a full network model of the CL. We hope that this study will become a cornerstone for many future studies.