Evaluation of Energy Recovery Potential by Anaerobic Digestion and Dark Fermentation of Residual Biomass in Colombia

This study provides the first overview in Colombia on energy recovery potential by anaerobic digestion (AD) and dark fermentation (DF) of three different residual biomasses: coffee mucilage (CFM), cocoa mucilage (CCM), and swine manure (SM). First, AD and DF models were developed based on the ADM1 model. Then, simulated biogas production yields were compared to experimental data to validate the models. The results of comparative simulations indicate that energy recovery potentials from biogas for the different Colombian departments range from 148 to 48,990 toe, according to the local production amounts of CFM, CCM, and SM in 2017. The study provides crucial information that can be used to assess the best design, operation mode, and locations of AD and DF plants in Colombia. The results indicate that biogas production performances and energy recovery yields improve by increasing CFM/SM and/or CCM/SM ratios of the feed, and by increasing organic load from 2 to 26 gCOD∙l−1.


INTRODUCTION
Biomass is considered a renewable source for several production processes such as biofuel and metabolite production (Osman et al., 2019;IEA, 2020), and, consequently, an increasing number of international studies have been investigating the valorization potential of different biomass types during the last two decades. Furthermore, legislation restrictions have led researchers to focus on developing sustainable biomass management and treatment (Gonzalez-Salazar et al., 2014;Rajesh Banu et al., 2020). Since 2010, in Colombia, the National Administrative Department of Statistics (DANE) has reported an increase in different agro-industrial activities (DANE, 2020), which has improved the income and social benefits of the national economy. According to DANE data (2018), coffee production increased from 0.64 to 0.86 Mton, cocoa production from 34.9 to 101 kton, and pork production from 0.23 to 0.41 Mton between 2011 and 2018 (Supplementary Figure S1), resulting in increased production of residual biomass from these agro-industrial sectors (e.g. coffee and cocoa mucilage, and swine manure). This residual biomass requires treatment processes for biological stabilization and/or final disposal, which represent an over cost for the agro-industry.
In recent years, biogas production has been identified as a sustainable process to treat and valorize residual biomass in Colombia (Hernández et al., 2014;Hernandez et al., 2018;Sagastume Gutiérrez et al., 2020). In this context, coffee and cocoa mucilage usually have high contents in carbohydrates, including dextrose, cellulose, hemicellulose, and others. Carbohydrates are fundamental compounds in the processes of dark fermentation (DF) and anaerobic digestion (AD) for biochemical conversions such as acidogenesis, acetogenesis, methanogenesis and fermentation, thus leading to the production of biogas and other byproducts (e.g. volatile fatty acids, alcohols) (Redgwell et al., 2003;Battista et al., 2016;Okiyama et al., 2017). Meanwhile, swine manure usually presents high alkalinity and high contents of proteins and nutrients (e.g. phosphorus and nitrogen), which may improve pH control and decrease the ratio of carbon to nitrogen (C/N) in co-digestion processes (Astals et al., 2011(Astals et al., , 2015Ye et al., 2013;Hernández et al., 2014;Neshat et al., 2017;Wang et al., 2018).
The simulation of bioprocesses presents new advantages due to the computational development that contributes to complex calculations. Indeed, process simulation represents a powerful tool for industry and researchers because it evaluates real scenarios with cheaper costs (Astals et al., 2011;Rajendran et al., 2014aRajendran et al., , 2014b. The Anaerobic Digestion Model No 1 (ADM1) is one of the most referenced models in the literature for the description of anaerobic processes (Batstone et al., 2002). It has been used as a base model for the development of most of the latest AD and DF models (Penumathsa et al., 2008;Antonopoulou et al., 2012;Saady, 2013;Rajendran et al., 2014a;Alexandropoulou et al., 2018). Although studies have been carried out to assess the energy potential of a large number of different substrates, in the literature to date there is no systematic and complete study of AD and DF for coffee and cocoa mucilage or possible blends with swine manure.
This study aims to provide the first overview in Colombia on energy recovery potential by AD and DF of coffee mucilage, cocoa mucilage, and swine manure. The main objectives are 1) to develop the AD and DF models for the simulation of the bioprocesses, 2) to identify the most suitable blends and process schemes for biogas production (e.g. AD, DF, and DF followed by AD), and 3) to evaluate the energy recovery potential of residual biomass in Colombia. First, AD and DF process models were developed based on the ADM1 model. Conversion reactions and yields, inhibition modules, and kinetic constants from the literature were included in the models, and simulations were performed at mesophilic temperature conditions (35°C) by using the commercial software Aspen Plus. The biogas production yields obtained by simulations were compared with experimental results of laboratory tests to validate the models. Finally, the developed models of AD and DF were used to assess the energy recovery potentials of the different Colombian departments based on the local availability of residual biomass. The energy recovery potentials using a comparative simulation approach between three biogas production processes were assessed to identify the most suitable scenarios for residual biomass valorization in Colombia.

