Climate resilient interconnected infrastructure: Co-optimization of energy systems and urban morphology

Co-optimization of urban morphology and distributed energy systems is key to curb energy consumption and optimally exploit renewable energy in cities. Currently available optimization techniques focus on either buildings or energy systems, mostly neglecting the impact of their interactions, which limits the renewable energy integration and robustness of the energy infrastructure; particularly in extreme weather conditions. To move beyond the current state-of-the-art, this study proposes a novel methodology to optimize urban energy systems as interconnected urban infrastructures affected by urban morphology. A set of urban morphologies representing twenty distinct neighborhoods is generated based on fifteen influencing parameters. The energy performance of each urban morphology is assessed and optimized for typical and extreme warm and cold weather datasets in three time periods from 2010 to 2039, 2040 to 2069, and 2070 to 2099 for Athens, Greece. Pareto optimization is conducted to generate an optimal energy system and urban morphology. The results show that a thus optimized urban morphology can reduce the levelized cost for energy infrastructure by up to 30%. The study reveals further that the current building form and urban density of the modelled neighborhoods will lead to an increase in the energy demand by 10% and 27% respectively. Furthermore, extreme climate conditions will increase energy demand by 20%, which will lead to an increment in the levelized cost of energy infrastructure by 40%. Finally, it is shown that co-optimization of both urban morphology and energy system will guarantee climate resilience of urban energy systems with a minimum investment.


Introduction
Due to rapid urbanization [1], cities are witnessing a drastic growth. The number of cities with a population of five hundred thousand to one million is currently estimated to increase from 598 (in 2016) to 710 (in 2030) [2]. This massive urban development has significantly changed the morphology of cities, as new multi-functional urban areas appear within and beyond the borderline of megacities. Cities and urban areas are characterized by their high energy density and heterogeneity in energy use profiles [3]. They accommodate around 50% of the world's population [4] and are responsible for two-thirds of the global primary energy consumption, inducing 71% of global direct energy-related greenhouse gas (GHG) emissions [5,6]. It is projected that due to population growth in urban areas, 68% of the world's population will live in urban areas by 2050 [4]. This, together with climate change and economic growth will place enormous pressure on material and energy resources [7].
It is widely known that energy demand in urban areas is highly affected by climate variations [8]]. Increasing the share of renewable energy generation based on distributed energy sources, such as solar and wind energy, will increase the dependence of energy supply on climate and can change the roles and relations in the energy market [9]. Climate change including extreme events can considerably influence the reliability and resilience of energy systems, affecting different aspects of the energy flow from generation to demand [10]. Reaching a win-win situation where both climate change mitigation and adaptation are fulfilled requires understanding the interactions between variations in climate, energy demand and supply. It is vital to consider the climate dependence of the energy demand and supply and to accommodate their variations in the energy system design when planning to increase the sustainability and resilience of energy systems [11]. The energy conditions are also strongly influenced by the sectoral and technical changes that take place in cities and urban areas, such as increased urbanization and electrification of transport [8]. This means that the energy system is an interconnected urban infrastructure with complex interactions with climate, urban morphology, urbanization, etc.
Urban morphology has a considerable impact both on the energy sustainability and the climate resilience of cities. In the process of urbanization, the main elements of urban morphology, such as urban form (e.g. density, shape, layout, height, etc.), urban function (e.g. functional needs of buildings, size, location, etc.) and urban pattern (e.g. streets, canopies, open spaces) have transformed into far more complex interconnected urban structures [12,13]. The complexity of urban morphology is one of the major drivers of climate variations at urban and local scale [14]. Urban climate can be defined as the interactions between urban morphology elements and climate variables [15]. More specifically, air temperature, wind flow patterns [16], relative humidity [17], and solar radiation [18] are considerably affected by morphological elements. Among other, urban morphology can modify the amount of shading [19] and desirable/undesirable solar radiation in urban areas [20]. This can increase the dependence of buildings on air conditioning systems to maintain thermal comfort [21] as well as the electricity required to provide desired lighting [22]. A clear link also exists between urban morphology and renewable energy potential in urban areas [23,24], which demonstrates the importance of urban morphology from the perspectives of energy efficiency and sustainability.
A number of recent studies have focused on improving the efficiency of urban forms. Natanian et al [25], proposed a framework to optimize generic urban forms in terms of energy demand and spatial daylight autonomy. Ye et al [26], assessed different methods to optimize urban forms to reduce CO 2 emissions. Martins et al [27], introduced a multiobjective optimization framework to optimize solar energy potential in a simplified urban form. Only a few studies focus on the impacts of urban morphology on urban energy system design and renewable energy integration [28]. Perera et al. [29,30] showed that urban density as well as building layout have a notable impact on investment and operation cost, renewable energy integration and autonomy level of urban energy systems. According to their studies, certain building layouts will facilitate higher integration levels of renewable energy technologies while minimizing the cost (of both operation and investment) of the energy systems. However, it is impossible to define an optimal urban layout that can be used universally since it depends on the local seasonal changes in climate, energy demand patterns as well as renewable energy generation potential. Instead, it makes sense to develop methodologies to optimize urban morphology considering energy as well as climate resilience.
Climate change complicates the design problem due to its multidimensional impacts on both the energy and the building sector as well as due to the considerable uncertainties that exist in climate change projections [31]. The sequential method applied in the present state of the art, in which urban planning is performed as a sequential process where decisions related to urban morphology are taken first and followed by energy system designis unsuitable to deal with this level of complexity [32,33]. The thermal performance of urban morphologies can significantly change due to extreme climate events such as heatwaves, ultimately collapsing the energy infrastructure [9,10]. This highlights the importance of harmonizing the urban planning and the energy system optimization process, which demands a major change in methods used for energy system design and urban planning.
Reaching a higher level of energy sustainability and climate resilience in cities requires developing frameworks that concurrently optimize the urban morphology and the energy system while accounting for their interconnected interactions. This study develops such a framework, considering the following objectives: • Developing a co-optimization algorithm that can facilitate layout planning for the building sector as well as optimize energy system design. Understanding the connectivity between building and energy infrastructure plays a vital role in this regard which has not been discussed in the present state of the art. Optimizing both building and energy infrastructure considering 1. Urban form 2. Urban density 3. Energy system design (capacity of energy system components) 4. Operation strategy implemented by the energy system is performed in this study.
• Evaluate the impact of future variations on energy optimized design to assess climate resilience, since future climate variations and extreme climate events can have a notable impact on the energy demand as well as on renewable energy generation, especially considering extreme conditions.
In this work, the focus is mainly on the interactions between the energy system (demand, supply and system design), outdoor climate and urban morphologies. The scope of the study is limited to these three major interacting factors, aiming to avoid increasing uncertainties in the assessment, understand the interactions, investigate the probable correlations and reach viable results. The research paper is organized as follows: Section 2 of the manuscript presents the workflow that links the future climate model with the urban simulation model and the energy system model. Section 3 presents the model used to generate climate files required for the building and energy system design process with the support of regional and global climate models. Section 4 presents an overview of the computational model used for the urban simulation, which generates the energy demand profiles for a set of different urban morphologies. The co-optimization of urban form and the energy system are explained in the latter part of Section 5. Finally, Section 6 is devoted to discussing the results of the study followed by conclusions to pre-sent the major findings of the study in Section 7.

