## Abstract

Systematic enumeration and identification of unique 3D spatial topologies (STs) of complex engineering systems (such as automotive cooling systems, electric power trains, satellites, and aero-engines) are essential to navigation of these expansive design spaces with the goal of identifying new spatial configurations that can satisfy challenging system requirements. However, efficient navigation through discrete 3D ST options is a very challenging problem due to its combinatorial nature and can quickly exceed human cognitive abilities at even moderate complexity levels. This article presents a new, efficient, and scalable design framework that leverages mathematical spatial graph theory to represent, enumerate, and identify distinctive 3D topological classes for a generic 3D engineering system, given its system architecture (SA)—its components and their interconnections. First, spatial graph diagrams (SGDs) are generated for a given SA from zero to a specified maximum number of interconnect crossings. Then, corresponding Yamada polynomials for all the planar SGDs are generated. SGDs are categorized into topological classes, each of which shares a unique Yamada polynomial. Finally, within each topological class, 3D geometric models are generated using the SGDs having different numbers of interconnect crossings. Selected case studies are presented to illustrate the different features of our proposed framework, including an industrial engineering design application: ST enumeration of a 3D automotive fuel cell cooling system (AFCS). Design guidelines are also provided for practicing engineers to aid the application of this framework to different types of real-world problems such as configuration design and spatial packaging optimization.

## 1 Introduction

Engineering systems design [2–4] often involves choosing the most suitable candidate among many alternative design solutions to meet specific system performance criteria and design constraints using techniques such as comparative design analysis and optimization [5–8], and configuration design [9–18] methods. In most engineering design efforts, the component technologies and the component-to-component connectivity (referred to here as *system architecture* (SA)) remain fixed to preserve the functionality of the system, while different feasible system spatial layouts are explored. More precisely, SA specifies what components comprise a system, and how ports on components are connected to specific ports on other components. Distinct SAs represent specific technology options to perform each desired function and define the flow paths of material, energy, and/or information from one component to another. System architecture enumeration problems have been studied previously for engineering design examples such as hybrid electric power trains [12,14,19–21], automotive vehicle suspension systems [22], design of structures [23], plate heat exchangers [24], planar robot configurations [25], and optimal cooling system layouts [26,27] for dynamic thermal management. However, for each system architecture, many spatial topologies (STs) may exist, each with its own optimal geometry.

Any 3D engineering system geometric optimization problem can be simplified by decomposing the problem into two stages: (1) identification of unique spatial topologies and then (2) geometrically optimizing each system configuration, accounting for physical interactions and subject to all relevant geometric, failure mode, and other constraints. For example, if an interconnect links ports *P*1 and *P*2, many options exist for how this interconnect passes around various other interconnects and components in the system. Two spatial topologies are *equivalent* when there is a continuous deformation of component locations and interconnect trajectories that takes one topology to the other. This continuous deformation must be possible without cutting and rejoining any interconnects. Figure 1(a) shows 3D systems A and B1 having two different system architectures because interconnect IC1 is connected between components {1, 2} in A but between {1, 3} in B1. In other words, component-to-component connectivity is different in A and B1, respectively. Systems B1 and B2 have the same system architecture as all the component-to-component interconnections are the same. However, B1 and B2 have different spatial topologies because the interconnect IC2 between the component 1 and 2 is topologically different (please see the crossing patterns in Fig. 1(a)). As both the ends of the interconnect IC2 are fixed, it cannot be continuously morphed between B1 and B2. Hence, B1 and B2 are topologically different systems. Similarly, Fig. 1(b) is another example. The scope of this article encompasses the challenge of enumerating and identifying such unique spatial topologies for each system architecture within a design problem. Other recent articles have focused on continuous optimization of systems once spatial topology is specified (see Refs. [28–31], for example). The example shown in Fig. 1 is kept simple for illustration purposes, but the framework can be used to generate STs for more complex architectures and larger interconnected systems with multiple crossings.

### 1.1 Need for Efficient Spatial Topology Exploration Methods: To Address Fundamental Engineering Design Research Questions.

SA design automation methods studied in Refs. [19,22,26,32,33] have largely ignored spatial topology decisions. For instance, Refs. [19,32] focused on simultaneous configuration and sizing optimization of automotive powertrains but not on the spatial topological invariance of the connections between the components. From an engineering design point of view, there is only a little or practically no emphasis on leveraging the “spatial topological” aspects of these 3D systems during the design exploration or optimization process. How are the components and interconnects oriented or interacting among themselves *in real 3D space*? How are they located in that 3D space relative to each other? For example, do two interconnects cross each other and how are their physical fields affecting components or other interconnects in their proximity? What engineering design representations would be suitable to capture the spatial interactions (or the spatial topology) between the different system elements? There are several fundamental research questions that need to be carefully and thoroughly addressed. In other words, studying spatial dimensions and spatial interactions is very important as it directly has practical implications on the amount of 3D volume occupied by the engineering system, the accessibility of its components and how each of the subsystems in turn get influenced by their multiphysics spatial fields such as heat losses, electromagnetic radiation, and thermal stress or mechanical fatigue acting within the system.

For example, there exists a lot of wiring/piping, and components such as a battery, electric motor, and power converter in a hybrid electric vehicle (HEV) system. The spatial locations and topologies of these components and interconnects contribute to the overall spatial packaging density (SPD) of the HEV. SPD directly influences the vehicle’s efficiency, road range, and total energy consumption. Two HEV system models might have the same system architecture and component connectivity; however, their performance will be completely different because of the spatial locations and arrangements of its components and interconnect network. If lesser volume is occupied, heat loss radiation, noise and vibration, and cooling requirements should be addressed more aggressively. Furthermore, compactness makes it challenging to address accessibility and maintenance issues. A maintenance engineer would be glad if there is lot of working space available between components and subsystems with lesser entangled wiring for performing efficient overhaul, repair, or replacement operations. But the customer and product designer may want a compact vehicle owing to better efficiency, road range, and energy consumption. Therefore, there are key design tradeoffs that exist between different vehicle performance attributes such as overall system efficiency, packaging density, and accessibility due to spatial topological decisions. Similarly, for designing automotive climate and thermal heating, ventilation, and air-conditioning (HVAC) systems, sometimes the hot and cold refrigerant pipelines must not be close to each other in terms of physical proximity to avoid heat radiation effects. In such cases, spatial topological exploration is very helpful. Would it be possible to find other configurations where the pipelines are far apart while the system connectivity remains the same? It might be possible to take such decisions for one or two components or pipes. But for a large vehicle with hundreds of system elements, it is a humanly impossible task to come up with an optimal arrangement that takes into consideration different design choices. Hence, from an engineering design perspective, spatial topological decisions are vital, and there is a need to efficiently navigate through the discrete 3D design space to achieve optimal system configurations that balance between spatial packaging density, physics performance, and design for accessibility and maintenance.

As described earlier, there are many pain points in real industry practice due to 3D spatial constraints, and such design challenges must be thoroughly addressed to develop efficient, sustainable, and easily maintainable systems. This article aims to arouse in-depth discussion among engineering design researchers and help the larger design community realize the need for new representations and design automation frameworks to address such complex design problems and make relatively faster and better decisions. In summary, this article is trying to address some of the following fundamental research questions as listed below:

Does there exist a design representation that can capture 3D spatial topology of engineering systems? There has been a lot of work done in system architecture exploration. However, we need a design representation that can capture the “spatial” relationships between components and interconnections within the 3D system.

For systems having components immersed in multiphysics spatial fields, there might be different operating conditions that need to be met. For example, hot and cold fluid lines may need to stay apart, and components may be separated by bays or walls to avoid electromagnetic radiation or noise/vibration issues.