Model Development and Simulation Design
Simulation of the biological processes of AD and DF was performed in Aspen Plus V9 software, the Soave-Redlich-Kwong (SRK) method was chosen as it correlates and calculates the mole fractions and activity coefficients. In addition, the SRK model is suitable for vapor-liquid equilibrium (VLE) (Madeira et al., 2017). The starting point for the simulation of AD was the ADM1 model, which relates the stoichiometric equations for the stages of hydrolysis, acidogenesis, acetogenesis, and methanogenesis (Serrano, 2010). The same methodology was also used to simulate the stages of hydrolysis, amino acid degradation, and fermentation for DF, as proposed by Serrano (2011).
For each simulation in Aspen Plus, carbohydrates (e.g. dextrose), complex carbohydrates (e.g. cellulose, hemicellulose), fats (triolein, tripalmitin), proteins, and water mass fractions of the feed substrate were included in the reactor inlet stream. A maximum of 3 U were used in Aspen Plus to simulate all stages of the biological processes ( Figure 1).
The first unit is a stoichiometric reactor (RStoic) that simulates the hydrolysis stage for both the AD and DF processes. It calculates the mass flows for carbohydrates, proteins, and fats based on the hydrolysis reactions and fractional conversions (Supplementary Table S1) (Rajendran et al., 2014a).
The second unit is a continuous stirred tank reactor (CSTR) for both the AD and DF processes. The CSTR was selected following the most common operating mode of AD and DF reactors in the literature (Barca et al., 2015;Adarme et al., 2017;Neshat et al., 2017;Rajesh Banu et al., 2020). The CSTR simulates the processes of amino acid degradation, acidogenesis, acetogenesis and methanogenesis for AD, and amino acid degradation and fermentation for DF ( Figure 1). A set of 42 reactions was used to describe the main metabolic pathways associated with the bioconversions for the AD and DF processes (Supplementary Table S2) (Rajendran et al., 2014a). The compounds in Table S2 come from the feed substrates and/or from the hydrolysis stage (Supplementary Table S1). All the kinetic reactions are supposed to follow first-order kinetics (Angelidaki et al., 1999). Simulations were performed according to the process simulation model (PSM) proposed by (Rajendran et al., 2014a). PSM model is based on the ADM1, and it includes ammonia inhibition on the calculations of the kinetic constants for the conversion reactions, thus following the model proposed by Angelidaki et al. (1993) and Angelidaki et al. (1999). The model also takes into account variables such as pH and temperature to the calculations of the kinetic constants by correlation functions presented by (Siegrist et al., 1993;Vavilin et al., 1994;Serrano, 2010). All these correlation functions were introduced in the CSTR unit through calculator blocks available in Aspen Plus. Each calculator block had a FORTRAN code to determine the kinetic constants of the conversion reactions. In addition, for DF process a calculator block was integrated into the model to establish the restriction of hydrogen gas production yield from carbohydrates within the values 1-4 mol H 2 ·mol hexose −1 , where 4 is the maximum stoichiometric yield (Alexandropoulou et al., 2018).
The third unit was developed and applied for the DF process to simulate the homoacetogenic stage, which involves the consumption of hydrogen (H 2 ) and carbon dioxide (CO 2 ) gas (Eq. 1) (Saady, 2013).
H 2 consumption by homoacetogenesis usually increases by increasing the mass transfer limitation of H 2 from the liquid to the gas phase (Saady, 2013). In this study, the minimum conversion rates reported in the literature (Saady, 2013) were considered to simulate low mass transfer limitation conditions in the DF reactors. According to this, homoacetogenesis is carried out in an RStoic block with a conversion factor of 0.1 ± 0.1, as a function of the H 2 and CO 2 composition of the biogas produced during the fermentation step. The methanogenic stage controls the H 2 concentration as H 2 is consumed by hydrogentrophic methanogenesis (pathway n°39, Supplementary Table S2). Therefore, the homoacetogenic stage is discarded for the AD process.

