A Literature Review of Hybrid System Dynamics and Agent-Based Modeling in a Produced Water Management Context

This paper explores the possibility and plausibility of developing a hybrid simulation method combining agent-based (AB) and system dynamics (SD) modeling to address the case study of produced water management (PWM). In southeastern New Mexico, the oil and gas industry generates large volumes of produced water, while at the same time, freshwater resources are scarce. Single-method models are unable to capture the dynamic impacts of PWM on the water budget at both the local and regional levels, hence the need for a more complex hybrid approach. We used the literature, information characterizing produced water in New Mexico, and our preliminary interviews with subject matter experts to develop this framework. We then conducted a systematic literature review to summarize state-of-the-art of hybrid modeling methodologies and techniques. Our research revealed that there is a small but growing volume of hybrid modeling research that could provide some foundational support for modelers interested in hybrid modeling approaches for complex natural resource management issues. We categorized these efforts into four classes based on their approaches to hybrid modeling. It appears that, among these classes, PWM requires the most sophisticated approach, indicating that PWM modelers will need to face serious challenges and break new ground in this realm.


Introduction
Our contemporary managerial problems, such as produced water management, are becoming increasingly complex, and we need to equip ourselves with modern analytical tools and modeling approaches to address these problems. Produced water, the brine water in a geological formation and flow-back water from the hydraulic fracturing process that is a coincidental byproduct of oil and gas production, has high variability in volume and quality at the local scale, and has broader implications for freshwater availability at a regional scale. We report here our exploration of the necessity of developing a hybrid system dynamics (SD) and agent-based (AB) simulation approach for evaluating the regional water budget impacts of produced water policies and management decisions.
Hybrid simulation approaches are needed to advance our understanding of future problems and to solve them effectively while mitigating the unintended consequences of our solutions [1,2]. In this regard, the primary goals of this research were to (1) identify the main dynamic characteristics and complexities that a comprehensive produced water management (PWM) model should take into account, (2) review the current body of literature where such complexities could be addressed using the particular hybrid simulation modeling approaches, (3) assess the necessity and usefulness of hybrid modeling to PWM issues, and (4) provide recommendations for future hybrid modeling of PWM. A secondary goal of this research was to bring together common terminology used in the modeling literature in order to guide future hybrid modeling efforts. This paper is organized as follows. The literature about hybrid modeling is reviewed in Section 2. Section 3 presents the history and background of produced water in New Mexico. We discuss in Section 4 why a hybrid modeling approach could be effective for understanding the complexities around PWM. In Section 5, we develop and use a conceptual framework to explain the different levels of complexity in decision-making processes for managing the issues associated with produced water, and its impact on local and regional areas. Section 6 describes our literature review approach that identifies previous efforts of hybrid modeling in various contexts. Section 7 presents the results of this review and the insights we gained through this process. Based on the findings, we categorize distinct approaches for hybrid modeling that should be most applicable for future research. Section 8 then concludes the paper by summarizing the necessity of hybrid SD-AB modeling for produced water management issues. It also provides modelers with guidelines for proper integration of system dynamics, agent-based modeling, and geospatial data for the specific problems associated with produced water.