The workflow of the computational model
Urban energy planning is too complex to be handled by a single computational model. Therefore, computational platforms that link several models are often developed to achieve this task. For example, Guan et al. [34] and Mohajeri et al. [28] developed a computational platform linking a GIS module, an urban simulation module and an energy system design module to support energy planning for a Swiss village. A similar workflow is not necessarily adequate to assess the impacts of climate change. To do so, the workflow needs to be able to  Workflow of the study in five major steps including generating urban morphologies (Step 1), and future weather datasets (Step 2), energy simulations (Step 3), assessing renewable energy potential (Step 4), collecting techno-economic data (Step5), and achieving optimal system design and urban morphology (Step 6).
handle large weather data sets and climate uncertainties or to generate representative weather data sets from climate models [35]. Moreover, it should be able to take into account climate variations with fine temporal and spatial resolutions (e.g. hourly variations minimum at the mesoscale and preferably at the urban-scale) and their impacts on both energy demand and supply. These all mean that the workflow should be able to provide a seamless transfer of weather data from climate models to energy models, either its own or representative weather data sets that take into account climate uncertainties as well as both typical and extreme conditions [31]. Therefore, the workflow of the urban planning process is extended in this study to include the inputs from climate models as presented in [10] and shown in Fig. 1.
The workflow presented in this work, first generates representative future weather files as presented in Section 3 and uses them in urban climate modelling and energy simulations. The workflow begins with formulating urban archetypes as described in Step 1, urban simulation is performed to identify the energy demands in the future considering typical and extreme weather conditions. Step 2 presents climate data for the urban simulation conducted in Step 3. Future climate weather data is used to obtain renewable energy generation potentials as presented in Step 4. Technoeconomic data is collected in Step 5, and finally used to arrive at an optimal system design and urban morphology in Step 6. The representative weather data sets, including typical and extreme weather conditions, are used in the urban simulation model to generate energy demand for heating, cooling and electricity for different alternative urban morphologies. The urban energy simulation model considers the weather conditions, interactions among buildings (shading effect), occupancy in buildings, lighting conditions and equipment usage at the building level when computing the energy demand for each morphology. A comprehensive explanation of the urban simulation model is presented in Section 4.
Subsequently, a pool of time-series corresponding to different urban morphologies are used to find the optimal urban morphology. The generated time-series of energy demand for the considered urban morphologies are used to optimize the urban morphology (form and density) along with the energy systems. The energy system optimization model is used to size the energy system, i.e. to define the optimal capacities for dispatchable energy sources, non-dispatchable renewable energy technologies such as wind and solar PV, energy storage, and energy conversion devices such as heat pumps. Usually, this is formulated as a bilevel optimization problem considering dispatch optimization at the primary level and the system configuration optimization at the secondary level. Energy demand is considered as an input parameter to the energy design model, which does not depend on the decision vector. The main change introduced in this study is that we consider urban form and urban density as a decision vectors. As a result, energy demand is subject to vary with the decision vector, which enables a concurrent definition of optimal design and urban form. A comprehensive explanation of the formulation of objective functions is presented in Section 5.