Substrate Characterisation
The residual biomasses selected for this study were coffee mucilage (CFM), cocoa mucilage (CCM) and swine manure (SM), which are related to the most relevant agro-industrial activities in Colombia. Raw samples of CCM, CFM and SM were collected from the departments of Santander and Cundinamarca from the farms on a medium production scale. Then, the total volatile solids (TVS) and total organic load (COD) of the samples were determined by standard methods (APHA, 2005). Total carbohydrates, cellulose, hemicellulose, lignin, and protein were determined by high performance liquid chromatography (HPLC), by using a chromatographer Agilent Model 1200 (Agilent Technology, United States). Total organic carbon (C) was determined by the Walkley Black titration method according to the Colombian technical norm NTC 5167 (2003). Total organic nitrogen (N) and total organic sulfur (S) were determined after acid digestion according to the Kjeldahl method and to the barium sulfate precipitation method, respectively. Then, oxygen (O) and hydrogen (H) contents were estimated stoichiometrically. The results of biomass characterization are summarized in Table 1. These results are in good agreement with those provided by the online database of the Netherlands Energy Research Center (Phyllis2, 2020).
A new category for the substrates had to be created in Aspen Plus due to the database not having this specific kind of complex biomaterials. Therefore, CCM, CFM and SM were defined in the database of Aspen Plus through a set of specific properties, as summarized in Table 1. In the inlet streams for AD and DF models, total fats were divided in equal fractions between triolein and tripalmitin.

Model Validation
The validation of the AD and DF models was carried out by comparison of simulated results to those of batch laboratory tests that aimed at evaluating the biogas production potential of different blends of CFM, CCM, and SM (Cabeza et al., 2016;  Table 2), based on three independent parameters, was used to determine the number of blends for the comparative experiments: 1) the dried mass ratio between CFM and CCM (CFM:CCM ratio), 2) the mass ratio of C to N (C/N ratio), and 3) the total organic load. For each parameter, three levels were considered: the CFM:CCM ratios were 3:1, 1:1, and 1:3; the C/N ratios were 25, 35, and 45; and the organic loads, based on chemical oxygen demand (COD), were 2, 5, and 8 gCOD·l −1 . SM was used as N rich substrate to set the C/N ratio of each blend. Then, distilled water was added to set the COD concentrations of each blend.
The experiments were carried out in 250 ml bottles with a total working volume of 200 ml. The reactors were operated under batch conditions and at mesophilic temperature (35°C). The inoculum used for AD and DF experiments was obtained from the biodigester of the wastewater treatment plant of a dairy company (ALPINA, Sopó, Cundinamarca). Before DF tests, a heat shock treatment of the inoculum (30 min at 100°C) was applied to select sporeforming hydrogen-producing bacteria and to destroy hydrogen-consuming archaea that are usually non-sporeforming (Barca et al., 2015). The measurement of methane (CH 4 ) and hydrogen gas (H 2 ) was taken daily by volume displacement. NaOH solution (pH > 10) was used as a carbon dioxide gas (CO 2 ) trap and, therefore, it is assumed that the displaced gas only corresponds to CH 4 in process AD and H 2 in process DF. The duration of the tests was between 10 and 20 days for AD, and between 5 and 12 days for DF experiments. The experimental biogas yields (ml CH 4 ·g COD −1 and ml H 2 ·g COD −1 for AD and DF, respectively) were calculated based on the cumulated volume of biogas produced from each experiment. All the tests were performed in triplicate and average biogas yields were determined for each experimental condition. As the duration of the tests varied between the different experiments, it was used as the input value of hydraulic retention time for each simulation. All simulations were performed at 35°C and 0.1 MPa. The simulated biogas production yields (ml biogas·g COD −1 ) were calculated based on the mass amount of CH 4 , H 2 , and CO 2 , and on the density of each compound at standard conditions (25°C and 0.1 MPa).