Many real-world engineering systems contain hollow objects. Design decisions might include whether a pipe or duct should pass through or outside the hole.

Can new algorithms be created that can navigate 3D discrete space efficiently to explore possible design candidates that satisfy operational and maintenance constraints?

For a given system, how many spatial topological solutions exit for multi-pipe-routing or wiring?

Are there mathematical ways to map abstract designs into polynomials for faster evaluation and comparison with filter nonunique spatial topologies?

What would be running time complexity for such design automation algorithms or frameworks?

Would such representations support simultaneous geometry, physics, and topology optimization? Can gradient-based optimization methods be used?

Would such a spatial topology enumeration framework help maintenance engineers determine what is the optimal trajectory that components can take to be easily removed or replaced from the system. What is the system configuration that has maximum access ports for critical components?

What is best spatial configuration of piping or wiring to support tight packaging constraints?

Can these generic algorithms be applied to systems at various scales?

Can artificial intelligence and machine learning techniques be employed to make these methods efficient and scalable?

The aforementioned questions are generic and fundamental and applicable to any interconnected 3D engineering system. Some common design engineering questions have been outlined, and this is not considered to be an exhaustive list. However, the decisions involved in each question can have great impact on overall system performance, early design stage process, and several life-cycle consideration. Interesting case studies have been demonstrated in Sec. 5 as an attempt to answer some of these aforementioned research questions using the proposed enumeration framework. The inferences obtained from the detailed case studies have been presented in Sec. 6. However, it must be noted that this article is a preliminary attempt to answer these challenging questions, and answering all questions is definitely out of scope of this article. Further investigation is being done on these topics, and a lot of support is required from the larger design engineering community to make impact.

### 1.2 Design Process and the Power of Design Representations.

Any systematic design process [34,35] involves four key tasks: *representation*, *generation*, *evaluation*, and *design guidance*. *Representation* refers to the task of describing a system using a generic model that captures the functionality of the various system elements. Depending on the design analysis tools and application requirements, design representations can be mathematical, graphical, physics-based, or conceptual [36,37]. A formal representation of engineering systems is necessary for generation of their different spatial configurations. Previous works utilized different representations such as stick diagrams [38], canonical graphs [39,40], lever diagrams [41,42], design structure matrices [43], and bond graphs [19,44–48]. Each of these representations has advantages and disadvantages that impact the computational effort required to generate the feasible design space of engineering system configurations. The *generation* phase involves creating feasible design alternatives using the representation based on design rules. Configuration synthesis and automatic enumeration strategies [14,19] that meet design constraints have been developed previously by researchers for applications such as planetary gear sets for automotive hybrid transmission systems [49], split hybrid configurations for single planetary gears [50], robotic manipulator configurations [51], and vehicle suspension mechanisms [52]. *Evaluation* is the process of measuring the design quality in terms of the performance criteria. Finally, *design guidance* is providing feedback for the generation task based on the evaluation output to find better alternatives in the design space. In design automation methods, the generation, evaluation, and design guidance tasks are performed in an automated loop that, when successful, converges to a design solution. However, the element of this process that actually enables high design accuracy, comprehensive search space navigation, and computational efficiency is the design representation that is selected before the start of this iterative process. Hence, a design representation for 3D interconnected systems that captures the relevant problem attributes and aligns well with a generation computation method is critical for navigating through the discrete 3D spatial topology options efficiently.

The design elements of any three-dimensional interconnected engineering systems are its components, their 3D spatial locations and orientations, their port valencies, interconnections (diameters/sizes and trajectories), and the crossings between their interconnections. Unlike 2D system layout enumeration and optimization performed in the design of very large-scale integrated circuits [53–64], the 3D spatial layout enumeration problem is fundamentally different and an even more challenging problem as described in Refs. [31,65,66]. More specifically, even for a given system architecture, options exist for how the interconnects are routed relative to one another and to the components (e.g., say, should duct A go over cable B or under, or should a pipe be routed through the hole in the casing, or around the edge of the casing, etc.). These are topologically discrete design differences. To cater to that need, in this article, we have identified a viable mathematical design representation to describe 3D interconnected spatial layouts as spatial graph diagrams (SGDs).

### 1.3 Objectives and Contributions.

Complex engineering systems such as autonomous aerial vehicles [67], electric power trains [19,68], aero-jet engines [69], or vehicle thermal management and cooling systems [26,70] have different kinds of components connected together either through wires, ducts, or pipes entangled with one another in a tightly packed three-dimensional volume. Engineering design alternatives with distinct spatial topologies may exist, with different values of metrics such as efficiency, spatial packaging density, maintenance costs, and design complexity. Current practice for exploring different system spatial topologies relies largely upon human expertise, design rules, creating new designs derived from existing ones (resulting in incremental design advancements), and manual adjustments (limiting complexity of designable systems).

This approach bounds the pace and scale at which spatial topologies can be explored for practical application to typical complex systems. The realization of significantly enhanced functionality or performance is prevented by the profoundly incomplete design space coverage of current practice. More efficient automated methods are needed to unleash engineering system design efforts from the restraints of current design practice for spatial topological decisions. Enumeration methods in low-dimensional topology for knots are well known, which are a special class of spatial graphs, see Ref. [71] for an overview. By using additional techniques from hyperbolic geometry, it is possible to exactly enumerate all knot topologies with less than 20 crossings, of which there are more than 350 million [72]. Compared to such massive computations, the prior work on spatial graphs where the underlying graph is more complicated than just a loop is limited: mostly tabulations of less than 100 topologies [73–78]. For example, the authors in Ref. [76] generated two vertex bouquet spatial graphs with a maximum of seven crossings. We have taken inspiration from works in mathematical low-dimensional topology to create new algorithms that can generate 3D system spatial topologies using spatial graphs. As our cases studies in Sec. 5 demonstrate, the enumeration strategy in this article allows for much larger scale enumerations, with arbitrary specified system architectures. This is very suitable for representing large-scaled complex engineering systems easily, and for enumerating their 3D spatial topologies efficiently.

The main objectives of the work presented in this article are as follows: (1) to present fundamental engineering design research questions listed in Sec. 1.1, and (2) to develop an automated and systematic enumeration framework to both represent 3D engineering systems using SGDs and efficiently generate distinct STs for a given system architecture using a rigorous mathematical approach as an attempt to address the research questions outlined at the end of Sec. 1.1.

The major contributions of the proposed design framework include:

This is a new way to represent 3D engineering system configurations using spatial graph theory. This representation supports description of components as nodes, interconnects as edges, allows multiple interconnect crossings, and supports variable component valency.

Combinatorial enumeration of all SGDs for a given system architecture up a maximum number of interconnect crossings. The proposed method presented in this article supports enumeration of spatial topologies for moderately scaled engineering system architectures containing (approximately) up to 20 components, 60 interconnects, and 10 crossings. To enumerate spatial topologies for much larger scale systems, an effective decomposition-based method, such as the approach presented in case study in Sec. 5.5, can be utilized. As part of future work, we also plan to investigate alternative methods such as deep learning and pattern recognition as mentioned in Ref. [79] to efficiently explore newer topologies using large-scaled system datasets of 3D system architectures. This sort of strategy has been used previously to scale SA design to problems too large for enumeration using machine learning strategies in conjunction with enumerated data as articulated in Refs. [80–83]

Efficient and systematic identification of unique SGDs from the exhaustively enumerated SGD set using Yamada polynomial invariants. The Yamada polynomials are used as a tool to identify any duplicate spatial graph topologies, similar to isomorphisms for standard graphs. This serves as a foundation to explore the discrete 3D spatial topology design space thoroughly.

Topological equivalence between spatial graphs is tested here using Yamada polynomials rather than comparing graph diagrams directly.