Climate models and developing extreme climate conditions
Global climate models (GCMs) are used to simulate future climate. GCMs are forced by different concentrations of anthropogenic Greenhouse Gas (GHG), depending on the selected Representative Concentration Pathways (RCPs), four of which are defined based on the IPCC Fifth Assessment Report (AR5): RCP 2.6, RCP 4.5, RCP 6.5 and RCP 8.5, labeled after a possible range of radiative forcing values in the year 2100 (2.6, 4.5, 6, and 8.5 W/m2, respectively) [36,37]. GCMs have a coarse spatial resolution of 100-300 km and their outputs cannot be considered weather data [38], which is required for energy simulations (check Fig. 2). Therefore, GCMs are downscaled using two main approaches; dynamical and statistical downscaling [31]. The statistical version only reflects changes in the average weather conditions and underestimates extremes. Therefore, dynamical downscaling (using regional climate models or RCMs) is recommended, which simulates physically consistent weather data across different variables with the temporal and spatial resolutions suitable for energy simulations [39]. The spatial resolution of RCMs is mostly between 10 and 50 km (mesoscale), although there exist a few with finer resolutions, even down to 2.5 km. However, RCMs mostly do not consider the impacts of urban texture and energy flows in the urban scale since this usually needs downscaling to spatial resolutions of less than 1 km. Meanwhile, at the urban-and micro-climate scale, considerable variations may occur to mesoscale weather, amplifying or dampening extreme weather [14] and influencing the energy performance of buildings [9].
In this work, the microclimate around buildings has been simulated using the weather data simulated by RCA4 [40], which is the 4th generation of the Rossby Centre RCM. When doing an impact assessment of climate change, it is not considered appropriate to rely on only a few climate scenarios and short term periods [39]. Multiple climate scenarios (e.g. with different GCMs and RCPs) should be considered for study periods of 20-30 years, which increases the computational load extensively [31]. Nik [31] developed a method to overcome the challenges of big data sets and climate uncertainties. The method is based on synthesizing three sets of representative future weather data sets out of RCMs: Typical Downscaled Year (TDY), Extreme Cold Year (ECY), and Extreme Warm Year (EWY). The three synthesized one-year-long weather data sets represent typical and extreme weather conditions while considering multiple climate scenarios. This has the advantages of shortening the analysis time considerably while still accounting for climate uncertainties and extremes. The application of the method has been compared with other available approaches and weather data sets [41] and verified for several types of simulations and impact assessments [31]. The method has been developed further to quantify the impacts of climate change on urban energy systems [10]. In this work, two groups of data sets have been generated; one for calculating the energy demand based on the method in [31] and one for renewable generation potentials based on the method in [10]. Representative weather data sets have been synthesized for the three 30-year periods of 2010-2039, 2040-2069, and 2070-2099, considering 13 future climate scenarios (five different GCMs and three RCPs: RCP 2.6, RCP 4.5 and RCP 8.5). For more details about synthesizing the weather files, the reader is referred to [31,10].
As mentioned above, RCM outputs, which are meso-climate, can be considered as weather data, reflecting the conditions and variations that affect energy performance of buildings, energy systems and human comfort (Fig. 2). Using mesoscale climate data (with the horizontal resolution of 5-100 km) for energy simulations is widely practiced and accepted by the scientific community. Although buildings and infrastructures are affected by the urban/micro-climate around them, it is very rare to use such high-resolution climate for energy simulations due to lack of access to observed/simulated weather data and their strong case dependency (i.e. valid for the specific cases that the urban/microclimate have been observed or simulated). Therefore, so far the common approach in the energy simulation community has been using meso-climate data. For example, all the available and ready to use weather data sets for energy simulation tools (like the well-known TMY format) are meso-climate [42]. As is reviewed and discussed thoroughly in previous works of the authors [31,41], there are different approaches to generate weather data sets from RCMs to be used energy simulations. Being able to properly count for sub-meso features and variations that occur to weather due to the urban texture will increase the accuracy of energy simulations and better reflect the impact of extreme weather conditions, as is discussed in our previous works [14,43]. Since the most extreme scenarios have been used in the energy simulations and analysis (extremes considering outdoor temperature, solar radiation and wind speed with monthly and hourly temporal resolutions for energy demand and generation calculations) considering 13 future climate scenarios (390 single years) with five different GCMs and three RCPs (which RCP 8.5 is considered as the worst case scenario), we can be very confident that the most pessimistic and extreme scenarios are considered in the resilience assessment.

Case study and urban morphologies
Representing urban complexity through a set of alternative morphologies requires selected morphologies to be matched with the local conditions while presenting energy-efficient alternatives. Section 4.1 presents a brief overview of the considered case study elaborating the local conditions while Section 4.2 describes the process used to generate alternative urban morphologies considered in the study.

Case study
Athens, the capital and largest city in Greece, is selected as the case study in this research work. It has an overall land area of roughly 414 km 2 and an urban population of approximately 3.2 million [44]. The area selected for the purpose of this study is located in Central Athens (hereafter referred to as Athens) with a total area of 87.4 km 2 , a population of about 1 million and a population density of 12,000 per 1 km 2 , containing eight municipal districts [45]. The city has a subtropical Mediterranean climate with hot arid summers and fairly mild winters [46,47]. High urban density, narrow streets and scarcity of green spaces contribute to the development of a strong urban heat island (UHI) index [48]. This results in a very large cooling demand during the summer months with considerably high peak values [49]. To recognize the urban morphology of Athens, the most frequent physical characteristics of urban neighborhoods (e.g. layouts, built density, building height, etc.) were studied for each municipal district. The majority of buildings in Athens are multi-story apartments with flat roof shapes, and high thermal mass (cement/concrete as the main material) [50][51][52][53]. Fig. 3 shows the studied area and site plan of four urban areas with different physical characteristics.