Process Scheme Simulation
The biogas production performances of each substrate blend ( Table 2) were further simulated at three different process schemes: anaerobic digestion (AD), dark fermentation (DF), and dark fermentation followed by anaerobic digestion (ADF). The main objective was to identify the most suitable blends and process schemes, for the production of biogas or a specific gas (CH 4 , H 2 , CO 2 ), using a comparative simulation approach, thus taking into account relevant factors such as the C/N ratio and the organic load of the feed. Indeed, previous studies have already shown that these factors may significantly affect biogas production performance (Zhang et al., 2014;Xie et al., 2016;Kovalovszki et al., 2017;Mosquera et al., 2020). All the simulations were performed at 35°C and 0.1 MPa, and according to hydraulic retention times of 20 and 10 days for AD and DF, respectively. Then, the biogas yields (ml biogas·g COD −1 ) were calculated to assess and compare the performances of the different process scenarios.

Evaluation of Energy Recovery Potential in Colombia
The energy recovery potential in Colombia was based on residual biomass production data of 26 departments (out of the 32 Colombian departments) according to the agricultural production reports of the National Department of Statistics, the Colombian Institute of Agriculture, and the National Pork Producers Association (DANE, 2020). A factor of 0.75 was used to estimate the CFM and CCM effectively available (Şenol, 2019), thus assuming 25% of losses in their collection. In addition, the SM availability was calculated as 75% of the total SM from technical farms (ICA, 2018). Once the availability of CFM, CCM, and SM was established, the annual total was divided as daily flow and the ADF process was simulated for each department. The inlet organic load of 26 gCOD·l −1 was used as an input value for all the department simulations based on the COD concentration of the substrates ( Table 1) and their total available amounts for each department ( Table 3). Temperatures and heat duties, from the results of ADF simulations, were analyzed and gathered to perform the heat integration procedure, according to the pinch method, thus considering the minimum approach temperature ΔTmin of 3°C. Heat transfer from the effluents was integrated into the calculations by considering their cooling down to a temperature of 25°C (Abdou Alio et al., 2021). For each department, the total daily masses of inlet substrate M substrate (CCM, CFM, SM) (kg·day −1 ), outlet biogas M biogas (NH 3 , CO 2 , H 2 , CH 4 , H 2 S) (kg·day −1 ), and outlet digestate M digestate (propionate, butyrate, ethanol, benzene, acetate, furfural, residual carbohydrates, and proteins) (kg·day −1 ), were calculated based on the results from simulations. Then, total input power P input (MJ·day −1 ), substrate power P substrate (MJ·day −1 ), biogas power P biogas (MJ·day −1 ), and digestate power P digestate (MJ·day −1 ) were determined by Eqs 2-5, respectively, where LHV is the related lower heating value (MJ/kg) and P process (MJ·day −1 ) represents the energy flow required to heat the feedstock and maintain a constant process temperature in the reactors.  Frontiers in Energy Research | www.frontiersin.org June 2021 | Volume 9 | Article 690161 5 P biogas M biogas · LHV biogas (4) The energy recovery potentials E (toe) from biogas and digestate by ADF process were calculated using Eqs 6, 7, respectively, where d represents the days of the year (365) and f E represents the conversion factor ton of oil equivalent per Joule (2.388ˑ10 −5 toe/MJ).
The energy recovery yields for biogas η biogas (%), digestate η digestate (%), and total energy recovery η total (%) were calculated by Eqs 8-10, respectively. tη biogas P biogas P input (8) η total P biogas + P digestate P input Finally, a sensitivity analysis of the ADF model was made to assess the effect of nine different organic loads (from 2 to 26 gCOD·l −1 ) on the total energy recovery yields η total (%) for the five departments with the highest production of coffee, cocoa and pork (Antioquia, Boyacá, Cundinamarca, Meta, and Santander). Figure 2 shows the correlation between simulated and experimental biogas production yields for the different blends tested. Overall, the results show an average ratio between simulated and experimental yields (y sim /y exp ) of 1.05 ± 0.24 for AD and 1.17 ± 0.28 for DF (average values ± standard deviations), respectively, thus indicating a good general agreement between the data.