Practical demonstration of the proposed design framework on a real-world industry application: enumeration of unique spatial topologies of a 3D automotive fuel cell cooling system (AFCS).

Engineering design insights and guidelines help system design engineers apply this proposed framework to different kinds of practical applications. An in-depth discussion is presented in Sec. 6.

The remainder of this article is organized as follows. The terminology and notation of the proposed spatial graph representation are introduced in Sec. 2. Section 3 describes the characteristics of Yamada polynomials and how they are evaluated for an individual SGD. In summary, Secs. 2 and 3 provide a detailed discussion introducing readers to the mathematical theory of spatial graphs, Yamada polynomial invariants, and graph relations through a few illustrative examples. Section 4 demonstrates the proposed six-step design framework that utilizes spatial graphs to represent, enumerate, and identify distinctive 3D topological classes for an engineering system, given its SA. SGDs are generated for a given SA from zero to a specified maximum crossing number. Corresponding Yamada polynomials for all the enumerated SGDs are then generated. SGDs are categorized into topological classes, each of which shares a unique Yamada polynomial. Finally, for each topological class, 3D geometric models are generated that can be used for multiphysics configuration design optimization. Section 5 presents several interesting practical case studies based on the proposed framework. The results are discussed in Sec. 6 from an engineering system design perspective. Finally, the conclusion and future work items are presented in Sec. 7.

## 2 Spatial Graphs

The study of graphs in 3-space has been mathematically formalized using *spatial graphs* [84–86], which we will now describe. Suppose *G* is a graph, that is, a set of vertices and a set of edges, where an edge is just a pair of vertices. (Edges are undirected and multiple edges between the same pair of vertices are allowed.) A *spatial embedding* of a graph *G* is a set of points (nodes) in **R**^{3} corresponding to the vertices of *G*, and a set of smooth arcs (links) corresponding to the edges of *G* that join appropriate pairs of vertices; here, each arc meets the vertices only at its two endpoints, and it intersects other arcs only at these vertices. Collectively, these points and arcs form a *spatial graph* with underlying (abstract) graph *G*. More formally, the spatial embedding is a function *f* : *G* → **R**^{3}, whose image $G~:=f(G)$ is the spatial graph. See Fig. 2(a) for a sample spatial graph. The natural topological notion of equivalence for spatial graphs is *isotopy*, when two spatial graphs $G~1$ and $G~2$ can be continuously deformed from one to the other without any arc passing through another arc or itself.

Spatial graphs are a natural extension of the knot theory, which is the study of circles embedded in **R**^{3}, since we can put vertices on a knot to make it into a spatial graph. While the study of knot theory has its origin in the physics of the late 19th century [71], the spatial graph theory has its roots in chemistry [87,88] and is different from the graph theory because graph theory studies abstract graphs, while the spatial graph theory studies embeddings of graphs in **R**^{3} or even in other three-manifolds [89–91]. This theory was used in polymer stereochemistry [87,92], ecology research [93], and molecular biology (e.g., protein classification and protein structure prediction [94–97] using deep learning) to distinguish different topological isomers. A folded protein can be thought of as a spatial graph where residues are the nodes and edges connect the residues in close proximity. This article’s goal is to improve and tailor the immense potential of these representations in developing new engineering design automation methods.

If a spatial graph is projected onto a plane, then some arcs (edges) may appear to cross in the projection. If information about which arc is on top at the apparent crossings is omitted, the projection is called a *shadow* of the spatial graph, as shown in Fig. 2(b). If we keep track of which arc is on top at each apparent crossing, the projection or planar representation is called a *diagram* of the spatial graph, as shown in Fig. 2(c). In other words, diagrams are the images of embedded graphs under a projection **R**^{3} → **R**^{2} whose singularities are a finite number of crossings of edges equipped with over-under crossing information. Hence, many different spatial diagrams of a spatial graph *G* may have the same shadow. We can produce this family of SGDs by assigning all possible permutations of overstrand or understrand information.

### 2.1 Reidemeister Moves.

Given an abstract graph *G*, we can use the spatial graph diagrams above to begin enumerating spatial embeddings of *G*. The challenge is then to determine which of these SGDs actually describe isotopic spatial embeddings (i.e., are topologically equivalent), so that later steps in the design process consider each topological possibility only once. Fortunately, it has been shown that two diagrams represent isotopic embeddings if and only if they are related by a finite sequence of fundamental *Reidemeister moves* (*R*0 to *R*6) [98–100] as shown in Fig. 3. Figure 4 shows a simple illustration of three diagrams where SGDs A and B are *topologically equivalent* under the first Reidemeister move *R*1, whereas C is not equivalent to either A or B as its edges cannot be continuously deformed using the Reidemeister moves to attain A or B, so they represent *topologically distinct* spatial graphs.

### 2.2 Flat Vertex Graphs and Ribbon Graphs.

The topological formulation of spatial graphs is quite idealized in that each vertex has no local structure and the edges are infinitely thin. We can impose additional, but still purely topological, structure by considering *flat vertex graphs (FVGs)* and *ribbon graphs*, which may be more suitable for certain design applications. A flat vertex graph is a spatial graph where the vertices correspond to flat disks in **R**^{3} as shown in Fig. 5. In particular, this gives the edges coming in to each disk a cyclic order. A flat vertex graph can also be encoded by a SGD, with the convention that each disk is rotated parallel to the projection plane before projecting. Two SGDs represent isotopic flat vertex graphs if and only if they differ by a series of Reidemeister moves *R*0 to *R*5; here, *R*6 is no longer allowed since it would change the order of the edges coming into the vertex disk.

A ribbon graph is a spatial graph whose vertices have become flat disks and whose edges have become thin bands, depicted in Fig. 6. These too can be encoded as SGDs by using the blackboard framing convention (a way to view a knot diagram on a plane) illustrated in Fig. 6. This framing is obtained by converting each component to a ribbon lying flat on the plane. Two SGDs represent isotopic ribbon graphs if and only if they differ by a sequence of Reidemeister moves *R*0 and *R*2 to *R*5. The basic notation of a spatial graph introduced in Sec. 2 is sometimes referred to as a *pliable* spatial graph to contrast the notion with flat vertex and ribbon graphs. Here, we will focus on pliable and flat vertex graphs, but note that ribbon graphs would be useful for measuring twisting along interconnects in the final 3D system.

## 3 Yamada Polynomial Invariants

Reidemeister moves are valuable for identifying when two embeddings are isotopic (i.e., topologically equivalent); however, finding the specific sequence of moves between two equivalent spatial diagrams can be extremely challenging, especially when the spatial graphs have many nodes and edges. Even for knots, which are the simplest class of spatial graphs, it is unknown whether there exists a polynomial-time algorithm for determining when two knots are isotopic. (It is not impossible that such an algorithm exists: the question of whether a knot is equivalent to a round circle is in $NP\u2229coNP$ [101,102].) To show that two embeddings are not *isotopic* requires an *invariant*: a function of the embeddings whose output is not changed by isotopies and which takes different values on the two embeddings [103–106]. Mathematicians often use such invariants that are computable and yet powerful enough to detect some delicate differences of embeddings of the same graph. Over the last century, many *polynomial invariants* [104,107–109] were discovered by knot theorists, such as the Alexander-Conway [110], Jones [111], Kauffmann [112], and Yoshinaga [113] polynomials. Some of these have been extended to spatial graph theory [114–116] using similar constructions. These invariants satisfy nice *skein relations* that are mathematical tools that give linear relationships between the polynomials of closely related diagrams. Relevant skein relations are sufficient to calculate the polynomials recursively and are relatively convenient to use for this purpose. The proof of invariance then relies on using the skein relation to show the value of the invariant is unchanged by Reidemeister moves.