Generating urban morphologies
A method introduced by Javanroodi et al. [54], "Building Modular Cell" or "BMC", is adopted in this work to generate urban morphologies . BMC is based on the horizontal and vertical expansion of a basic module (4 × 4 × 4 m cube) to generate quadrilateral layouts in respect with several design-based constraints (e.g. connectivity between modules at least from one side, the generated configuration should be one integrated shape in each sub-site, the minimum number of floors for each building is two and maximum is twelve). A total number of fifteen influencing morphological parameters are used to generate the urban morphologies. The influencing parameters are categorized into two comprehensive groups of Urban Density (including plot to area ratio (PAR), volume to area ratio or (VAR), site coverage index (SCI), building density index (BDI), urban plan area density (λ p ), frontal area density (λ f ), relative compactness (Rc), building materials, and occupancy density (D oc ), and Urban Form (including building layout and building height distribution or 'BHD', neighborhood size, building function, streets and canopies geometry and materials) based on their morphological definitions. Table A1 in Appendix presents a description of each parameter and Fig. 4 illustrates the graphical definition of them.
A hypothetical site with a total area of 5120 m 2 is defined as a template to generate an urban neighborhood (20 × 16 grid). The defined site is fragmented into a matrix with three rows and columns resulting in nine equal sub-sites. Each sub-site is based on a 4 × 4 grid (covering 45% of the site) and a network of streets and canopies (covering 55% of the site). One sub-site (covering 5% of the site) is selected as open space, in which a total area of 64 m 2 (covering 0.01% of the whole site, and 25% of the sub-site) is considered as green space (there is scarcity of green spaces in the studied area). Two main streets (each 4 × 16 grid) and six canopies (each 2 × 4 grid) are defined between the nine equal sub-sites. To cover a wide range of urban areas in the selected case study, five ranges of urban densities are defined based on the adopted morphological parameters as follows: Each cell represents one thermal zone with a total area of 16 m 2 (a 4 × 4 m space). For urban form, four major architectural layouts are selected as the main design constraints in the developed version of the BMC technique, including Cube (C), Courtyard (CY) and L and U forms. As an attempt to define the vertical expansion of the developed layouts, five ranges of building heights (number of floors) are considered based on the defined UD ranges. However, due to the comparative approach of this study, at the central sub-site, a twelve-story building with a height of 48 m is placed in all twenty urban neighborhoods. Based on the location of each cell at the outer sides (building's elevation), a constant glazing ratio is defined: North elevation: 16%; West and East elevations: 8%; and South elevation: 25%.
According to this procedure, twenty distinct architectural layouts are selected out of 400 generated building layouts using the BMC technique, based on the most frequent building layouts observed in the studied area. Each UD category has a constant total area and number of floors with four different architectural form (C, CY, L, U). Moreover, two different building types including residential and commercial (office) are considered in all twenty cases as urban function parameters (with different materials, performance, heating/cooling systems, lighting/ equipment, etc.). The combination of two building types will lead to a more general and realistic representation of the energy system and the integration of renewable energy potentials based on complex urban neighborhoods. Fig. 5 illustrates the 3D model of the urban morphologies generated in this study with corresponding values for each parameter.
The models generated in Grasshopper are converted and exported into EnergyPlus models using Archsim, Diva for Rhino 4.0 plugin. Each cell is defined as a thermal zone with at least one connection with the adjacent cell at each architectural layout. The total number of thermal zones for cases in each UD category is 172, 144, 130, 109, and 88 respectively for UD 1 to UD 5 . The cooling/heating demand, equipment, and lighting energy demand profiles are calculated for all thermal zones on each floor of building blocks as the combination of latent and sensible loads. The impacts of the surrounding buildings on each building are calculated through the true view factor algorithm in EnergyPlus, which modifies the amount gained direct/diffuse solar radiation beams for external surface of each building. The total heat transfer through external surfaces (wall, roof, floors, and windows) as well as internal loads (e.g. people density, Equipment, Lightning, etc.) and infiltration through the building envelope and openings is defined based on the function of each building. To account for wind speed variations, the verified model in the previous works of authors [54,55,70] have been adopted based on BMC technique.

Co-optimization of urban morphologies and energy system
A typical formulation of energy system optimization maps the decision space variables that correspond to the energy system design to the objective space through a time series simulation. This can be formulated in a deterministic, stochastic, robust, or stochastic-robust way. Irrespective of the method used, design parameters of the energy system and dispatch variables are considered in the decision space. Both single and multi-level optimization algorithms have been introduced to optimize distributed energy systems. For example, Wang and Perera [56] introduced a multi-level optimization algorithm that can design distributed energy systems along with the connected network while maintaining n-1 security. Evins [57] introduced a bi-level optimization algorithm where the dispatch strategy of the distributed energy system is optimized at the primary level while system design is optimized at the secondary level. In general, multilevel optimization algorithms are computationally extensive. For example, the optimization algorithm introduced by Evins [57] takes several days to arrive at the optimal solution. Therefore, a single level optimization algorithm is used in this study. Single level optimization algorithms have been used previously to optimize the design of energy systems. Further, Wu et al. [58] used a single level optimization algorithm to determine optimal building energy system design along with the building renovation level. However, considering urban morphology adds a level of complication since urban morphology changes bring notably affect the demand profile, which cannot be simply presented using linear programming.

Optimization algorithm
Bringing energy system design and urban form into the optimization process entails formulating non-linear and non-convex objective functions, for which linear/mixed-integer linear programming techniques or a gradient-based technique are unsuitable. The direct influence of urban form on the energy demand makes it difficult to decouple the optimization process into two levels where heuristic methods with linear programming techniques can be used as introduced by Wang and Perera [56]. Hence, heuristic methods are used in this study as the optimization technique. The decision space vector presents the variables related to the energy system configuration, dispatch strategy and energy demand. System configuration defines the size of system components such as PV panels, wind turbines, battery bank, etc. Dispatch strategy defines the fuzzy and finite-state transition rules. Energy demand is determined by two decision space variables corresponding to urban form and density. The entire set of variables amounts to 30, which is a considerable number to be handled by a heuristic algorithm. Net Present Value of the system and grid integration levels are considered as objective functions while power supply reliability is used as a constraint in the optimization. The constraint tournament method [59] is used to consider the constraint function during the optimization. A co-operative co-evolutionary algorithm (COCE) often used to consider very extensive decision spaces is used to conduct the optimization. The co-evolutionarily algorithm enables partitioning the decision space while cooperativeness enables considering non-variable separable problems along with the partitioning of the decision space variables. An extended explanation is presented in [60].

Outline of the energy system
An energy hub model introduced by Guidl and Anderson [61], which has been widely used during recent studies, is adopted in this work to represent the energy system. Solar PV panels and wind turbines are considered as non-dispatchable energy technologies. An internal combustion generator (ICG) is used as the dispatchable energy source. A battery bank is considered as energy storage. The heating and cooling demand of the buildings is assumed to be catered by the heat pumps and air conditioners. The energy hub interacts with the grid when catering the fluctuating energy demand of the buildings. Grid curtailments are considered when injecting the excess generated as well as buying electricity from the grid. The hourly pricing scheme is considered for both selling and purchasing electricity to and from the grid.

Dispatch strategy
Dispatch strategy stands for the energy management strategy implemented in the energy system to withstand the variations in energy demand, renewable energy generation and grid conditions. The dispatch strategy determines the optimal use of energy storage and can be presented in different ways when conducting energy system optimization. A number of different methods have been used to implement a dispatch strategy in energy system optimization. Bi-level optimization has often been used to link dispatch and system optimization problems. Rulebased models, grey-box models and reinforcement-learning algorithms have been used to implement the dispatch strategy in the energy system optimization problem. Although reinforcement learning can perform exceptionally well concerning complex energy systems, it demands much higher computational time and resources [60]. Therefore, the grey box model is considered in this study. A bi-level dispatch strategy is used in the grey box model [62]. The primary level is based on fuzzy automata, where operating conditions of the ICG are determined based on the mismatch between energy demand and renewable energy generation, state of charge in the battery bank and price of electricity in the grid. The secondary level uses finite state automata to derive the energy interactions between grid and energy storage. The energy interactions are determined considering the price of electricity in the grid, the state of charge of the battery bank, and grid curtailments. A more comprehensive definition of the state space of the dispatch strategy is presented in Ref. [63].

Co-optimization of the energy system and urban form
A detailed energy flow model is used in this study, linking the urban simulation model with the energy system model. Energy demand on an hourly scale is fed into the energy system model, which depends on the urban morphology. Hourly power generation using Solar PV panels and wind turbines is computed subsequently. Power generation from the PV panels and wind turbines depends on hourly global solar irradiation on the tilted solar PV panel surface and wind speed at the wind turbine hub level as well as installed capacities of these technologies. For example, renewable power generation using solar PV panels (P PV t ) for time step t is computed using Eq. (1). ] and x pv denotes the area of the solar panel and number of PV panels installed, which is optimized in the optimization algorithm. Complexities inherent to city morphologies make it difficult to consider the impact of shading, especially when performing energy system sizing. We introduce a shading factor ς PV to account for the losses due to shading. A similar approach is used to compute the power generation using wind turbines Eq. (2).
In Eq. (2), P W t (v t ) [kWh] denotes the power generation from a single wind turbine. x w represents the number of wind turbines which is optimized using the computational model. Wind speed will be dropped due to the urban boundary layer and losses in power converters are represented byς w in the system (which is obtained using the optimization algorithm) and the power losses. Power generation from the dispatchable sources, energy interactions with the grid, and energy storage are determined by the dispatch strategy. The state of the charge model is used to determine the charge level of the battery bank while the polynomial power curve is used to determine the fuel consumption of the dispatchable energy source.

Formulation of objective functions and constraints
Interactions with the grid are maintained within bounds. Grid curtailments are introduced for injecting and purchasing electricity to and from the grid. Grid curtailments act as the primary barrier, which guarantees smooth operation of the grid while accommodating renewable energy technologies at a large scale. As the study focuses on autonomous energy districts, dependence on the grid is be minimized. Less dependence on the grid encourages integrating more renewable energy technologies within the district while helping to maintain the stability of the grid. In most of the instances, the grid acts as virtual storage when integrating renewable energy technologies, taking the excess generation within the micro-grid while supporting the catering of the energy demand whenever renewable energy generation is not sufficient. Support from energy storage and dispatchable sources is required to balance the difference between demand and generation, which will add an additional cost to the system, creating a Pareto front between cost and the autonomy level of the district. Therefore, cost and autonomy level are considered as the objective functions. Grid integration level is used as a representative of the autonomy level. Lower grid integration levels will lead to a higher autonomy level. Grid integration level is defined in this study according to Eq. (3) [62].
In this equation, E IG t [kWh] and E D t [kWh] denote energy imported from the grid and energy demand of the hub. The energy demand curve is strongly influenced by the urban morphology. Hence, certain morphologies would lead to maintaining minimum interactions with the grid. This makes it important to optimize the archetype along with the system components.
The net present value of the system represents the financial feasibility. Initial capital cost (ICC) of the system components as well as the operation and maintenance costs are considered under the net present value. Initial capital cost represents the acquisition and installation cost of system components such as PV panels, wind turbines, battery bank, etc. The expenditure of constructing the building stock is not considered in this study. Operation and maintenance costs consist of variables as well as fixed operation and maintenance costs.
In Eq. (4), CRF c denotes the Capital Recovery Factor for the c th component. PRI denotes the real interest rate calculated using both interest rates for investment and the local market annual inflation ratio.
The year considered is represented by h.
In addition to the two objective functions considered for the Pareto optimization, Loss of Load Probability (LOLP) is considered as a constraint in the optimization process. LOLP has often been considered to represent the power supply reliability [65][66][67][68].
The loss of power supply (LPS) of the system is computed according to Eq. (5).

Results and discussion
Energy requirements in the urban context can be influenced by both urban morphology and climate conditions, which can notably influence the design as well as the operation of distributed energy systems. A comprehensive quantification of the influences brought by urban morphology and climate variations can lead to improving energy efficiency as well as sustainability from the generation perspective. Towards this objective, Section 6.1 presents the influence of urban climate and urban morphology on the energy demand while Section 6.2 quantifies the impact of both these sectors while simultaneously optimizing the energy system and urban morphology.

Impacts of urban morphology on energy demand
Results in this section help to understand how the energy demand profiles (summation of heating, cooling, lighting, and appliance demands) alter based on the considered urban morphology. Distributions of the hourly values of energy demand over one year (for typical and extreme weather conditions) are shown in Fig. 6  As is visible in Fig. 6, a downward trend can be noticed, in which a higher urban density results in higher energy demand. The difference in the average annual energy demand for typical and extreme weather conditions varies for each building form and urban density (UD 1 to UD 5 ). For example, for the C-form buildings, the highest energy demand on an annual scale is observed in ECY; where the average absolute difference between ECY and EWY is +20, +19, +17, +14, and +13 percent for UD 1 to UD 5 , respectively. The lowest difference between average annual energy demand for ECY and EWY is observed in the L-form buildings in UD 5 with +9 percent; which clearly indicates the impacts of urban morphology on the energy performance of buildings with similar building forms in extreme weather conditions. An important part of the assessment is to quantify the impacts of each UD (with similar built density, conditioned volume, and the number of floors) on the annual distribution of energy demand. Fig. 7 compares the cumulative distribution of the annual energy demand for all the building forms over three time periods with typical and extreme weather data sets. The CY-form buildings have the best performance in ECY for all the urban morphologies, except for UD 2 and UD 3 where U-and C-form buildings have the lowest energy demand on an annual basis. Similar distributions can be observed in EWY for CY-form buildings in UD 1 , UD 2 , and UD 4 ; while in UD 3 and UD 5 the U-and L-form buildings have a notably lower energy demand. For example, in UD 1 (A a = 13824 m 2 ), CY-form have the best performance with average energy demand of 211. 8, 197.4 and 175.1 kWh in ECY, TDY and EWY; while L-form buildings showed about 9%, 8% and 5% higher energy demand respectively. Table 1 shows the average energy demand for all urban morphologies based on building form and UDs. There is an upward trend during To assess the impacts of urban density and building form on the energy demand profiles, the annual energy demand values per square meter are compared in Fig. 8. The lowest annual energy demand for ECY and EWY conditions is observed in UD 5 for all building forms. For example, in C-form buildings, the annual energy demand for UD 1 to UD 4 is 10, 19, 10, and 11 percent higher than UD 5 with 100.4 kWh/m 2 energy demand in EWY. A similar condition is observed for ECY; where annual energy demand difference for UD 1 to UD 4 is +18, +27, +14, and +6 percent compared to UD 5 with 117.3 kWh energy demand. For TDY, these numbers are +14, +23, + 11, +0.1 percent indicating a similar trend. The highest annual energy demand in ECY for CY-, L-, and U-form buildings is observed in UD 3 with 149.9, 146.9, 145.8 kWh/m 2 respectively. The highest annual energy demand overall observed in UD 4 (138.8 kWh/m 2 ) using EWY weather dataset. Similarly, the energy demand of all urban morphologies shows the highest increment by time with up to 17 kWh/m 2 . It is interesting to mention that this increment is even higher in L-form urban morphologies. This is why C-and CY-form buildings with more compact geometries and semi-courtyard spaces show a better performance; particularly during warm seasons.
Another important part of the assessment is comparing the impact of urban morphology on the peak energy demands during typical and extreme weather conditions (Fig. 9). The highest peak demand in all UDs and building forms occurs in EWY, except for the U-form buildings in UD 1 whose highest peak demand is observed in ECY. The magnitude of peak demand increases notably by time in all UDs for EWY. It reaches over 1073 kWh during 2070-2099 for the U-form buildings in UD 1 , which is 25% more than the value for 2010-2039. Urban morphologies with C-and CY-form buildings showed a lower magnitude of peak energy demand in all UDs. According to the results, energy demand in EWY  conditions will increase notably in all urban morphologies due to the impacts of climate change on the variations of air temperature. This will considerably affect the energy system design for the future conditions and mitigate the larger and more frequent average and peak demands for each UD.

Effect of urban morphology on the energy system
The impact of urban morphology can easily go beyond the energy demand and influence energy system operation, which might be vital in case of extreme climate events. To quantify the impact of urban morphology on the energy system, the energy system is optimized considering the sensitivity of building form, urban density, and both these aspects together, as presented in Sections 6.2.1 and 6.2.2 respectively. Subsequently, the impact of extreme climate conditions on energy system design is investigated in Section 6.2.3.

Sensitivity of building form
To understand the impact of building form on energy systems, a Pareto optimization is conducted for each building form considering the density of UD 3 , taking demand profiles for the 2010-2039 period. Although the Pareto fronts do not show a significant difference in NPV (Net Present Value), a noticeable split is observed in two groups (Fig. 10). L and CY can be grouped into one class, which shows a relatively lower cost compared to the groups of U and C. The difference in cost is not distinguishable when the grid interactions are low. The Pareto fronts almost coincide when reaching the standalone operation mode. The difference in cost is kept within 10% throughout the Pareto front, which indicates that it closely follows the results of the building simulation where annual energy demands are kept within a 10% bound. This indicates that the changes observed in annual demand does not reflected in the Pareto fronts.

Sensitivity of urban density
The influence of urban density is clearly visible in both annual and peak energy demands as discussed in Section 6.1. To understand the influence of urban morphology on energy systems in a more holistic manner, the sensitivity of energy systems to urban density must be considered. Towards this objective, the energy system is optimized along with the urban form while maintaining the building form (building form C). As is visible in Fig. 11 (a), a notable difference in cost is observed when changing the urban density. For example, cost per unit  area can increase by more than 50% when moving from UD 5 to UD 1 . Such a significant change in the objective function values easily results in a significant change in the required energy system design. When analyzing the five Pareto fronts, a higher density will lead to a higher cost per unit area except for UD 3 . A significant shift in the objective function values can be observed when moving from one density class to another. UD 1-3 shows a gradual reduction in cost with the increase of grid integration level while UD 5 presents a notable reduction in cost when the grid integration levels are low in contrast to other urban densities. However, the drop-in cost gradually diminishes with the increase of grid integration level. By contrast, UD 1-3 maintains the drop in cost even after the grid integration level reaches beyond the 1250 kWh/ m 2 limit as shown in Region B. More importantly, UD 1-5 reach lower NPV values when compared to UD 5 in Region B especially towards the end of the Pareto fronts. This analysis makes it easier to understand the complete Pareto front considering both system configuration and urban morphology as presented in Fig. 11 (b). It is interesting to investigate the reasons for the cost per area increasing with urban density. The higher cost may occur due to1) higher energy demand in buildings, and 2) higher cost in catering to the demand. Table 1 provides a clear explanation for the first hypothesis. It is clear that certain urban forms are more energy-efficient than the others, both for annual and peak energy demands. This clearly explains why UD 5 leads to a lower cost per unit area when compared to the others. However, the analysis of the energy demand does not explain the behavior of the Pareto fronts in Region B, which leads us to the second point (performance of the energy system and higher cost in demand catering). In Fig. 12 the Levelized Energy Cost (LEC) is plotted against the levelized grid integration, which provides an overview of the cost of generating energy units for each Pareto solution for the five Pareto fronts. When analyzing the LEC of the Pareto solutions it is clear that all the design solutions show similar LEC when grid interactions are at a minimum level. However, the gap between the Pareto fronts becomes noticeable when increasing the grid interaction levels. The morphologies with higher urban density become more economical when compared to the morphologies with a lower urban density. For example, the LEC decreases from 0.17 to 0.12 Euros when moving from UD 5 to UD 1 within Set A, reducing the cost by 45%. UD 1 caters a much larger energy demand when compared to UD 5 , which leads to generating energy at a much cheaper cost due to the scale of the economy. This explains that LEC decreases with an increasing grid integration level. However, when reaching standalone conditions the energy system that caters denser morphologies needs to work harder, since morphologies with higher urban density need to cater higher peak demands when compared to the morphologies with a lower urban density. There are significant drops in part-load efficiencies at lower operating lead factors, diminishing the advantage obtained from the scale of economy. As a result, LEC values coincide with all the Pareto solutions when grid integration levels are low. As the LEC is similar for all the densities and UD 1 and UD 2 have a much higher demand, the cost of UD 1 and UD 2 is high when compared to UD 5 in Region A. However, the levelized generation cost notably decreases for both UD 1 and UD 2 compared to UD 5 , which compensates the higher energy demand. As a result, the NPV per unit area is low for UD 5 when grid interactions are lower (Region A), while UD 2 presents a lower NPV per unit area when the grid interactions are high. This leads us to understand the Pareto front obtained, considering the full decision space including energy system configuration, urban density, and building form as presented in Fig. 11 (b).