Model Validation
With the only exceptions of blends 1, 2, and 13, the ratio y sim / y exp for AD varied from 0.79 to 1.26 (Figure 2A), thus indicating that the AD model gives a good estimation of CH 4 production yields. The lowest differences were observed for blends 3, 7, and 10 (y sim /y exp ratios of 0.92, 1.00, and 1.09, respectively), whereas the highest differences were observed for blends 2 and 13 (y sim / y exp ratios of 0.70 and 1.79, respectively). The results indicate that the differences between simulated and experimental data for AD increase by decreasing the organic loads and/or by increasing the SM content of the blends. Higher organic loads probably resulted in dominant organic compounds in the blends used for the experiments, thus leading to more stable and more predictable conversion pathways. Moreover, variability in the single substrate composition used for the different experiments, especially for SM, may further explain the differences between simulated and experimental yields.
With the only exceptions of blends 4, 7, and 12, the ratio y sim / y exp for DF varied from 0.77 to 1.28 ( Figure 2B), thus confirming the validity of the DF model. The lowest differences were observed for blends 3 and 6 (y sim /y exp ratios of 1.03 and 1.06, respectively), whereas the highest differences were observed for blends 4, 7, and 12 (y sim /y exp ratios of 1.72, 1,94, and 1.43, respectively), thus indicating an overestimation of H 2 production yields compared to the experimental data. The results of DF seem to indicate that the differences between simulated and experimental yields increase by increasing the CCM and/or CFM loads of the blends. Most probably, incomplete hydrolysis of complex carbohydrates (e.g. cellulose, hemicellulose) may account for the lower experimental yields for blends at high loads of CCM and/or CFM (blends 4 and 7, respectively).
The results of the validation model suggest that variability on single substrate composition between the different experiments and incomplete hydrolysis of carbohydrates may account for the differences between simulated and experimental yields. However, using the same AD and DF models for evaluating the potential of biogas production potential from model substrates may overcome these technical limitations of laboratory experiments and thus allow comparison of the potential of energy recovery potentials from different blends and different process schemes.