Prior Work
Hybrid modeling has many different forms and types, and there is currently no clear and cohesive definition for it [3]. Two examples of differing hybrid modeling definitions found in the literature include "an approach that merges recent advancements in nonparametric analysis with standard parametric methods" [4], and "mathematical models that can handle various types of information and combine diverse theoretical methods on multiple temporal and spatial scales" [5]. Eldabi and others [6] attribute this lack of consensus to "the very nature of hybridization where models are based on mixing several paradigms, making it difficult to be housed within one." Here, we simply use the term "hybrid modeling" to refer to a process of combining two or more dynamic simulation methods, in particular, system dynamics (SD) and agent-based (AB) modeling. In a broader context, hybrid simulation, also known as multi-paradigm simulation, is usually defined as any combination of the three main simulation paradigms, i.e., SD, agent-based modeling (ABM), and discrete event simulation (DES) [7].
To address the dynamic complexity of PWM, we can take different dynamic simulation approaches. The main approaches to consider are system dynamics (SD) and agent-based (AB) modeling. In theory, pure SD or AB models could be applied to any dynamic problems. Each of these approaches has its own strengths and weaknesses. SD models are efficient computationally, have great clarity of exposition, and provide easily tractable analysis [8]. AB models, on the other hand, have an advantage with respect to expressing and characterizing heterogeneity, and can also include spatial interactions within and between agents and their environment [9]. SD models could be designed to take heterogeneity into account by the use of subscripts or arrays. Ruth [10] provides one of the earliest examples of this kind of modeling. However, this approach is inflexible in terms of interactions between agents as explained in detail in BenDor and Kaza [11]. As an alternative approach, Rahmandad and Sterman [12] show how system dynamics could be used to represent an approximate AB model. Roach and Tidwell [13] applied a spatial variation of the Compartmental SD (CSD) approach to groundwater resources management. However, it was shown that the simulation results of an AB model that fully accounts for network structures differ significantly and substantially from a CSD that does not [8].
AB models could also be applied to any dynamic problem. Like SD models, they can take feedback mechanisms and nonlinearities into account. However, unlike SD models that treat feedback loops as the main unit of analysis [14], AB models focus on agents as the unit of analysis [15]. Also, compared to SD models, they are more difficult to validate and verify, and lack effective architectures and protocols to represent agents and their interactions [9].

Contextual Background
Oil and gas production is increasing exponentially in New Mexico (Figure 1), meaning that the produced water management issue will become increasingly significant in the future. Oil and gas production in New Mexico began as early as the 1920s using conventional drilling techniques and remains a substantial source of the state's revenue. For instance, oil and gas production accounted for 23.8 percent of the state's revenue in 2008 [16]; this percentage varies by year and fluctuates with the price of oil and gas. New Mexico oil production peaked temporarily during the 1960s; however, a new boom starting in 2011 has made the state the third largest oil producing state in the U.S., yielding 329.4 million barrels of oil and 1.8 billion MCF of natural gas in 2019. Approximately 97 percent of oil and gas production in New Mexico occurs in the southeastern corner of the state. The discovery of what is currently considered the world's largest unconventional oil play within the Permian Basin [18], and advancements in drilling and production techniques, have renewed the importance of Southeast New Mexico for energy production in the national arena. Particularly, Lea and Eddy counties ( Figure 2) are two of the top oil and gas producing counties in the United States. In 2016, there were 46,232 oil wells and 8045 gas wells operating in these two counties.
A recent assessment by the United States Geological Survey estimated 46.3 billion barrels of oil and 281 trillion cubic feet of gas are recoverable in the Permian Basin in Southeastern New Mexico and western Texas [19]. Increases in oil and gas production causes increases in produced water, such as over 42 billion gallons of produced water generated in New Mexico in 2018 [20], that imposes significant economic and environmental challenges to oil producers and to society generally. Various end uses for treated produced water have been suggested, such as agriculture, potash mining, energy production, surface water discharge, and managed aquifer recharge [21]; however, the current regulations and public concern prohibit produced water use outside the oil and gas industry. The large volumes of available oil and gas in the region solidify the future of the oil and gas industry in the region for the foreseeable future.
The volumes of produced water commonly range between a produced water to oil ratio of 3:1 and 13:1 [18], and contain varying levels and composition of dissolved solids [22], making it expensive to treat and dispose of [23], and thus it has remained a major challenge for policy makers, industry, communities, and environmental protection agencies [24]. Total dissolved solids in produced water in the western United States ranges from 1000 mg/L to greater than 400,000 mg/L [25] and samples taken at different times from the same well can vary more than 100,000 mg/L [22], which have direct implications for treatment costs and availability of technological treatment solutions. The six common disposal options for produced water are: discharge, underground injection for disposal, underground injection for reservoir pressure maintenance to increase oil recovery, evaporation ponds, offsite commercial disposal, and beneficial reuse [26]. Management of the freshwater supply at the regional level in southeastern New Mexico is under increasing scrutiny as demand increased due in part to increase in of the oil and gas production, which requires water for enhanced recovery and hydraulic fracturing. In Eddy and Lea Counties, the overall trend of annual precipitation is decreasing and abnormally dry to exceptional drought continues to be a frequent occurrence. A lack of a reliable surface water supply and dwindling groundwater supplies in the region require a detailed accounting of the regional water budget for water planning.
Policy decisions and management choices of produced water disposal can have both local and regional impacts on the volume and quality of local and regional freshwater supplies, as well as on seismicity levels, transportation infrastructure, oil and gas production costs, soil quality, and ecosystem health. These impacts also link to societal effects such as employment, quality of life, environmental advocacy, and agricultural sustainability.