### 3.1 Yamada Polynomial Properties.

The specific polynomial invariant used here is the Yamada polynomial, which associates to each SGD a polynomial in an indeterminate *A*, which is an arbitrary independent variable. For example, it turns out that the Yamada polynomial for the SGD C in Fig. 4 is −*A*^{−6} − *A*^{−5} − *A*^{−4} − *A*^{−3} − *A*^{−2} − 1 − *A*^{2} + *A*^{6}. The Yamada polynomial is defined in terms of a polynomial invariant *H* of ordinary (nonspatial) graphs. Continuing the example, the abstract theta graph that underlies all the SGDs in Fig. 4 has *H*(*A*) = −*A*^{−2} − *A*^{−1} − 2 − *A* − *A*^{2}. We first discuss the polynomial invariant *H* and its properties P1–5, which can be used to recursively compute *H* for any graph. We then turn to the Yamada polynomial and its fundamental properties S1–S4, which similarly allow for recursive computation of it for any SGD. After giving some helpful additional properties Y1–Y8 of the Yamada polynomial, it is computed for four different example SGDs in Sec. 3.2.

Let *G* = (*V*, *E*) be an abstract graph, where *V* is the vertex set and *E* is the edge set of *G*. For two graphs *G*_{1} and *G*_{2}, $G1\u2294G2$ denotes the disjoint union of *G*_{1} and *G*_{2}, and *G*_{1}∨*G*_{2} denotes a wedge at a vertex of *G*_{1} and *G*_{2}, that is, $G1\u2228G2=G1\u222aG2$ and $G1\u2229G2={vertex}$. In addition, a graph *G* has a cut-edge *e* (also known as bridge or isthmus) if *G* − *e* has more connected components than *G*. First, following Ref. [117], a polynomial invariant *H*(*G*)(*A*) of an abstract graph *G* is described, where *A* is an indeterminate (arbitrary independent variable); precisely, our *H*(*G*)(*A*) is Yamada’s *h*(*x*, *y*) with *x* = −1 and *y* = −*A* − 2 − *A*^{−1}. The polynomial *h*(*G*)(*A*) is characterized by the following properties:

P1. $H(emptygraph)=1$ and $H(simpleloop)=A+1+A\u22121$.

P2. $H(G1\u2294G2)=H(G1)H(G2).$

P3.

*H*(*G*_{1}∨*G*_{2}) = −*H*(*G*_{1})*H*(*G*_{2}).P4. If

*G*has a cut edge, then*H*(*G*) = 0.P5. Let