Simulated Biogas Production Performances
The total biogas production yields and compositions obtained by the simulations of AD, DF, and DF followed by AD (ADF) processes, for the different substrate blends ( Table 2), are summarized in Figure 3.
The biogas production yields for DF varied between 147 and 405 ml biogas·gCOD −1 , and the H 2 and CO 2 content of biogas varied between 65-70% and 30-35%, respectively ( Figure 3C). These results are in good agreement with the theoretical compositions of biogas from the fermentative pathways used for DF simulation (n°39 and 40, Supplementary Table S2) (Gomes et al., 2015). The lowest biogas production yield for DF was observed for blend 9, where the dominant substrate is SM (more than 80% based on total organic load) and COD concentration is 2 gCOD·l −1 . Instead, the highest biogas production yield for DF was observed for blend 7, where the dominant substrate is CFM and COD concentration is 5 gCOD·l −1 . These differences in biogas production yields are probably primarily related to the different compositions of proteins, fats, and carbohydrates considered for CFM and SM ( Table 1). Higher carbohydrate content of CFM resulted in higher H 2 production yields for blend 7, as carbohydrates represent the most suitable substrate for DF (pathways n°39 and 40, Supplementary Table S2), whereas higher protein content of SM resulted to lower C/N ratios that may inhibit DF as described by the ammonia inhibition modules integrated into the model (Madeira et al., 2017;Serrano, 2011). These results are in good agreement with the findings of previous experimental studies (Zhu et al., 2008;Ghimire et al., 2015), which have indicated that low C/N ratio and high concentration of ammonium ion (NH 4 + ) inhibit H 2 production during DF when SM is used as the only substrate, with maximum production yields of 3.65 ml H 2 ·VS −1 . Indeed, higher production yields (from 68 to 652 ml H 2 ·gCOD −1 ) were obtained by co-fermentation of substrate blends (e.g. swine manure, sewage sludge, food waste, agro-industrial waste) (Kim et al., 2004;Perera et al., 2012;Hernández et al., 2014).
The biogas production yields for AD varied between 295 and 459 ml biogas·gCOD −1 , thus showing a lower variability according to the different blends rather than DF. The CH 4 , CO 2 and H 2 biogas components were 57-70%, 29-42%, and 0.3-0.4%, respectively ( Figure 3B). Higher yields of biogas production were observed for blends containing higher amounts of CFM and CCM rather than SM. Overall, biogas production yields obtained by simulations are within the ranges reported in the literature (61-650 ml CH 4 ·gCOD −1 ). This large variability in experimental yields is mainly attributed to the different organic substrates and operating parameters (e.g. retention time, temperature, inoculum) between the different studies (Garfí et al., 2011;Hernández and Rodríguez, 2013;Astals et al., 2015). The blends 12 and 13, with equal CFM and CCM contents, showed the highest biogas yields (433 and 459 ml biogas·gCOD −1 , respectively) with the highest CH 4 content (70 and 69%, respectively). Therefore, they appear to be the most suitable blends for AD. These results are consistent with the experimental results of Wang et al. (Wang et al., 2018), who indicated that biogas production performances of AD are significantly affected by the ratios of swine manure, corn stove, and cucumber residues in the co-digested mixture.
The biogas production yields for DF followed by AD (ADF) varied between 602 and 864 ml biogas·gCOD −1 , and the CH 4 , CO 2 , and H 2 content of biogas varied between 25-51%, 30-35%, and 17-40%, respectively ( Figure 3D). The highest total biogas production yields were observed for blends 7 and 8 (863 and 864 ml biogas·gCOD −1 , respectively), where the CFM and CCM are above 65% of the total organic load. Additionally, one of the substrates, CFM or CCM, dominates over the others. The lowest yields (<160 ml biogas·gCOD −1 ) were obtained for blends 5, 6, 9, and 10, where SM is greater than 80% over the total organic load. However, despite the lowest biogas production yields, these blends also show the highest CH 4 contents in biogas. This appears to confirm that SM is more suitable for CH 4 than for H 2 production, most probably because it has a lower content in carbohydrates compared to the other substrates. While, blends with an SM load equal or lower than those of CFM and CCM (blends 7, 8, 11, 12, and 13) generate biogas with CH 4 and H 2 contents comprised between 28 and 40%. ADF showed higher yields of total biogas production compared to AD and DF, this for all the substrate blends studied. As shown in Figure 3A, total biogas production yields for ADF were 1.3-2.3 times higher than DF, and 1.5-2.4 times higher than AD, thus indicating a greater conversion efficiency of organic substrate to biogas. This is coherent with the assumptions of the ADM1 model that was used for the simulations, as organic products from DF (e.g. volatile fatty acids) represent viable substrates for methanogenesis in the AD process (Hernández et al., 2014), thus increasing the total biogas production yields (Yeshanew et al., 2018).
The results of the simulations confirm the effectiveness of codigestion of CFM, CCM, and SM, in a two-step process of DF followed by AD compared to the individual process steps. The higher content in proteins, fats, and hemicellulose of SM supports the production of fatty acids, amino acids, and carbohydrates by hydrolysis in the DF step (Supplementary Table S1). These byproducts are further converted and used as substrates for CH 4 production by acidogenesis, acetogenesis, and methanogenesis in the AD step (Supplementary Table S2). Instead, the higher content in carbohydrates of CFM and CCM conducted to the production of acetate and H 2 by fermentation in the DF step (Supplementary Table S2), which are used as substrates for CH 4 production by acetotrophic and hydrogenotrophic methanogenesis in the AD step (pathways 37  Table S2). Moreover, the simulation results show that different organic loadings do not seem to have a significant impact on the production yields and biogas composition for the different blends and processes studied. This behavior was related to all blends and processes simulated with the same stoichiometric ratios, considering residence times long enough to reach conversion equilibrium and without considering excess substrate inhibition. However, despite showing similar production yields, higher organic loads may give higher flow rates of biogas production in real operating conditions, as more organic substrate is available for biogas production per unit of time (Barca et al., 2015).