Sensitivity of energy systems to climate variations
As discussed in Section 6.1, climate change notably influences the energy demand of buildings. Typical demand profiles consider the gradual changes in the demand profile but not extremes while ECY and EWY do consider extremes. An increased average temperature will result in a higher cooling demand on average while stronger extreme events will increase the peak energy demand. In this section, the sensitivity of the energy system to climate variations (considering both long-and short-term variations) is discussed. Pareto fronts obtained in previous sections are derived considering the typical energy demand as explained in Section 3, which is the usual practice. Rising temperatures put building HVAC systems under strain by increasing the energy demand especially during extreme climate events such as heatwaves. Considering such extreme climate conditions will be vital in the future, especially with a larger penetration level of renewable energy sources and limited availability of backup dispatchable generators based on fossil fuels. Hence, Pareto optimization is performed for urban density UD 5 under typical and extreme climate conditions. A significant increase in NPV by up to 40% (Fig. 13) is observed when moving from a typical scenario to extreme conditions. The increase in cost that is observed in energy infrastructure goes well beyond the annual energy demand increase for extreme climate events. This clearly reflects the magnifying effect of climate impact with regard to extreme climate events when shifting focus from the building sector to the energy sector.

Conclusions
The present study extends the scope of the present state of the art by considering the influences of urban morphology on energy generation, and simultaneously optimize urban morphology and energy systems. As a first step, three groups of future weather datasets were synthesized based on thirteen climate scenarios for Athens, Greece, for three time periods (2010-2039, 2040-2069, 2070-2099). Each group contained three weather datasets to represent typical (Typical Downscaled Year or "TDY"), and extreme (Extreme Cold Year or "ECY", and Extreme Warm Year or "EWY") weather conditions (nine weather datasets in total). Then, a total number of twenty urban morphologies were generated based on physical characteristics of urban areas in eight municipal districts in Central Athens using a technique called Building Modular Cells or "BMC". The energy performance of each urban morphology in typical and extreme conditions was thoroughly discussed based on major morphological parameters in three time periods. Subsequently, the urban morphology was optimized along with the energy system configuration to optimize energy system design, building form and urban density. A clear link between peak energy demands with building form and urban density was observed, where the C-and CY-form urban areas with higher than 60% urban density had the lowest peak demand. Results also showed that the average annual energy demand increases with time (from period 2010-2039 to 2070-2099) by up to 16 percent in extreme warm conditions and by up to 20 percent in the case of peak energy demand. This is due to the air temperature increment, both on average and hourly (extreme values) basis, particularly in extreme warm conditions such as occur in Athens. Although the influences of urban density and building form on the energy demand profiles were visible, the latter showed only trivial impacts on the energy system.
Urban density presents a clear influence on the energy system. The impact of urban density on the NPV per unit area is different from the levelized energy cost, which clearly reflects that the problem can be looked at from two different perspectives. Nonetheless, both NPV per unit area as well as levelized energy cost can be reduced by above 45% by appropriately selecting the most suitable urban form where the cooptimization algorithm used in this study can be helpful. The study reveals that the energy system can withstand changes in energy demand up to a threshold without demonstrating much change in performance indicators. This leads to trivial changes in Pareto fronts concerning the building morphology. However, a notable change in the energy system can be observed when moving beyond this threshold, considerably exceeding the changes in energy demand observed. This clearly reflects the importance of considering buildings and distributed energy generation together as interconnected. Urban morphology plays a significant role when linking these two sectors since there is a notable change in energy demand when comparing alternative designs. The same principles apply when studying the impacts of climate variations and extreme events. The considered time span does not affect the energy system substantially when only typical climate conditions are taken into account, using the deterministic framework considered for the study. However, extreme climate events lead to significant variations in the energy system, as there is a notable increase in energy demand during extreme events. With the increasingly frequent occurrence of extreme climate events, it is essential to quantify the combined effect of buildings and energy infrastructure, which enable improving the climate resilience of cities. The influence of urban climate or urban microclimate on the energy demand was not considered in this work, which might lead to significant changes. This is going to be investigated as a future extension of the developed framework. The proposed framework will be immensely helpful in such contexts to design climate-resilient cities considering the close interactions between buildings and energy systems.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgement
The Swedish contribution was supported by the framework of the joint programming initiative 'ERA-Net Smart Energy Systems' focus initiative on Integrated, Regional Energy Systems, with support from the European Union's Horizon 2020 research and innovation programme [775970].

Appendix A
See Table A1 and Figs. A1 and A2.

Table A1
The definition of fifteen morphological parameters used in this study to generate urban morphologies based on Urban density and Urban Form.

Category Parameters Descriptions
Urban Density (UD 1 , UD 2 , UD 3 , UD 4 , UD 5 ) VAR (volume to area ratio) Total volume of all buildings divided by the total site area of site PAR (plot to area ratio) Total floor area of all buildings divided by the total site area of site SCI (site coverage index) The area of the ground floor of the building divided by the total area of site BDI (building density index) The area of the ground floor of the building divided by the area of each sub-site λ p (urban plan area density) The built area projected onto the ground surface divided by horizontal section A h λ f (frontal area density) The surface area of three central buildings divided by their plot area of a horizontal section A h Rc (relative compactness) The volume of each building divided by the total area of external surface, R Ca the average R C of buildings in an urban morphology Building materials The major inputs used for the material for all urban morphologies D oc (