*e*be a nonloop edge of a graph*G*. Then,*H*(*G*) =*H*(*G*/*e*) +*H*(*G*−*e*), where*G*/*e*is the graph obtained from*G*by contracting*e*to a point and*G*−*e*is*G*with*e*deleted.

*Yamada polynomial*[116–120]. Let

*g*be the a spatial graph diagram. For a crossing

*c*of

*g*, three

*graph reductions*are defined as follows:

*s*

_{+},

*s*

_{−}, and

*s*

_{0}(denotes a vertex), with the class of spin +1, −1, and 0, respectively, as shown in Fig. 7. These graph reductions are used to replace crossings in a spatial graph for Yamada polynomial calculation. Let

*S*be the planar graph obtained from

*g*by replacing each crossing with a spin.

*S*is called a

*state*on

*g*and

*U*(

*g*) denotes the set of states on

*g*obtained by applying all possible reductions in its crossings. Set ${g|S}=An1\u2212n2$, where

*n*

_{1}and

*n*

_{2}are the numbers of crossings with spin of +1 and spin of −1, respectively, and

*A*is an indeterminate. The Yamada polynomial

*R*{

*g*}(

*A*) ∈

**Z**[

*A*,

*A*

^{−1}] is defined as follows:

*g*does not have crossings, then

*R*(

*g*) =

*H*(

*g*). This Yamada polynomial for a spatial graph can be computed recursively using the following skein relations and the properties of

*H*:

S1.

*R*() =*AR*() +*A*^{−1}*R*() +*R*(),S2.

*R*() =*R*() +*R*(), where*e*is a nonloop edge.S3. $R(g1\u2294g2)=R(g1)R(g2)$,

S4.

*R*(*g*_{1}∨*g*_{2}) = −*R*(*g*_{1})*R*(*g*_{2}).

So far, the Yamada polynomial is a function of the given diagram *g*, and we need an invariant of the spatial graph $G~$ it describes. Yamada showed:

(I1) Any two diagrams

*g*and*g*′ whose flat vertex graphs $G~$ and $G~\u2032$ are isotopic have*R*(*g*′) = (−*A*)^{n}*R*(*g*) for some integer*n*.(I2) If every vertex has valence at most three, then two diagrams

*g*and*g*′ whose spatial graphs $G~$ and $G~\u2032$ are isotopic have*R*(*g*′) = (−*A*)^{n}*R*(*g*) for some integer*n*.(I3) Any two diagrams

*g*and*g*′ whose associated ribbon graphs $G~$ and $G~\u2032$ are isotopic have*R*(*g*′) =*R*(*g*).

The next set of relations for the Yamada polynomial can be derived from the previous ones and are very useful aides for its calculation. Detailed proofs for these relations (Y1–Y8) are provided in Ref. [117]. They are as follows:

Y1.

*R*() =*B*, where*B*=*A*+ 1 +*A*^{−1},Y2.

*R*() = −*BR*(),Y3.

*R*() = −*AR*() − (*A*^{2}+*A*)*R*(),Y4.

*R*() = −*A*^{−1}*R*() − (*A*^{2}+*A*^{−1})*R*(),Y5.

*R*() = −*AR*(),*R*() = −*A*^{−1}*R*(),Y6.

*R*() =*A*^{2}*R*(),*R*() =*A*^{−2}*R*(),Y7. Edge subdivision does not change the polynomial:

Y8. Petals to concentric self-loops:

= −

*B*^{2}= − (*A*+ 1*A*^{−1})^{2}.

### 3.2 Illustrative Examples.

Yamada polynomials for a few spatial graphs are calculated by reducing the spatial graph diagram into a linear combination of smaller elements based on the skein relations stated earlier.

Theta (*θ*_{1}) graph: The Yamada polynomial for a standard theta graph is calculated as follows:

A spatial graph that is isotopic (by the R6 move) to the standard theta graph. Its Yamada polynomial is calculated as follows, though one could instead use property Y5 as a shortcut.

The spatial graph is . This example involves extensive use of Yamada polynomial skein relations (Y1–Y8). Its Yamada polynomial is calculated as follows:

The spatial graph is . This example involves extensive use of skein relations (Y1–Y8). Its Yamada polynomial is calculated as follows:

## 4 Proposed Spatial Graph Diagram Enumeration Framework

Figure 8 shows the steps of the proposed design framework to represent, enumerate, and categorize unique spatial topologies of a 3D system (assuming intercomponent connectivity is fixed). The detailed steps of the enumeration design framework are as follows:

Define system architecture: Provide the specific 3D SA for which spatial topologies must be enumerated. From the SA, extract the number of nodes (components), their valencies, system interconnectivity, and the corresponding edges (interconnects in the system).

Enumerate spatial graph diagrams: Combinatorially enumerate all possible SGDs for the SA from zero crossings to the maximum crossing number (

*k*= 0, 1, …,*k*_{m}) using their corresponding shadows.Check graph planarity: Planar diagrams of spatial graphs are used for the calculation of the Yamada polynomials. However, before calculation, each enumerated graph must be checked to determine whether it is planar, for which there are linear-time algorithms [121]. The graphs start with a circular order of the edges at each vertex, making the planarity check (PC) even easier. The algorithm shown in Fig. 9 recursively contracts the edges of a graph until the diagram is a bouquet of circles, and then use the fact that if the diagram is planar, there must exist at least one loop edge whose endpoints come consecutively in the cyclic ordering around the vertex, i.e., a self-loop. This self-loop can be removed without altering planarity. The recursive steps for the PC algorithm are enumerated as follows:

PC1. Convert all vertices to crossings, as it does not affect planarity.

PC2. Contract all nonloop edges (edges shared between two vertices) to a vertex, as it does not affect planarity.

PC3. Remove all planar self-loops at a vertex.

PC4. Empty vertex does not affect planarity, so remove it. By doing all these steps recursively, if the result is an empty diagram, then the original diagram is planar. If not, the diagram is nonplanar.

Evaluate Yamada polynomials: The Yamada polynomials for all the valid planar SGDs are evaluated using the Yamada polynomial properties detailed in Sec. 3.

Categorize different spatial topologies: Cluster SGDs into topological classes so that the SGDs belonging to each class have the same Yamada polynomial or differ by a factor (−

*A*)^{n}based on property I2. Consequently, no two classes have the same Yamada polynomial. 3D geometric models can then be generated for SGDs under different spatial topology classes as required.3D model generation: SGDs from step 5 are utilized as underlying skeleton structures for generating various 3D system geometric models. As discussed in step 5, any pair of SGDs that represent isotopic spatial graphs will be in the same Yamada class. 3D geometric models for SGDs within the same topology class can be generated based on two different design requirements as follows:

In some engineering applications, additional crossings beyond what is required between a pair of interconnects may be undesirable. For example, there may be no advantage to two wires intertwining several times unless they intentionally function as a twisted pair. Similarly, crossings between pair of pipes may be difficult to fabricate and more complex to install and remove. Having fewer crossings helps reduce complexity and helps in easier system diagnosis and maintenance. In such cases, an SGD having the fewest crossings is selected to generate a 3D geometric model.

However, in some applications, this assumption that the simplest SGD is desirable may not always be applicable. Sometimes due to interface requirements, interconnects may require to exhibit complex topology (e.g., pipes or electrical lines passing through a fire-rated interface in a wall). Other examples include cases where multiple crossings enable easier assembly, disassembly, design effective wire harnesses, or simpler instantiating of other SGD elements. Quite often a collection of intertwining wires or ducts can help increase mechanical strength to provide support and conserve space.

## 5 Case Studies

In this section, a number of case studies are provided to demonstrate the proposed enumeration framework discussed in Sec. 4. All computations (Yamada polynomial calculation, planarity check, etc.) in the case studies were performed using wolfram Mathematica 11.3 software with an Intel Xeon E5-2660 CPU @ 2.00 GHz, 64 GB DDR4-2400 RAM, Windows 10 64-bit workstation.

### 5.1 Case Study 1: Components With Equal Valencies.

In this case study, we consider a 3D system with architecture as shown in Fig. 10. This system contains four identical trivalent components (nodes) and six interconnects (edges). We find the unique 3D STs of the system for crossing numbers varying from zero to three. The notation used to indicate each SGD is given by SGD_*k*, where *k* is the crossing number of that diagram, and the letter refers to the specific SGD. SGD_0 is the original system architecture without any crossings. By using the proposed framework, we combinatorially enumerate all the SGDs and pass them through a planarity check procedure. Yamada polynomials were then calculated for 3, 31, 118, and 231 valid planar SGDs having 0, 1, 2, and 3 crossing numbers, respectively. We group the SGDs having the same Yamada polynomial or differing by a factor (− *A*)^{n} under the same topological class based on property I2. Through this, we attain a total of four unique Yamada classes as shown in Fig. 10. The Yamada polynomials for these different classes of topologies are shown in Table 1. For class 1, it is observed that three distinct Yamada polynomials exist for different crossing numbers, but they all differ by a factor (− *A*)^{n}. This strongly suggests that the SGDs shown under class 1 are isotopic to SGD_0, and indeed, this can be verified using the R6 move. In contrast, the two SGDs in class 3 can be shown to be nonisotopic using other tools. A sample 3D model of class 2 spatial topology candidate design solution is also shown in Fig. 10. The computational time for this entire case study is 78.3 s.

Classes | Yamada polynomials |
---|---|

Class 1 SGD_0: | A^{−3} + A^{−2} + 3A^{−1} + 2 + 3A + A^{2} + A^{3} |

Class 1 SGD_1: | ( − A)^{1}(A^{−3} + A^{−2} + 3A^{−1} + 2 + 3A + A^{2} + A^{3}) |

Class 1 SGD_2: | ( − A)^{2}(A^{−3} + A^{−2} + 3A^{−1} + 2 + 3A + A^{2} + A^{3}) |

Class 1 SGD_3: | ( − A)^{3}(A^{−3} + A^{−2} + 3A^{−1} + 2 + 3A + A^{2} + A^{3}) |

Class 2 | A^{−4} + A^{−3} + A^{−2} + 2A^{−2} + 2A + 2A^{3} + A^{4} + A^{5} + A^{6} |

Class 3 | A^{−7} + A^{−6} + 3A^{−5} + 3A^{−4} + 4A^{−3} + 3A^{−2} + 3A^{−1} − 3A^{2} − 2A^{3} − 2A^{4} − A^{5} + A^{6} + A^{8} |

Class 4 | A^{−7} + A^{−6} + A^{−5} + 2A^{−4} + 2A^{−3} + A^{−2} + 3A^{−1} + 2A − A^{2} − A^{5} + A^{6} |

Classes | Yamada polynomials |
---|---|

Class 1 SGD_0: | A^{−3} + A^{−2} + 3A^{−1} + 2 + 3A + A^{2} + A^{3} |

Class 1 SGD_1: | ( − A)^{1}(A^{−3} + A^{−2} + 3A^{−1} + 2 + 3A + A^{2} + A^{3}) |

Class 1 SGD_2: | ( − A)^{2}(A^{−3} + A^{−2} + 3A^{−1} + 2 + 3A + A^{2} + A^{3}) |

Class 1 SGD_3: | ( − A)^{3}(A^{−3} + A^{−2} + 3A^{−1} + 2 + 3A + A^{2} + A^{3}) |

Class 2 | A^{−4} + A^{−3} + A^{−2} + 2A^{−2} + 2A + 2A^{3} + A^{4} + A^{5} + A^{6} |

Class 3 | A^{−7} + A^{−6} + 3A^{−5} + 3A^{−4} + 4A^{−3} + 3A^{−2} + 3A^{−1} − 3A^{2} − 2A^{3} − 2A^{4} − A^{5} + A^{6} + A^{8} |

Class 4 | A^{−7} + A^{−6} + A^{−5} + 2A^{−4} + 2A^{−3} + A^{−2} + 3A^{−1} + 2A − A^{2} − A^{5} + A^{6} |

### 5.2 Case Study 2: Enumeration of Spatial Topologies of an Automotive Fuel Cell System.

There are two main goals of this case study: first, to demonstrate the proposed spatial topology enumeration framework on a practical industry application, the Ford AFCS. Second, to provide insight into how the number of unique spatial topologies varies with increasing the number of components and the maximum number of crossings in a 3D system.

*Automotive Fuel Cell System Description:* An AFC system contains fuel cell stacks; the required cooling system components such as heat exchangers, pumps, radiators, cooling fans, particle filters, and cabin heat generation unit; and a cooling hose interconnect network as shown in Fig. 11(a). Underhood 3D spatial packaging of the AFCS is essential for efficient thermal management in fuel cell electric vehicles. However, this requires finding a suitable AFCS 3D spatial topology that can be optimized to minimize underhood volume while delivering the required vehicle capability and physics-based system performance. Unique AFCS spatial topologies can be obtained using the proposed enumeration framework to thoroughly navigate through the discrete 3D design space for a given AFCS architecture. These unique spatial configurations can be utilized as initial 3D designs to perform parametric gradient-based multiphysics optimization to reduce overall AFC system volume subject to geometric and physics-based constraints. Finally, an AFCS spatial topology configuration is selected that satisfies all the performance criteria and constraints. However, this optimization process is beyond the scope of this article and has been separately addressed in detail in our previous articles [28–30,122–124]. In this study, we are only interested in demonstrating an important use case scenario of this proposed enumeration framework.

Here, we consider only AFCS system architectures where all nodes (components) are trivalent. The original AFCS system shown in Fig. 11(a) contains 14 trivalent components. To test our framework on *highly connected* AFCS architectures containing different number of components, we group together some of the closer or smaller components of the 14-component 3D AFCS layout to create architectures having 6, 8, 10, and 12 components. For example, Fig. 11(b) shows an AFCS architecture containing only 10 trivalent components where the 2 WEG (Werner, Eggon, and Gerald company) heater units and heater core are combined to a single larger component, the small particle filter is merged into the full cell stack, and the 2 cooling fans are combined with the radiators, thus reducing a 14-component system to a 10-component system architecture. Then, unique spatial topologies are enumerated for such AFCS system architectures containing 6, 8, 10, 12, and 14 trivalent components with interconnect crossing numbers varying between 0 and 6. Table 4 shows the number of topologically distinct AFCS spatial configurations (with distinct Yamada polynomials) generated by this study. Table 5 shows the Yamada polynomials for the 6 different classes of spatial graph ttoplogies illustrated in Fig. 12.

*Running Times:* Tables 2 and 3 show the enumeration running times in seconds, and in log-base-10, respectively, for different system architectures with varying number of nodes and crossings in the spatial graph diagrams. Computations were performed in parallel on a 10-core Xeon workstation using code written in python. Note that increasing either the number of nodes by 2 or the number of crossings by 1 increases the running time by a factor of 10.

Crossings | |||||||
---|---|---|---|---|---|---|---|

Components | 0 | 1 | 2 | 3 | 4 | 5 | 6 |

6 | 0.03 | 0.06 | 0.11 | 1.03 | 9.46 | 98.44 | 1103.10 |

8 | 0.04 | 0.09 | 0.89 | 14.35 | 224.73 | 3087.12 | |

10 | 0.05 | 0.37 | 9.26 | 213.36 | 4202.83 | ||

12 | 0.10 | 2.37 | 98.08 | 2879.99 | |||

14 | 0.40 | 25.17 | 1258.59 | 45363.99 |

Crossings | |||||||
---|---|---|---|---|---|---|---|

Components | 0 | 1 | 2 | 3 | 4 | 5 | 6 |

6 | 0.03 | 0.06 | 0.11 | 1.03 | 9.46 | 98.44 | 1103.10 |

8 | 0.04 | 0.09 | 0.89 | 14.35 | 224.73 | 3087.12 | |

10 | 0.05 | 0.37 | 9.26 | 213.36 | 4202.83 | ||

12 | 0.10 | 2.37 | 98.08 | 2879.99 | |||

14 | 0.40 | 25.17 | 1258.59 | 45363.99 |

Crossings | |||||||
---|---|---|---|---|---|---|---|

Components | 0 | 1 | 2 | 3 | 4 | 5 | 6 |

6 | −1.47 | −1.25 | −0.96 | 0.01 | 0.98 | 1.99 | 3.04 |

8 | −1.43 | −1.06 | −0.05 | 1.16 | 2.35 | 3.49 | |

10 | −1.30 | −0.44 | 0.97 | 2.33 | 3.62 | ||

12 | −1.01 | 0.37 | 1.99 | 3.46 | |||

14 | −0.39 | 1.40 | 3.10 | 4.66 |

Crossings | |||||||
---|---|---|---|---|---|---|---|

Components | 0 | 1 | 2 | 3 | 4 | 5 | 6 |

6 | −1.47 | −1.25 | −0.96 | 0.01 | 0.98 | 1.99 | 3.04 |

8 | −1.43 | −1.06 | −0.05 | 1.16 | 2.35 | 3.49 | |

10 | −1.30 | −0.44 | 0.97 | 2.33 | 3.62 | ||

12 | −1.01 | 0.37 | 1.99 | 3.46 | |||

14 | −0.39 | 1.40 | 3.10 | 4.66 |

Figure 12 shows a sampling (6 out of 37) of the 3D geometric models that have been automatically generated for a 10-component system with 3 crossings. Similarly, Fig. 13 shows 3 out of 283 spatial topologies for a 10-component system with 4 interconnect crossings. The top and orthogonal views along with their corresponding spatial graphs were generated using autodesk fusion 360 software. The top view helps in visualizing the component–component connectivity and the crossings between the interconnects. Please note that the components and interconnects in Fig. 12 are laid out much farther apart from each other compared to the original AFCS layout in Fig. 11(b) to help visualize the different spatial diagrams clearly. Figure 14 shows 2 of 69 possible unique spatial topologies of the full-scale AFCS shown in Fig. 11(a) containing 14 components and 3 crossings.

Crossings | |||||||
---|---|---|---|---|---|---|---|

Components | 0 | 1 | 2 | 3 | 4 | 5 | 6 |

6 | 0 | 1 | 0 | 3 | 5 | 35 | 211 |

8 | 0 | 1 | 1 | 7 | 36 | 294 | |

10 | 0 | 1 | 3 | 37 | 283 | ||

12 | 0 | 1 | 5 | 74 | |||

14 | 0 | 0 | 5 | 69 |

Crossings | |||||||
---|---|---|---|---|---|---|---|

Components | 0 | 1 | 2 | 3 | 4 | 5 | 6 |

6 | 0 | 1 | 0 | 3 | 5 | 35 | 211 |

8 | 0 | 1 | 1 | 7 | 36 | 294 | |

10 | 0 | 1 | 3 | 37 | 283 | ||

12 | 0 | 1 | 5 | 74 | |||

14 | 0 | 0 | 5 | 69 |

Classes | Yamada polynomials |
---|---|

Class 1: | −1 − A + 2A^{2} − 3A^{3} − 2A^{4} + 13A^{5} − 24A^{6} + 3*A^{7} − 39A^{8} |

+ 42A^{9} − 30A^{10} + 23A^{11} − 5A^{12} + 7A^{14} − 3A^{15} + 2A^{16} | |

Class 2: | 2 − 4A + 2A^{2} + 2A^{3} − 18A^{4} + 24A^{5} − 41A^{6} + 37A^{7} |

−41A^{8} + 28A^{9} − 18A^{10} + 6A^{11} + 3A^{12} − 3A^{13} + 3A^{14} | |

Class 3: | −1 − A^{2} − 2A^{3} − 3A^{4} + 6A^{5} − 19A^{6} + 26A^{7} − 36A^{8} + 36A^{9} |

−35A^{10} + 25A^{11} − 16A^{12} + 4A^{13} + A^{14} − 4A^{15} + 2A^{16} − A^{17} | |

Class 4: | −1 + 2A − 5A^{2} − A^{3} − 6A^{5} − 8A^{6} + 18A^{7} − 36A^{8} + 45A^{9} − 47A^{10} |

+ 42A^{11} − 28A^{12} + 10A^{13} + 5A^{14} − 14A^{15} + 11A^{16} − 6A^{17} + A^{18} | |

Class 5: | −2 − A^{2} − 4A^{3} + 7A^{5} − 19A^{6} + 33A^{7} − 40A^{8} + 42A^{9} − 39A^{10} |

+ 24A^{11} − 15A^{12} − 3A^{13} + 5A^{14} − 8A^{15} + 3A^{16} − A^{17} | |

Class 6: | −1 + 3A − 6A^{2} + 4A^{3} − 2A^{4} − 3A^{5} − 5A^{6} + 14A^{7} − 33A^{8} + 41A^{9} |

− 49A^{10} + 42A^{11} − 33A^{12} + 13A^{13} + 3A^{14} − 15A^{15} + 15A^{16} − 9A^{17} + 3A^{18} |

Classes | Yamada polynomials |
---|---|

Class 1: | −1 − A + 2A^{2} − 3A^{3} − 2A^{4} + 13A^{5} − 24A^{6} + 3*A^{7} − 39A^{8} |

+ 42A^{9} − 30A^{10} + 23A^{11} − 5A^{12} + 7A^{14} − 3A^{15} + 2A^{16} | |

Class 2: | 2 − 4A + 2A^{2} + 2A^{3} − 18A^{4} + 24A^{5} − 41A^{6} + 37A^{7} |

−41A^{8} + 28A^{9} − 18A^{10} + 6A^{11} + 3A^{12} − 3A^{13} + 3A^{14} | |

Class 3: | −1 − A^{2} − 2A^{3} − 3A^{4} + 6A^{5} − 19A^{6} + 26A^{7} − 36A^{8} + 36A^{9} |

−35A^{10} + 25A^{11} − 16A^{12} + 4A^{13} + A^{14} − 4A^{15} + 2A^{16} − A^{17} | |

Class 4: | −1 + 2A − 5A^{2} − A^{3} − 6A^{5} − 8A^{6} + 18A^{7} − 36A^{8} + 45A^{9} − 47A^{10} |

+ 42A^{11} − 28A^{12} + 10A^{13} + 5A^{14} − 14A^{15} + 11A^{16} − 6A^{17} + A^{18} | |

Class 5: | −2 − A^{2} − 4A^{3} + 7A^{5} − 19A^{6} + 33A^{7} − 40A^{8} + 42A^{9} − 39A^{10} |

+ 24A^{11} − 15A^{12} − 3A^{13} + 5A^{14} − 8A^{15} + 3A^{16} − A^{17} | |

Class 6: | −1 + 3A − 6A^{2} + 4A^{3} − 2A^{4} − 3A^{5} − 5A^{6} + 14A^{7} − 33A^{8} + 41A^{9} |

− 49A^{10} + 42A^{11} − 33A^{12} + 13A^{13} + 3A^{14} − 15A^{15} + 15A^{16} − 9A^{17} + 3A^{18} |

### 5.3 Case Study 3: Components With Different Valencies.

In most real-world systems, every component may not have the same valency (number of ports). Unlike case study in Sec. 5.1, in this case study, we investigate a system with four components having different valencies. In addition, FVGs as described in Sec. 2.2 are used. As FVGs have local structures at nodes, the edge connectivity order around the nodes is preserved, and thus, FVG representations are highly suitable for design applications where nodes have a specific cyclic ordering of ports. Here, R0 to R5 moves are valid but not R6. Figure 15 shows some of the results obtained in this study. After computing the Yamada polynomials of hundreds of planar SGDs, a total of 27 unique Yamada classes are obtained. For illustration purposes, we show some isotopes of SGD_0 (original system architecture) as class 1 isotopes. Furthermore, unique SGDs belonging to some unique Yamada classes are shown for crossing numbers one, two, and three. Two final 3D system geometric models (referred as S1 and S2) are also shown in Fig. 15. The total computational time taken for study B was 211.4 s. It can be observed from this study that with components with different valencies, we get more unique Yamada classes than those with identical components. Thus, generating such designs manually is very challenging, and the automated enumeration framework we proposed here is very valuable.

### 5.4 Case Study 4: Circular Graph Representation.

While filtering out isotopic spatial graph diagrams in the previous case studies, we observed that in a few occasions where system topologies have many crossings, two edges in that diagram twist around each other multiple times. Although a higher crossing number is satisfied, unnecessary intertwining between edges can often be reduced by Reidemeister moves to a smaller crossing number, so essentially, no unique spatial topology is attained. In some cases, such intertwining is practically not observed or desirable between pipes or ducts in most complex systems (e.g., aero-engine externals, hydraulic systems). Some syntactic constraints need to be imposed to prevent more than a simple crossing between any two edges. This requires a representation that implicitly forbids twisting of two edges multiple times around each other. One way to get different spatial embeddings of an input abstract graph *G* as shown in Fig. 16 is to (1) pick an ordering of the nodes and use that to arrange them along a circle on the plane, (2) connect the nodes by straight lines corresponding to the edges of *G*. This gives the “shadow,” and (3) resolve the intersections lines of the shadow into over or under crossings. Figure 16 shows the shadow of graph *G* based on a particular cyclic order of nodes and one spatial graph embedding. As there are five crossings, a total of 2^{5} = 32 spatial embeddings are possible. The unique ones can be identified using the proposed design framework.

### 5.5 Case Study 5: Large-Scale System—Spatial Graph Decomposition Approach.

From the observations made in the previous case studies, it is evident that enumerating spatial topologies for most real-world systems containing many components and approximately hundreds of crossings is intractable with manual processes and can become computationally expensive with automated methods such as the one presented in this article. In contrast, enumerating spatial topologies of each subsystem of components can be a simple and efficient process. A complex spatial graph can thus be converted to a set of subgraphs, and the unique spatial topologies of these subgraphs can then be enumerated separately. The subgraphs can be decoupled and can be considered as super-nodes. This decouples the task into two subasks: (1) Enumerate STs of the system graph with only the subsystems as super-nodes, and (2) enumerate unique STs within each subgraph. This presents fewer design candidates about which to make decisions, which greatly reduces the overall computational expense. Figure 17 shows a random complex spatial graph with 14 nodes, 20 edges and allowing at most 10 edge crossings. Approximately 1.134 × 10^{4} SGDs are attained for this entire system that fall under 434 unique Yamada polynomial categories. As this is a very large set, decomposition of the graph into subgraphs (as super-nodes) is appropriate. First, a unique spatial topology of the super-nodes graph is found. Case study 1 in Sec. 5.1 is utilized as a subgraph for demonstration purposes. Note that while enumerating STs for the spatial subgraph, the rest of the system is condensed as an extra node in the subgraph to preserve spatial connectivity information. Finally, by using the proposed design framework, unique STs of the subgraph can be plugged into the original system to attain system configurations. The scope of this article only deals with enumerating unique STs, so we plan to show how each of these unique topologies affect overall system performance in future work.

To explain this decomposition concept using a concrete engineering design example, suppose that the complex network represents the spatial topology of a hybrid electric vehicle powertrain; one possible subsystem could be a fluid-thermal cooling circuit. Each distinct circuit topology can be geometrically optimized for fair comparison, revealing how the topological features contribute to the overall system efficiency, fuel economy, thermal loss management, and other figures of merit due to physical interactions between components, interconnections, and the environment. The best candidate ST can then be chosen according to the desired performance requirements, as done in Refs. [22,26], where the same procedure was followed but with the goal of ranking different SAs. As a part of future work on this topic, we plan to investigate alternative methods, such as deep learning and pattern recognition as mentioned in Refs. [79,95,125], to efficiently explore new topologies of large-scale systems where exhaustive enumeration of all possible topologies may not be tractable.

## 6 Discussion

As described in Sec. 1.1, this article attempts to answer several fundamental research questions (not considered thoroughly in existing literature) on the design automation of spatial topological decisions and their implications on engineering system design. Based on the aforementioned case studies, several design insights were obtained. This section summarizes a list of important observations made from the five case studies as follows from an engineering design perspective and how those inferences help answer some research questions to a reasonable extent:

First, the enumeration framework support representation of 3D systems as spatial graphs and allows an abstract representation of systems as graphs with nodes and edges. Especially from case studies presented in Secs. 5.1 and 5.2, it can be observed that the number of unique spatial topologies attained for a given interconnect crossing complexity is much smaller than the combinatorially enumerated set of spatial graph diagrams, as most of them are isotopic to each other under the Reidemeister moves. Furthermore, implementation of the framework for enumeration of fuel cell system topologies in 5.2 helps to identify nonobvious designs, thus providing value to creative design efforts for real-world engineering systems. This is a very promising result as navigation of 3D topological space efficiently is now possible and initial design filtering would help explore only useful candidate solutions.

Another contribution of this work is embodied in case study in Sec. 5.3, where STs are enumerated for components with different valencies, in contrast to existing work [73–77] that is mostly limited to two or three equivalent vertices. This would be very useful in extending the framework to systems of various scales and requirements, including hollow objects and void regions. Artificial intelligence (AI) and machine learning (ML) techniques can be easily used through the data generated from such diverse design solutions.

The circular graph representation method, presented in case study presented in Sec. 5.4, is a simple way to enumerate and realize SGDs and avoid edge intertwining, although Yamada polynomials should still be used for identifying unique STs. Furthermore, specific syntactic constraints can be added to significantly reduce the initial set of SGDs obtained for planarity checking and Yamada polynomial evaluation. For example, by adding constraints on total crossings allowed between two edges of a system, there is a greater control on the type of spatial topologies finally obtained. This will be studied more in future work. Imposing syntactic constraints will be very useful to address the issue of interconnects having different physics spatial fields and optimal design trajectories for accessibility. Such constraints would be very useful to eliminate infeasible designs and either introduce or remove crossings between interconnects to satisfy operational, physics, or accessibility requirements.

As seen in case study in Sec. 5.5, for large-scaled systems, the best way so far to achieve different STs and search effectively is by graph decomposition. The spatial graph of the subsystem, which plays a critical role in performance impact, can be extracted to find its unique topologies. This avoids the need to enumerate thousands of diagrams of a complex network, compute their polynomials, and compare them. Moreover, subgraph designs can be optimized for performance independently and then combined with the remaining system. In the decomposition-based approach, the isotopy equivalence check is first performed on subgraphs and then on the system graph containing the super-nodes of each subgraph. This case study is very helpful is addressing the question if it is possible to address subsystem level design challenges. For example, if certain portion of a larger system needs to be replaced or accessed for repair. Then graph decomposition would be very helpful to see accessibility scenarios for that subsystem with respect to the entire larger system, and similarly, module-based repair is possible, and the cost for overhaul becomes less expensive.

Another impactful aspect of this framework is that for one system architecture, there can be a range of nonobvious spatial topologies from zero to many crossings. The spatial embeddings with fewer interconnect crossings are generally useful for practical engineering purposes where reduction of system geometric complexity is desired. Therefore, for existing, complex real-world system designs with 10s of crossings, using the proposed design framework, a simpler spatial topology may be found for that network with a much lower crossing number, but still keeping the same system connectivity. However, sometimes there might exist an implicit design requirement where having multiple interconnect crossings supports a system function. For example, applications involving wire harnesses, twisted cable bundles, etc., may benefit from additional crossings to support tighter space constraints or provide mechanical support. Hence, this framework can be utilized to satisfy different design requirements.

Several research questions such as a new design representation for 3D systems, scalability, time complexity, operational and maintenance considerations, simple geometry, crossing numbers, and design automation for graph generation have been addressed in this article, and the scope for further research on this topic is now expanded. In summary, advantages of using a spatial graph representation are as follows: (1) its simplicity, while capturing necessary system elements and topological features, (2) ease of visualization, (3) flexibility to add new geometric features like size and shape, (4) the ability to detect distinct topologies using polynomial invariants, (5) scalable or even decomposable into a set of smaller graphs, (6) supports automated 3D model generation, and (7) accommodates features such as node locations, edge diameters, edge trajectory shape functions, port locations, crossing information, and other elements that can be parameterized for performing continuous numerical optimization. For example, items 6 and 7 can be very useful for using different 3D models as initial start points for multiphysics component placement and routing optimization. The features of the spatial graph design representation address several fundamental research questions; however, further investigation on this topic is needed to address the limitations of this framework and find out solutions to complex design issues.

## 7 Conclusion

The design representation presented in this article greatly enhances the study of unique 3D engineering system spatial topologies in a systematic manner and is supported by rigorous mathematical foundations in the spatial graph theory. Topologies of complex engineering systems, designed for particular applications, are conventionally created manually. But for more effective performance and design process efficiency, systematic identification, enumeration, and classification of possible system topologies can aid thorough navigation through challenging 3D discrete design spaces. A framework for representing three-dimensional interconnected engineering systems using spatial graph embeddings is presented. Initially, all the combinatorial spatial graph descriptions up to some fixed topological complexity are enumerated for an input system architecture. A polynomial invariant, the Yamada polynomial, is then calculated for the set of all the spatial graphs attained from the combinatorial permutations. The Yamada polynomial helps identify duplicate spatial graph topologies from the exhaustive set and a smaller set of unique spatial embeddings (equivalent topological classes) is obtained. This smaller set of spatial graphs can be used for generating three-dimensional geometric system models. Five case studies have been demonstrated using the proposed enumeration strategy, including implementation of an industry application, the AFCS. The results show that this method is efficient, scalable, and applicable to all general 3D interconnected system networks; allows comprehensive exploration of the design space; and greatly aids in the design and development of unprecedented system topologies.

Future work includes adding more geometric features to these spatial graph embeddings, such as representing nodes with geometric shapes and ports. Investigation of braid-based representations of interconnect networks is also anticipated. As designed systems become larger, evaluating Yamada polynomials for many SGDs is very time consuming. This can be overcome by implementing a mix of Reidemeister moves to eliminate isotopic diagrams quickly to produce a smaller set of diagrams that require Yamada calculations. Other application aspects include utilizing the unique spatial topologies obtained here as starting points for physics-based component placement and routing optimization of 3D systems. Furthermore, research areas that can benefit from SGD representations are 3D pipe routing, topological 3D path planning for robotic operations, aerial drone navigation, generation of new automotive cooling system configurations, 3D integrated circuit interconnect technology, and many others. In this article, the general concept of 3D spatial topology enumeration using spatial graphs is discussed. We hope that this initial work serves as a strong foundation to bridging the gap between engineering design and mathematical low-dimensional topology. There are many interesting aspects that are yet to be explored and can have a great impact when applied to practical engineering design problems.

## Acknowledgment

This material was based upon work supported by the National Science Foundation Engineering Research Center (NSF ERC) for Power Optimization of Electro-Thermal Systems (POETS) with cooperative agreement EEC-1449548. The opinions, findings, and conclusions or recommendations expressed are those of the author(s) and do not necessarily reflect the views of the National Science Foundation. Dunfield was partially supported by NSF (Grant No. DMS-1811156) and the Simons Foundation. The authors also acknowledge industry partners from Ford Motor Company who provided the representative data for the automotive fuel cell system (AFCS). The authors also thank Kyle Miller from the University of California, Berkeley for his helpful discussions on topics related to spatial graphs and Yamada polynomials.

## Conflict of Interest

There are no conflicts of interest.

## Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.

## References

_{2})

*Introduction to Design Processes*, 2nd ed., Department of Mechanical Engineering, Monash University Clayton, Victoria, https://nla.gov.au/nla.cat-vn556118