Evaluation of Energy Recovery Potential in Colombia
Among the 26 departments considered in this study, 17 report coffee production more than twice higher than cocoa production (DANE, 2020). Only eight departments (Bolívar, Córdoba, La Guajira, Meta, Arauca, Putumayo, and Vichada) show higher cocoa production than coffee production, and their total coffee production amounts are relatively low (<210 tons of coffee/year) compared to the other departments (>2000 tons of coffee/year). According to DANE data (2017), the departments of Huila, Cauca, Tolima, Caldas, and Antioquia generates 73% of the Colombian coffee production and they are therefore considered as the main Colombian coffee region. Additionally, their coffee production is more than 15 times higher than their cocoa production, thus indicating a much lower availability of CCM rather than CFM (Supplementary Figures S2, S3). Table 3 summarizes the CFM, CCM, and SM availability data for 26 Colombian departments in 2017. Overall, the data indicates that CFM production dominates over CCM production in the center and south of the Andean region (Antioquia, Huila, Tolima, Caldas, Risaralda, Quindío), in the Pacific region (Chocó, Valle del Cauca, Nariño, Cauca), and more generally in the south-west of Colombia, whereas CCM production dominates over CFM production in the Orinoquía region (Arauca, Casanare, Meta, Vichada), in the Caribbean region (Atlántico, Bolívar, Córdoba, Sucre, La Guajira Cesar, Magdalena), in the north of the Andean region (Boyacá, Norte de Santander, Cundinamarca, Santander), and more generally in the north-east of Colombia. The Amazon region does not show significant productions of CFM and CCM, with the only exception of the department of Caquetá. SM availability per department in 2017 varies from 3,347 ton (Vichada) to 1084000 ton (Antioquia). The mass ratio of available SM over the sum of available CFM, CCM, and SM (SM/total) varies from 0.878 (Santander) to 1.000 (Atlántico), thus indicating a larger availability of SM compared to CFM and CCM for all the departments.
For each department, the potential recovery of biogas and digestate was evaluated by comparative ADF simulations according to the total available amounts of CFM, CCM, and SM in 2017, and considering for all the simulations an inlet COD concentration of 26 g COD·l −1 . As shown in Table 4, the total mass of biogas recovered per department varied from 307.6 ton (Vichada) to 102200 ton (Antioquia), whereas the total mass of digestate recovered per department varied from 146.3 ton (Vichada) to 47,580 ton (Antioquia). The main components of biogas were CO 2 and CH 4 , with an evaluated production potential ranging between 176.6-58,960 ton and 105.7-34,950 ton, respectively, whereas the main components of digestate were ethyl cianoacetate and propionate, with an evaluated production potential ranging between 11.19-3,693 ton and 10.08-3,295 ton, respectively. Table 5 summarizes the energy recovery potential (toe) from the biogas and digestate obtained for each department based on production data of 2017, as well as the energy recovery yields for biogas η biogas (%), digestate η digestate (%), and the total energy recovery yields η total (%) (Equations 8-10). The highest energy recovery potentials from biogas were observed for the departments of Antioquia, Córdoba,and Cundinamarca (48,990,10,160,and 13,460 toe,respectively), which are the departments with the highest cumulated productions of CFM, CCM, and SM. Instead, the lowest energy recovery potentials from biogas were observed for the departments with the lowest cumulated productions of CFM, CCM, and SM (Putumayo and Vichada, 866 and 148 toe, respectively). According to the energy consumption data in 2017 (IEA, 2020), the total energy recovery potential from the biogas (151.6 ktoe) could supply up to 70% of the energy demand for agriculture/forestry activities in Colombia (215 ktoe for 2017).
As shown in Table 5, the energy recovery yields for biogas η biogas (%) varied between 34.3 and 39.3%. The highest η biogas (%) were observed for the departments of Cauca,Huila,Santander,and Tolima (36.3,36.3,39.3,and 36.5%,respectively), which are the departments with the highest CFM/SM and/or CCM/SM production ratios. On the contrary, the departments with lower CFM/SM and/or CCM/SM production ratios (Atlántico, Córdoba, Sucre, and Vichada) usually present lower energy recovery yields for biogas η biogas (%) and higher energy recovery yields for digestate η digestate (%). This seems to confirm that although SM is less effective than CFM and CCM in biogas production, it is more effective for recovering other by-products from digestate (e.g. propionate, butyrate, ethanol, benzene, acetate, and furfural). As shown in Table 5, the total energy recovery yields η total (%), which include biogas and digestate energy (Equation 10), varied between 48.6% (Atlántico, Sucre, and Vichada) to 52.5% (Santander).
The results show higher total energy recovery yields η total (%) for the departments with higher CFM/SM and/or CCM/SM production ratios, confirming that CFM and CCM are more effective than SM for energy recovery by ADF processes. It should be noted, however, that since this is a theoretical study of energy recovery potential, factors such as seasonal variations in agroindustrial production and variability in substrate composition between different production sites due to breeding and cropping techniques, as well as soil and climate characteristics, are not considered.