Why a Hybrid Modeling Approach?
Ultimately, the goal of the research question drives the necessity of a particular modeling approach. To justify a hybrid modeling effort, certain aspects should be important for the research question [27]: • Explaining how relationships emerge and evolve among agents (e.g., the geospatial distribution of disposal or injection wells). • Explaining how these relationships affect the state of the system (total cost, oil production, water levels, and so on, as influenced by the dynamics of the geospatial distribution of oil and water wells).
• Explaining how the state of the system affects the relationships (e.g., the distribution of oil and water wells as influenced by oil production and groundwater levels).
Our research questions include all of these three aspects, and thereby necessitate a hybrid modeling approach. As described in Section 3, the dynamic interactions of PWM decisions and regional water budget are difficult to model across space and time because of the dynamic nature of many of the variables involved. For example, water demand for hydraulic fracturing occurs only in the beginning of the well life cycle and is on the order of days, whereas produced water volumes typically follow a logarithmic curve throughout the well life on the order of years. Injection wells are not evenly distributed, and the geochemistry of produced water must be compatible with the geologic formation the saltwater disposal well is drilled into. Produced water treatment facilities can be centrally located, but options for mobile treatment units are becoming more prevalent. Similarly, transporting the produced water, either to a treatment facility or injection well, is done either by trucking or pipelines-with each form of transportation having its own set of feedbacks into the larger system. Selection of treatment options for produced water is primarily driven by the treatment type, feedwater quality and volume, energy cost, and the intended water use.
Analysis of these individual factors alone does not lead to solutions that account for the multiple levels of interaction either driving or affecting the outcomes of produced water management. In fact, the body of scientific literature is rich with disciplinary research addressing many facets of produced water management: policy [24], treatment technology selection [23], geochemical composition [22], risk assessment [28], and increases in seismicity [29]. The literature is also rich with examples of a growing interest in utilizing a hybrid modeling approach to support natural resource management decisions-for example [30]. However, applications that consider all key attributes of a typical social-ecological system such as feedback, nonlinearity, cross-scale dynamics, and heterogeneity in a single package for water management problems are nonexistent [31].
To address the cross-scale dynamic complexities of PWM, we require research methods and tools that can characterize and represent nonlinear system-level interactions as well as heterogeneous and spatial interactions over time. Single-method analytical solutions are not adequate for this purpose because they cannot seamlessly integrate these system and individual levels of analysis. Advanced dynamic simulation approaches are needed to fill the gap [2]. In this paper, we explore the potential for applying a hybrid dynamic simulation approach to PWM. We ask what the minimum boundary is for a comprehensive model of produced water that aims to capture its important dynamic complexities. By using the literature, produced water data from New Mexico, and our preliminary interviews with subject matter experts, we develop a conceptual framework of the problem to guide us through this inquiry. Then, we explore our methodological options to examine whether single-method) dynamic simulation approaches, such as ABM or SD, are sufficient to tackle the issue. We then ask if a hybrid modeling approach would add any net value (benefits minus costs) to this area of research. We carry out a systematic literature review to answer these questions and to facilitate our exploration of deploying a hybrid model, using water-scarce southeastern New Mexico as a case study.

Conceptual Framework
The conceptual framework presented in this section emphasizes the multiple aspects of the PWM problem that are key to our understanding of important issues such as the impact of PW on the dynamics of regional water budgets. To understand how PW can change dynamics of water budgets in a region, we need to take full complexity of the problem into account. Otherwise, disjointed information, even though in the same language, may not provide much insight. Produced water management decisions are made and being influenced at multiple levels of complexity. Environmental regulations and water quality requirements change as public perception of PW environmental risks change. These changes impose new constraints on PW management options by altering the cost functions of PW treatment, disposal, transportation, etc., thereby affecting oil and gas production strategies, including the geographical location of wells, capacity utilization, and so on. Changes in production patterns will affect future trajectories of PW volumes, where they are generated, where they are disposed, and how they affect the dynamics of quantity and quality of water resources. Dynamics of the water budget then feed back to the system to drive both regulatory and management decisions that further drive the changes in the system. The illustrative interactions described involve three key characteristics (described below), which require a hybrid SD-AB modeling to address the question of how PWM will affect water budgets. (a) Being dynamic: Dynamic complexity emerges from the interactions among the agents over time [32]. In our framework, the key variables of the system such as oil production, water use for hydraulic fracturing, produced water used in secondary recovery and reservoir pressure maintenance, and decisions for treatment and disposal, all interact with each other and all change over time.
Being spatial: It relies heavily on spatial information and data, as managers must make decisions on where to drill a new oil well (where the produced water is generated), where to dispose of the produced water, and where to inject treated produced water. (c) Being heterogeneous: Agents representing stakeholders, wells, or well owners act differently based on their different input, analysis, and interests; this heterogeneity also adds to the complexity of the problem.
In addition, our framework should accommodate and clarify the mechanisms by which particular stakeholders make their decisions and how they are impacted by those decisions. In order to achieve this integrated model, we find there is a need to study these decision-making processes with respect to their impact on at least three distinct areas of interest:

1.
Individual Companies:At this level, there could potentially be two types of agents ( Figure 3). First, there are oil companies that impact the system by making decisions such as where to drill new wells and how to deal with the produced water at each location. Second, there are oil wells with different specifications regarding the volume of produced water and associated geological formations. The main attribute that differentiates oil companies as different groups of agents in our model is their size. According to our interviews with research and industry experts, major oil companies are socially driven to explore options for using less freshwater, to treat and reuse produced water, and to reduce impacts to the environment. Compared to large players in the system, independent oil companies usually have more immediate considerations and fewer resources for long-term investments such as large-scale produced water treatment. At this level, although many factors are considered for PWM, cost is the ultimate driving force behind management decisions, followed by the need to maintain a positive public perception. Specifically, when the profit margin of production drops and remains below a certain threshold for long enough, production ceases and the well is closed [26]. On the other hand, for agents representing the wells, geospatial attributes such as their distance from proper disposal areas and the transportation cost of the produced water will drive the management decisions.

2.
Local Community: The impact of produced water management decisions can be seen mostly at the local level where resulting pollution directly impacts the environment ( Figure 3). For example, some produced water may spill during transportation, or it may be partly responsible for an increase in the seismic activity in nearby areas; such examples may impose significant challenges to the local population. One of the main drivers of the changes in regulations regarding produced water is the public pressure on regulatory institutions. Each new regulation requires the oil companies to modify their decisions toward better environmental outcomes. These decisions change other spatial variables such as quantity and quality of available water, seismicity risks, environmental pollution, and so on. These changes drive economic and system-level changes such as the water budget, environmental regulations and policies, and carry societal costs. These factors then feed back into the decision functions of produced water managers as informational inputs for their cost-benefit analyses.

3.
Aggregate region: Treatment of the produced water can potentially reduce the amount of water available for other activities such as agriculture or industry (Figure 3). For example, based on the quality of the treated produced water, it can be reused for fracking. This process could reduce the need for freshwater, and therefore, reduce the extraction of water from almost all non-renewable groundwater aquifers in the region.
Because of the level of aggregation, more research is needed in order to connect the cause and effect processes. This goal can only be achieved by using integrated tools such as hybrid modeling.
An overview of our PWM conceptual framework can be seen in Figure 3, which also shows the potential inputs and outputs of the model.