Effect of Organic Load on Energy Recovery Yield
One of the impact factors on biogas production yields is the input load (Jiang et al., 2020). Figure 4 shows the evolution of Frontiers in Energy Research | www.frontiersin.org June 2021 | Volume 9 | Article 690161 8 total energy recovery yield η total (%) as a function of inlet organic load for the five departments with the highest production of CFM, CCM, and SM (Antioquia, Boyacá, Cundinamarca, Meta, and Santander). The η total (%) values in Figure 4 were obtained by comparative ADF simulations according to the CFM, CCM, and SM production data in 2017. The results indicate that η total (%) increases from 10.2-11.4% to 48.4-51.4% by increasing inlet COD from 2 to 26 gCOD·l −1 for all the departments. This was primarily related to the lower energy consumption for the initial heating of the feed, since the higher the COD concentration, the lower the feed flow to be heated. However, the rate of η total (%) increase appears to decrease by increasing organic load, thus suggesting an optimal inlet COD concentration to obtain higher η total (%) while maintaining acceptable yields of substrate conversion. Moreover, increasing the inlet COD concentration over the average value of the available raw mixtures (26 gCOD·l −1 ) would require further energy consumption (e.g. drying and membrane processes) that may affect the energy recovery efficiency of the full process. Among the departments in Figure 4, the department of Santander shows the highest η total (%) for all the inlet COD concentrations. This was probably because Santander has the highest ratios (CFM + CCM)/SM and C/N, thus resulting in higher biogas production performances compared to the other departments.

Closing Remarks and Perspectives
In recent years the Colombian economy has come to strongly depend on hydrocarbon production from non-renewable natural resources. In April 2020, 16.6 million barrels of crude oil were exported, which represented a 19.7% drop compared to April 2019 (DANE, 2020). Biofuel production from alternative renewable sources such as agro-industrial wastes may represent a sustainable way to overcome the shortage of natural deposits of fossil carbon. During the last 5 yr, the increase in pork breeding, coffee, and cocoa production (DANE, 2020) suggests an interesting potential market of SM, CFM, and CCM valorization as organic substrates for biogas production. However, the fluctuation in production of these residual biomasses may represent a weakness for perspective evaluations, and therefore a regular updating of inlet data for the models is recommended to better follow future production dynamics. It should also be noticed that most of the coffee, cocoa, and pork production activities are located in rural areas, and, currently, more than 20% of transport roads in Colombia are without pavement (DANE, 2020). This may significantly hinder the collection of co-substrates from different production areas, and it may also generate supplementary transport costs.  The bioprocess simulation approach developed in this study may be followed to evaluate energy recovery potential from other types of locally available substrates (IEA, 2020), thus showing flexibility of adaptation to local contexts. Furthermore, the integration of thermochemical (e.g. gasification, pyrolysis, hydrothermal liquefaction) and separation (e.g. membrane, extraction, distillation) processes may allow the recovery of further by-products (e.g. nutrients, biochar, amino acids, organic acids) from the effluents of AD and DF processes (Cebreiros et al., 2017;Arora et al., 2018;Rajesh Banu et al., 2020;Uddin et al., 2021). Therefore, an accurate life cycle analysis (LCA), integrating energy balances and operating costs for each unit of the process, would be required to identify the best design, location, and size of the biorefinery plants in Colombia.

CONCLUSION
This study contributes to the development of AD and DF models to assess biogas production performances from residual biomass. The results indicate that the local availability of different types of residual biomass represents the most influential parameter in assessing the energy recovery potential in Colombia. The energy recovery potential in Colombia for 2017 is 155.1 ktoe according to the local production amounts of CFM, CCM, and SM. Overall, energy recovery yields improve by increasing CFM/SM and/or CCM/SM ratios of the feed, and by increasing organic load from 2 to 26 gCOD·l −1 .

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
MA and J-HF contributed to process design and simulation and PINCH analysis. MH, contributed to experimental design of the study, J-HF contributed design of the energy analysis. MA organized the database and wrote the first draft of the manuscript. MA and CB wrote sections of the manuscript. All authors contributed to manuscript revision, read, and approved the submitted version.