The Systematic Review Method
In order to support a dynamic simulation approach, we conducted a literature review focused on hybrid modeling beginning with a series of trial and error searches to identify which methods or combinations thereof, among all the dynamic simulation approaches, would be suitable for our problem. We presumed that we needed a dynamic simulation approach (e.g., system dynamics) that could explicitly take into account feedback at the system level, and at different geographical locations. Therefore, our initial review started with several potential search terms for locating relevant literature for hybrid modeling approaches. We paid particular attention to the literature of spatial dynamic modeling, for example, Roach and Tidwell [13], BenDor and Kaza [11], and Neuwirth, Hofer, and Schaumberger [33], which led to the idea of combining Cellular Automata (CA) with SD in order to capture spatial dynamics [34]. However, as we described earlier, we had another layer of complexity to consider and that was individual decision-making processes. Since AB modeling could be used as an advanced platform for CA modeling [35], we came to the conclusion that a SD-AB hybrid modeling approach would most likely provide the minimum technical complexity that we needed to deploy in order to achieve our goal. Consequently, we focused on these two dynamic simulation methods in our next round of literature review.
The method and application papers were assessed for usefulness based on the criteria of containing methodological conceptualization, practical technical guidelines, or model codes or equations. The goal was to identify the current state-of-the-art of hybrid SD-AB modeling and to provide a useful guideline for those who want to model produced water management issues. The literature reviewed for this paper is the result of searching the Web of Science for publications that contained both the terms "system dynamics" and "agent based." The initial resultant 212 papers were reviewed to determine if the papers were describing a hybrid modeling approach, and if so, they were first sorted into one of three categories: review paper, method paper, or application paper. A full listing of the 212 initially selected papers is provided as a supplementary spreadsheet data that accompanies this paper. Among these papers, we identified 77 papers as relevant to the purpose of our research. These papers provide useful information for how a hybrid SD-AB model can be developed including guidelines for identifying the kinds of problems that could benefit from a hybrid approach, the conceptualization of generic structures that could be applied to some specific problems, and example applications including codes or equations that could inform modelers as to how they might implement the method for their problem.

A Hybrid Dynamic Simulation for Produced Water Management
Here, we analyze the outcome of our literature review by focusing on the 77 relevant papers identified in the previous step. The temporal distribution of these papers is presented in Figure 4. The historical trend shows an exponential growth of SD-AB hybrid modeling efforts. A meaningful interpretation of this growth, however, requires a comparison with the general trend of scientific publications, and possibly an analysis of impact factors for the journals in which these publications appeared. In general, entering this area of research (hybrid simulation modeling) is considered challenging, and sometimes, daunting [36] mainly due to a lack of formal training or educational material or textbooks available to participants at the outset [3,37], the need to acquire sufficient computer programming skills by subject matter experts [38], and the limited availability of software packages that can adequately and easily handle the integration of multiple approaches [27]. The majority of the relevant papers (i.e., 39 papers) focus on the application of hybrid modeling in different contexts; 28 papers investigate the hybrid modeling methodology; and the rest, 11 papers, review the literature. The relatively large amount of method and application papers was a promising sign that we might be able to find some practical instructions for how a hybrid model could be effectively and efficiently developed for the case of PWM, especially because the majority of these applications were in the realm of natural resource management. Figure 5 shows that 15 papers were related to natural resources and environmental issues; 11 papers to energy issues; and 7 papers to water. Of these papers, only one investigated a water-energy nexus problem that is potentially more relevant to PWM. Breakout of papers by application area Figure 5. A breakdown of reviewed papers according to their research application area.
The reviewed papers mainly provide comparisons between the dynamic simulation methods and how and when each of these methods or their combination will be more useful. For example, Lattila, Hilletofth, and Lin [39] discussed the possible ways of combining SD and AB modeling. The paper identifies five different situations where it will be useful to do so. Depending on the characteristics of a given problem, the paper tries to suggest the most suitable approach. For another example, Brailsford and others [40] explored the possibility and plausibility of hybrid modeling applications in the realm of Operations Research (OR). They investigated the challenges of application and hybridization and provide a conceptual framework for how to integrate these methods for OR cases.
The methodology papers are useful for providing theoretical foundations for hybrid modeling. For example, Anderson, Lewis, and Ozer [27] developed a framework for how SD, AB, and network modeling and analysis could be combined using Vensim TM software, a system dynamics modeling platform [41]. Another example is Duggan [42], which introduced a method for integration optimization in an agent-oriented SD framework in the context of supply chain management. The language used for this work is XMILE [43].
The application papers are useful by providing real-world examples of hybrid modeling so users can learn about the practical challenges involved in the process. An example of applied hybrid modeling papers is Kieckhaefer, Volling, and Spengler [44], which presents a SD-AB hybrid model to analyze electric vehicle markets in Germany. The model is implemented in AnyLogic [45], which is a hybrid modeling platform. Swinerd and Mc-Naught [46] provided another example of applied SD-AB hybrid modeling. They created their model in NetLogo [47], a modeling platform primarily used for agent-based modeling (ABM), to analyze the problem of diffusion of innovation in an international setting.
The 77 papers on our list provide useful information for developing hybrid SD-AB models. However, the level and type of usefulness of each paper varies depending on how deep it digs into the actual modeling processes. We previously discussed how different types of papers (literature reviews, methodological, and application-based) can provide different kinds of guidance for practical modeling. The papers could also be broken down by the level of detail they provide for practical modeling. Out of the total 77 papers reviewed, 53 papers provide some sort of technical guidance, usually in the form of a conceptual framework. Among these papers, 48 provide partial technical help such as detailed diagrams or some coding or equations. Out of these 48 papers, only 31 provide a complete listing of model codes or equations. We consider this last group of papers to be the most useful for practical hybrid SD-AB modeling.
We paid particular attention to the papers that provide detailed technical instructions, because these papers are the most useful for those who are seeking to develop hybrid SD-AB models. However, the models presented in these papers are not all similar in terms of structural and methodological design.
Inspired by Shanthikumar and Sargent [48], Swinerd and McNaught [36] classified yhybrid simulation models into three broad classes: interfaced, sequential, and integrated. The interfaced models are those that have modules of different methods that work in separate environments. The only connection between these modules is an interface that integrates their outputs. Venkateswaran and Son [49] provided an example of this kind by combining SD and discrete event simulation (DES) modules. The sequential models are those in which a set of modules run first by one method to provide input for a set of modules that run by another method. There is no feedback from the latter to the former. Mazzoleni and Massheder [50] presented an example of sequential hybridization by introducing a platform that connects a system dynamics software (Simile) to GIS. Ahmad and Simonovic [51] also provided a similar approach. The integrated models are the only group of hybrid models in which feedback exists between modules of two or more methods. Our work is focused mainly on this group of hybrid models because we believe this is what we need in order to address the full complexity of PWM issues.
The integrated models can be categorized using a finer level of classification to help the modelers better understand how the current modeling approaches work. Swinerd and McNaught [36] suggested three subclasses for integrated models: models that include: (1) agents with rich internal structure, (2) stocked agents, and (3) parameters with emergent behavior. In the first subclass of models, some agents of the AB model contain some sort of stock and flow structure. In the second subclass, some stocks in a system dynamics module contain agents that follow some specific behavioral rules. In the third subclass, one variable or parameter of a system dynamics module influenced by agents that follow some specific behavioral rules.
The classification of integrated models by Swinerd and McNaught [36] was not meant to be exhaustive but illustrative. As a result, it is too narrow with respect to some aspects, and too broad with respect to others. This classification is too narrow because it omits some important modeling approaches that need to be distinguished from others. It is also too broad, because for practical modeling purposes, separating subclasses 2 and 3, despite being different technically, does not add meaningful value. Therefore, we propose a modification to the Swinerd and McNaught [36] classification of integrated hybrid models as follows.
Class A: This class is missing from Swinerd and McNaught [36]. In this class, there are two separate sets of AB and SD modules that work in parallel. AB modules have architectural design and philosophy that are independent of the SD modules, but they can talk to each other through a protocol.

•
An example in the literature: Schieritz and Größler [52] provide an example of this kind that can be captured only partially by Swinerd and McNaught's subclass 1.
Their hybrid model addresses a supply chain management issue. The model's agents (companies within the supply chain) have two system dynamics modules: ordering and evaluation. The AB model is written in eM-Plant while the SD modules are developed in Vensim. A third module stores and processes the input and output data of the system dynamics modules and regulates the communication between Vensim and eM-Plant via Dynamic Data Exchange (DDE), which is a communication system. In this model, eM-Plant connects to Vensim (the DDE server), and the input and output data and commands are transferred via the established channel. • An example in PWM: A SD module of oil production in a SD platform such as Stella or Vensim, and an AB module of spatial dynamics of produced water injection wells in an ABM platform such as NetLogo communicate through an external protocol that runs through input/output spreadsheets. The SD oil production module simulates the volume of produced water over time. This simulation output will be exported to the SD output spreadsheet, which will serve as the AB input. The AB model reads the simulated produced water and distributes it between the injection wells based on a decision rule, e.g., a cost minimization function that determines which injection well is closer to each produced water disposal pond.

Class B:
This class is similar to subclass 1 in Swinerd and McNaught [36]. That is, agents of the AB model contain SD structure. There is a subtle difference between Classes A and B. In Class A, SD and AB modules are separate in terms of design and structure. In Class B, in contrast, there is no real separation between the modules. SD model codes or equations are written within the AB modules. This requires the use of the same modeling platform for both sets of modules.
• Example in the literature: An example of this kind is Duggan's model [42], although the model is fully developed using SD tools and called "agent-oriented SD" by the author. In this model, each player (agent) within the supply chain has a stock and flow structure. The output of these SD models then drives the behavior of the rest of the AB model. • Example in PWM: An AB module of spatially distributed injection wells where each well, as an individual agent, contains a SD structure. The SD structure could be a stock representing the capacity of the well that is available for additional injection of produced water. The produced water could be allocated to an injection well that has the greatest remaining capacity and that is closest to the point of distribution, thereby having the minimum transportation cost.

Class C: this class is a combination of subclasses 2 and 3 in Swinerd and
McNaught [36]. That is, agents of the AB model are part of a SD structure. Similar to what we described in Class B, there is a subtle difference between Classes A and C, which has important practical modeling implications. Class C has more flexibility than Class A as changes to the model structure will not require alteration of communication protocols between AB and SD modules. However, the AB and SD modules need to be written in the same language or at least closely compatible platforms.

•
An example in the literature: An example of this kind is Anderson, Lewis, and Ozer [27], which investigates the dynamics of team performance in knowledge-based organization. In their model, expertise is modeled as stocks while interactions between members and diversity-based subgroups are agent-based. In general, each variable in the SD model is subscribed to work as a small AB module. • An example in PWM: A SD oil production module connects directly to an AB module of spatially distributed injection wells. This example is similar to the Class A example with one major difference: the SD and AB modules are connected in a single environment and the interaction between these two modules occurs directly at each time step without any mediation.

Class D:
In this class, some SD variables are driven by AB interactions, while some AB variables receive information from SD variables. This is, in our view, the most sophisticated approach to hybrid modeling as it involves a natural and fluid hybridization that follows a unified modeling philosophy, architecture, design, and implementation. For the same reason, this is also the most difficult modeling approach, as it takes a lot of preparation in terms of thinking and design before the modeling begins.

•
An example in the literature: Swinerd and McNaught [36] do not mention this class explicitly. However, the model they provide in their other works [46,53], reveal that they acknowledge the existence of this class which combines Classes B and C. The only other instance of this class in our review is Alfaris and others [54] that present a model for national energy planning in Saudi Arabia. • An example in PWM: The AB module of injection wells explained in the Class B example is connected to a SD module of oil production in a feedback loop. Institutional dynamics driving regulation changes could also be modeled as a SD module. Regulations then will feed into the decision functions in the AB module and affect total costs of produced water management options. These decisions will impact the regional water budget, which could be represented by another SD module. The outputs of the water budget (e.g., water availability index), in turn, feed back into the other AB and SD modules to determine new institutional regulations and PWM decisions.
We believe that PWM issues require a Class D hybrid modeling, as suggested by our conceptual framework (Figure 3), because of the complex feedback structure that connects different levels of the system through irregular sequences. For example, oil companies' decision-making that could be an AB module is part of a system-level feedback, that is, public perception of risks. This part of the model would be Class B. The hydrologic dynamics, which is probably a system dynamics module, would need to be replicated in different locations. This part of the model would be Class C. Therefore, the whole model that combines these classes will be Class D. This will make our future modeling practice very challenging, as very few practical examples are available from the literature to guide us and none of them are related to water management issues. Our preliminary analysis of the reviewed literature reveals that only three papers present models that could classify as D while the majority of papers (14) offer a Class A model ( Figure 6). As mentioned earlier, an important factor that hinders the application and use of hybrid modeling is the lack of software packages that could implement such models in a user-friendly environment. The majority of system dynamics models are created in platforms such as Vensim and Stella. While such platforms contain various analytical tools such as subscription arrays that can represent agents, they do not fully support the object-oriented nature of ABM. This argument is also valid for ABM platforms such as Netlogo where the integration of SD and ABM is almost impractical. The majority of current hybrid applications are developed using AnyLogic software [45], arguably because it provides a relatively user-friendly multimethod simulation modeling environment that supports agent-based, discrete event, and system dynamics modeling (Figure 7). Despite its dominance in the hybrid modeling arena, AnyLogic has important limitations in accounting for different types of agent decision-making processes associated with optimizations (e.g., resource allocation mechanisms) in a hybrid model with multistage uncertainties [55]. This sort of optimization could be critical for a hybrid model of PWM where modeled company agents are to make decisions based upon optimizations that allocate their resources to different investment options (PW treatment technologies, PW-fresh water injection ratios, etc.) at different locations under dynamic uncertainties. On the other hand, there are many programming languages with preexisting libraries, particularly Python for integrating geospatial components and dynamic simulation modeling. However, the computational efficiency, especially in systems as large as PWM systems, is not trivial and still needs to be addressed. How efficiently different libraries or programming languages can process hybrid models differs on the context of the model and still needs more exploration. Current literature (e.g., Anderson et al. [27]) suggests that programming a customized hybrid ABM and SD model is also not trivial; thus, there is still a need for further research.

Conclusions
In this paper, we have explored the necessity and applicability of a hybrid dynamic simulation model that would address produced water management (PWM) issues and that could be applied to other complex natural resource management issues. We first developed a conceptual diagram to identify cross-scale feedback mechanisms that are in play in an integrated PWM system to see if the cost of hybrid modeling could be justified by the level of complexity involved. Our conceptualization was based on the literature, the formal data of produced water in New Mexico, and some preliminary interviews with subject matter experts. This conceptualization revealed that hybrid modeling could add value to our understanding, and as a consequence, probably would provide better policy advice. However, the amount of value that this effort provides compared to its costs is still an open question. Although the hybrid modeling exercise could yield theoretical advantages, the potential benefits from increased understanding and superior policy advice may not compensate for the additional costs of introducing more complexity into the modeling effort [56].
To select the best modeling approach for PWM modeling, we reviewed the literature of hybrid modeling in the second phase of the project. The goal was to provide some useful guidelines for modelers who would like to work in this area. Our initial exploratory review revealed that a combination of system dynamics (SD) and agent-based (AB) modeling could be necessary and sufficient for the purpose of comprehensive PWM modeling. Therefore, in our next step, we focused exclusively on a systematic review of the SD-AB hybrid modeling literature. We used the Web of Science for our systematic search.
Our literature review indicated that despite its currently small size, the SD-AB hybrid modeling realm is a growing area of research. Seventy-seven papers were found to be useful with respect to the explication and development of hybrid modeling. Among these papers, 31 provided detailed explanations of how this kind of modeling could be performed. However, only one paper was related to coupled water-energy issues. We also found that the majority of the models presented in these papers were developed using the AnyLogic modeling environment. Although it is considered as the most powerful hybrid modeling software, AnyLogic has its own limitations, which underscores the fact that hybrid modeling is still a very challenging practice and in an embryonic stage of development.
To provide a more meaningful guideline for hybrid modelers of PWM, we classified the current state-of-the-art hybrid modeling practices into four classes, A, B, C, and D. Class A is the simplest form of modeling wherein a set of SD modules talk to a set of AB modules using a communication protocol. Class B involves AB models with agents that consist of SD models. Class C involves SD models with variables that are driven by AB rules. Finally, Class D is a combination of Classes B and C, where the structure is very flexible with mixed hierarchical design. We tentatively concluded that a comprehensive PWM problem is likely to require a Class D modeling approach.
We found confirmation that the selection of a modeling approach depends strongly on the purpose of the modeling. Here, we assumed the modeling goal is to provide a complete picture of a comprehensive PWM effect on a regional water budget. Projects with a narrower focus should first consider using simpler approaches, such as traditional SD or AB modeling. In this context, perhaps a minor modification of the standard models might capture sufficient richness without the difficulties and expenses of a highly complex numerical simulation approach.
The findings of this paper will inform our future research, which will include comprehensive systems modeling of produced water management and its impacts on regional water budgets. This paper will also guide hybrid SD-AB modeling in other domains where multiple levels of dynamic interactions are of significance.

Conflicts of Interest:
The authors declare no conflict of interest.