Optimisation of curvilinear external shading of windows in cellular offices

Shading of windows influences building cooling and heating loads through control of solar heat gains, and lighting load through access to available daylight. Shading shape thus presents an important factor both in building energy analysis and building aesthetics. Curvilinearity of solar paths suggests that the optimal shading shape may be curvilinear as well, and our aim here is to test this expectation. To accommodate curvilinearity of shading shape, outer edges of shading, which consists of overhang, western and eastern fins, are modeled as non-uniform rational basis spline (NURBS) curves, a widely accepted representation standard for curves in design industry. As a case study, a cellular office is considered in the Pacific Northwest National Laboratory (PNNL) office building model, with its overhang lined up by seven control points, and the fins lined up by five control points each, with two ending control points joint for the overhang and the fins. With control points allowed to take on nine different alternative depths, genetic optimisation is employed for 16 representative USA climates with respect to total equivalent source energy for heating, cooling and lighting loads. The main finding is that in a very close proximity to optimal shadings found by genetic optimisation there exist shadings with much simpler control point structure, obtained by identifying depths of successive control points, that have nearly rectangular overhangs. Since the difference between these simpler shadings and the optimal ones is less than 0.24%, this partially rejects the expectation that the optimal shading shape should be curvilinear. Structure of these near-optimal shadings also suggests a new way to partition shadings into independent regions: the lower and the upper parts of the western fin, joints of the overhang with the western and the eastern fin, the interior part of the overhang and the rest of the eastern fin.


Introduction
A large part of energy in buildings is used to provide comfortable thermal conditions to their occupants. Shading of windows is a well accepted way of passive reduction of cooling loads, popularised by Le Corbusier in the form of brise-soleils in his Unité d'habitation buildings in the 1950s [1]. Le Corbusier's goal was that brise-soleils fully shade windows at noon in summer months, but not to block the sun during winter months. This approach led to development of a number of methods for shading design that are based on solar path projections and cut-off a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 days and hours, during which complete shading of windows is required. Early examples include methods of Olgyay and Olgyay [2] and Mazria [3], where Olgyay and Olgyay use a projection of the sun onto a horizontal plane, while Mazria uses a projection onto the vertical cylinder. In both cases the overheating period is plotted onto the solar path diagram which helps to define the shading mask that blocks direct sun radiation during this period. Shaviv [4,5] developed a computer program that from a solar path generates a set of possible shadings which prevent direct radiation for each month. Arumí-Noé's algorithm [6] first finds a winter solar funnel surface that ensures full insolation, and then clips it subject to the summer shading conditions. Dubois [7] developed a chart, complementary to Mazria's solar path diagram, by taking into account additional information on the window solar heat gain coefficient and the incidence angle between the window and the sun beam. Marsh [8] proposed a method that obtains shading shape by projecting solar position at required cut-off dates and hours onto a plane of the shading. Cheung and Chung [9] divided the sky map into small 5˚× 5˚patches in order to be able to predict probable sunlight duration on windows in a densely packed building environment, which is then used in the design of exterior shading.
Another approach present in the literature consists of dividing shading support surface into smaller cells and processing each cell separately. In his cellular method Kaftan [10] divides the shading surface into two-dimensional array of cells and for each cell calculates the amount of direct solar gains it prevents during time periods for which shading is required. In their Shaderade method, which builds upon the Kaftan's method, Sargent et al. [11] calculate for each cell and each time period the desired fraction of solar beam energy transmitted through the cell at that time period and then find the optimal transmittance value for a cell from an annual calculation. A welcome characteristics of these two methods is that the cells may be sorted according to calculated values which allows one to identify most and least effective areas of the shading support surface and to collect cells in decreasing order of effectiveness in order to obtain shading with required surface area. However, the shapes suggested by these methods tend to be serrated, so that architect's interference is necessary to produce aesthetical design.
Nowadays, more than 50 years since Le Corbusier's seminal buildings, heating loads have become significantly lowered due to standard use of high insulation and tight sealing in construction. In office buildings in milder climates they are even similar in magnitude to lighting loads, due to both high internal gains and minimal illuminance level requirements. Hence it becomes necessary for office shading design methods to take into account all of the heating, cooling and lighting loads. Kaftan and Marsh [12] combined Kaftan's cellular method with Ecotect in order to predict necessity of shading based on both direct solar gains and thermal comfort indicators. Calculations of desired cell transmittances in Shaderade [11] already include EnergyPlus predictions of heating and cooling loads, but for the base case of a building without shading only. These methods, however, cannot account for lighting load, due to imponderable influence of each particular cell of the shading support structure.
On the other hand, recent shading design methods usually employ genetic algorithms to overcome the difficulty of handling lighting load together with heating and cooling loads. Castorina [13] encodes a shading covering the whole façade with a particle-spring system and employs genetic algorithms and EnergyPlus to optimise a single objective function representing weighted combination of illuminance ratio, lighting load, and the ratio of winter and summer solar gains. Ercan and Elias-Ozkan [14] also consider the whole façade shading by encoding the depth and angle of each shading device and use genetic algorithms and Radiance to optimise daylighting levels, but they do not consider interior thermal conditions. Manzan [15] uses DAYSIM to estimate lighting load and ESP-r to calculate heating and cooling loads, and then use genetic algorithm to determine optimal angle and reveal of a full width rectangular overhang.
Here we are interested in optimal shape of external shading for a south oriented window in a cellular office. As solar path is curvilinear, it is natural to expect that the optimal shading shape will also be curvilinear. To accommodate this curvilinearity, it is assumed that the shading, consisting of an overhang, eastern and western fin, tightly placed around the window, has its outer edges modeled with NURBS curves. Due to their rather general definition, NURBS curves can be used to approximate arbitrary curves, as can be evidenced from a number of earlier works [16][17][18][19][20][21]. Using NURBS curves here enables sampling the huge space of smooth shading designs by a controllable number of NURBS-based designs, as the shape of NURBS curves is determined by the number of control points and their feasible positions. A genetic algorithm is employed to optimize positions of these control points for representative USA climates with respect to a weighted sum of heating, cooling and lighting loads as a single objective. Discussion is then oriented toward understanding optimisation results that in some cases seem counterintuitive, and toward producing simplified designs with loads close to the optimal ones. Results reported here are continuation of studies initiated in [22,23].

Cellular office model
The prototype large office building models [24], developed by PNNL and derived from the Department of Energy (DOE) Commercial Reference Building Models, are used as a starting point for this study. These models, whose development and definitions are closely described in [25,26], represent realistic building characteristics and construction practices, aim to cover more than 70% of the commercial building floor area in the United States for new construction and conform to the American Society of Heating, Refrigerating and Air-Conditioning Engineers (ASHRAE) 90.1-2013 standard. The models are prepared for simulations with Energy-Plus, with version 8.4 used here. They are provided for all main USA climate zones according to the climate zone classification system developed by Briggs et al. [27]. Locations of the PNNL building models, together with a few selected parameters, are listed in Table 1. In addition to USA locations, a model is provided for Vancouver, Canada, as well, which represents the cool, marine climate zone 5C.
The present study aims to study effects of shading a single office window, so that a cellular office has been set up as a single zone, using settings, materials, constructions and schedules of the underlying PNNL models. The office has width 3.60m, height 2.80m and depth 5.16m, due to the 200 ft 2 /person requirement from [25]. The external wall is oriented toward south, while other walls are assumed to be adiabatic. In order to keep the same windows-to-wall ratio as in the PNNL office building model, the window of the cellular office has been set to be of width 2.99m and height 1.99m, with the window sill at 0.6m from the floor. The office model is illustrated in Fig 1(a). Window glazing properties depend on climate zone, and their U-values and solar heat gain coefficients are listed in Table 1. Lighting power density is set at 1 W/ft 2 . The office has daylighting control with two sensors placed at desk level at one third and two thirds of the office depth, and illuminance setpoint of 375 lux. Since the simulations are run for a single office instead of a whole building, detailed heating, ventilation and air-conditioning (HVAC) from the PNNL models had been replaced with IdealLoadsAirSystem.
Two EnergyPlus simulation parameters are important for proper calculation of shading effects. The calculation method field of the ShadowCalculation object is set to TimeStepFrequency to perform solar path, shadowing and diffuse sky modeling calculations at each timestep (set at 15 minutes). The solar distribution field of the Building object is then set to FullInteriorAndExterior to compute shadow patterns on exterior surfaces by the window shading and to calculate amounts of transmitted beam radiation falling on each internal surface by projecting the sun's rays through the window.
Archive with simulation files for a cellular office model for all locations is available at [28].

NURBS curves and epnurbs package
NURBS is a widely accepted standard in computer-aided design, engineering and manufacturing for describing and generating smooth curves and surfaces [29,30]. A NURBS curve is defined by a sequence of control points P i , i 2 I for some index set I, that act as if P i were connected to the curve by a spring of strength w i . Each point of the NURBS curve C(t), 0 t 1,   is actually a convex combination of control points: where N i (t) are suitably calculated basis functions. The basis functions are determined by a degree d and a knot vector which partitions the interval [0, 1] into knot spans, in such a way that ∑ i2I N i (t) = 1 holds for each t 2 [0, 1] and that each basis function has d + 1 consecutive knot spans on each of which it reduces to a polynomial of degree d while it is equal to zero outside these knot spans. These conditions ensure that each curve point is determined by d + 1 closest control points. Details of computation of basis functions may be found in [29,30]. The window shading for the cellular office model consists of three parts: western fin, overhang and eastern fin. They are placed tightly around the window in vertical planes (fins) or in a horizontal plane (overhang), with their outer edges modeled as NURBS curves. There is a total of 15 control points: whose coordinates are given with respect to the lower left corner of the outside surface of the exterior wall. The NURBS curve of the western fin is determined by control points P 0 , . . ., P 4 , the one of the overhang by control points P 4 , . . ., P 10 , and that of the eastern fin by control points P 10 , . . ., P 14 . All control points are of unit weight. For each curve the basis functions are of degree three and are determined by a clamped uniform knot vector, that starts and ends with three empty knot spans, which ensures that the curve starts with its first control point and ends with its last control point. Due to the joint control points P 4 and P 10 , ends of the overhang coincide with the upper ends of fins, giving the shading shape a continuous look. It should be noted, however, that internal control points of curves do not necessarily belong to them, so that the y i values do not represent actual shading depths, but only the (negative) distances of control points from the external wall. An example of window shading for the cellular office model, with control points shown as small boxes, is illustrated in Fig 1(b). In addition to smoothness, another benefit of using NURBS to model outer edge of shading is possibility to control the size of the search space. Here it is prescribed that y i takes its value from the set of alternatives {0, 0.25, 0.5, 0.75, 1, 1.25, 1.5, 1.75, 2} for each i = 0, . . ., 14, so that the search space contains a total of 9 15 feasible curves. Although large, this set of NURBS curves with uniformly distributed control points, with clamped uniform vector and unit weights certainly cannot serve to approximate all imaginable smooth curves that could be used as outer edge of shading of depth at most 2m. Nevertheless it still provides a representative finite sample of the latter infinite set. The fact that NURBS curves are determined by a relatively small number of parameters, in this case by the y i values for i = 0, . . ., 14, further makes it possible to search for the optimal shading shape by using genetic optimisation. The number of control points (15) and their feasible positions (9 for each) were chosen here so that the resulting size of the search space (9 15 candidate curves) is, in our opinion, neither too small to give a poor sample of all possible smooth shading designs nor too big for genetic optimisation.
Before that, however, another obstacle had to be overcome as EnergyPlus can model building geometry using rectilinear surfaces (with at most four vertices) only and cannot handle NURBS curves directly. A general solution for this problem is to uniformly divide the domain interval [0, 1] by points t i = i/k for i = 0, . . ., k and a selected positive integer k, and then to approximate NURBS lined shading with a number of adjacent trapezoidal shadings, where the i-th shading, for i = 0, . . ., k − 1, has as vertices the curve points C(t i ), C(t i+1 ) and their projections on the wall surface. In preliminary studies [22,23] calculation of NURBS curve points was hard coded in the cellular office model using the EnergyPlus macro language EPMacro. Due to the absence of programming loops in EPMacro, calculation code had to be rewritten for each point, which made this a rather cumbersome solution that is not easy to update. Coding was further complicated by the fact that EnergyPlus does not allow surfaces that have a side shorter than 0.01m. This is hard to control in advance as the domain interval [0, 1] is not mapped uniformly to a NURBS curve: parts with higher curvature require a larger portion of the domain, while a smaller portion of the domain is sufficient to represent parts with lower curvature.
In order to have a solution that is easily applicable to different models, we developed a Python3 package epnurbs that creates a NURBS lined shading for a general wall surface in EnergyPlus. Its source code is openly available at [31], while the simplest way to install it is by issuing terminal command pip install epnurbs, which also installs necessary dependencies: eppy [32,33] to handle reading, updating and writing EnergyPlus files and NURBS-Python [34,35] to calculate NURBS curve points.
The main method of epnurbs package is createnurbsshading with signature createnurbsshading(idd_filename, idf_filename, base_surface, shading_str, ctrl_points, evaluated_points = 20). Arguments have the following meaning: • idd_filename is a path to the EnergyPlus idd file, which contains definitions of all objects in a particular version of EnergyPlus; • idf_filename is a path to the EnergyPlus idf or imf file, which contains definition of the building model; • base_surface is the name of the surface in the building model to which NURBS lined shading should be attached; • shading_str is the string containing definition of shading objects that will be created as an approximation of NURBS lined shading. This string may contain the following placeholders that will be replaced with actual values: • <IDX> is replaced by the ordinal number of the shading object used in approximation; • <BASESURFACE> is replaced by the name of the base surface; • <VERTICES> is replaced by a list of vertex coordinates in counterclockwise order; • <COUNTERVERTICES> is replaced by a list of vertex coordinates in clockwise order.
A sample shading string used in this study is 'Shading:Zone:Detailed, Shading <IDX>, <BASESURFACE>, , , <VERTICES>;'; • ctrl_points is a list of coordinates of control points for the NURBS curve defining outer edge, where each control point is given as a triplet [X, Y, Z] of its coordinates calculated with respect to the zone that contains the base surface; • evaluated_points is the number of points to be calculated on the NURBS curve, which is also equal to the number of trapezoids that will be created to approximate the NURBS shading. If not given, its default value is 20.
Method createnurbsshading first loads the EnergyPlus idf file and finds the base surface, and then calculates NURBS curve points and their feet of perpendiculars to the base surface. For each pair of consecutive curve points it then creates a trapezoidal shading object with the pair of curve points and their feet of perpendiculars as vertices and adds it to the Ener-gyPlus idf file, provided that all sides of such trapezoid are at least 0.01m.
Package epnurbs currently contains two more methods: createnurbsopening that approximates a NURBS lined window with a set of rectangles, and createrectshading that creates a sequence of rectangular shadings with given depths. These methods are not used in the present study, so further details about them may be found in [31].

Simulation management with jEPlus and jEPlus+EA
EnergyPlus simulations for locations mentioned in Table 1 were managed with jEPlus [36][37][38]. jEPlus enables one to perform parametric EnergyPlus simulations by describing a search space with sets of alternative values for specified simulation parameters and running simulations either for the whole search space or its representative sample. Parameters in this study are the positions of the control points P 0 , . . ., P 14 , encoded with integers from {0, 1, . . ., 8} which are multiplied with -0.25m to produce the respective y i coordinate. jEPlus also provides the ability to call Python methods to preprocess simulation files: it supplies the method with the names of three folders that contain project files, simulation results and EnergyPlus idd file, and other arguments specified in the parameter definition, passed in as a comma-delimitted string. Python preprocessing method then creates a list of alternative positions of control points, selects the appropriate positions based on current parameter values and calls createnurbsshading from epnurbs to add to the cellular office model the western fin approximated with 10 trapezoids, the overhang approximated with 15 trapezoids and the eastern fin approximated with 10 trapezoids, after which the model is simulated with EnergyPlus v8.4.
However, jEPlus cannot be used directly in the search for optimal shading shape due to the prohibitively large search space and instead an optimisation method has to be applied. Coupling of building energy simulation tools with optimisation methods has become mainstream in the study of energy and buildings after Caldas and Norford [39] used it prominently to facilitate performance-based façade design. A number of reviews on this topic are available: Machairas et al. [40] review methods and tools used for the building design optimisation while, more specifically, Kheiri [41] reviews optimisation methods for building geometry and envelope design and Stevanović [42] reviews work on optimisation of passive solar design of buildings. Although different optimisation methods have been used in building design optimisation, such as direct search, simulated annealing, particle swarm optimisation, harmony search or ant colony optimisation, they appear rather sporadically in literature, while large majority of building design optimisation studies rely on genetic algorithms, which tend to work well for problems considered in this area an are used in this study as well. In this aspect, one might object that building design optimisation community is somewhat lagging behind other engineering communities, in which newer optimization methods, such as gravitational search [43,44] or Jaya [45], are more easily embraced and applied in research (see, e.g., [46][47][48][49]).
Genetic algorithms, inspired by biological evolution, work by evolving a population of candidate solutions for the optimisation problem over a number of generations by repeated application of selection, reproduction, mutation and recombination, with the goal of improving candidates' fitness, which is given by the problem's objective function. The most well known genetic algorithm variant is the non-dominated sorting genetic algorithm NSGA-II [50], which is implemented in the jEPlus-based optimisation tool jEPlus+EA [51] used for this study. Note that Python preprocessing in jEPlus+EA is not available in versions prior to v1.7.7 beta. For each model location, population size was set to 50, somewhat larger than the number of bits needed to describe parameter values. Crossover rate was set at 90% to quicker explore part of the search space generated by current population. Mutation rate was set at relatively high p mut = 10%, aiming to make possible to explore wider part of the search space, but without losing good solutions since NSGA-II is an elitist strategy. The population could be expected to become mature after 3/p mut = 30 generations [52], so that the population was set to evolve for at most 100 generations here.
For each individual variant of the cellular office model, described by the parameters (P 0 , . . ., P 14 ), EnergyPlus reports after simulation: • annual energy used for district heating, which is the heating load H(P 0 , . . ., P 14 ), • annual energy used for district cooling, which is the cooling load C(P 0 , . . ., P 14 ), and • annual amount of electricity used for interior lights, which is the lighting load L(P 0 , . . ., P 14 ).
Since these loads use different types of end energy, they were converted to equivalent source energy. According to [53], district heating has efficiency 0.3 and uses natural gas with source energy factor 1.092, so that the heating load is multiplied by to obtain equivalent source energy. District cooling has COP of 3.0 and uses electricity with source energy factor 3.317, so that the cooling load is multiplied by while the lighting load is multiplied just by the source energy factor c l ¼ 3:317: Thus genetic algorithm was set to minimize the objective function which represents equivalent source energy for heating, cooling and lighting loads, under the constraint that P 0 , . . ., P 14 2 {0, 1, . . ., 8}. Following this optimisation, exhaustive search was performed in a small neighborhood of the best solution found by jEPlus+EA in order to finetune optimal positions of control points.

Convergence of solutions
With the search space consisting of 9 15 % 2.06 × 10 14 shading variants, the first question to tackle is whether the process of genetic optimisation had converged to a local optimum (which need not be a global optimum). Since NSGA-II keeps fittest candidates over different generations, total population created over 100 generations for each model location consisted on average of 4300 shadings that were simulated in EnergyPlus. Positions of control points of the best 200 candidates found during optimisation for each model location are shown in Fig 2. It can be seen from these diagrams that positions of the overhang control points P 4 , . . ., P 10 converge in ten model locations (Miami, Houston, Phoenix, El Paso, San Francisco, Baltimore, Salem, Boise, Vancouver and Fairbanks), but do not converge in the remaining six locations (Memphis, Albuquerque, Chicago, Burlington, Helena and Duluth). While this divergence is in most cases due to oscillation of control point positions between two adjacent values, one can see that P 4 in Memphis, Albuquerque and Burlington, P 9 in Albuquerque, Chicago, Burlington and Helena, and P 10 in Memphis, Burlington and Duluth take on wider range of values in solutions closer to the optimal one found.
Situation is rather opposite with the fin control points P 0 , . . ., P 3 and P 11 , . . ., P 14 : while some of them clearly converge like P 3 in San Francisco, P 11 in Miami or P 12 in Memphis, most of them do not seem to be determined in candidates close to the optimal one. One feasible explanation for this situation is that positions of fin control points have less influence on total equivalent source energy than positions of overhang control points, which forces the genetic algorithm to decide on values of more important parameters first, and leave fine-tuning of less important parameters for future generations. It is thus possible that positions of the fin control points would converge if genetic algorithm were let to run for additional number of generations. However, one should note that running genetic optimisation in jEPlus+EA for additional 100 generations requires approximately 6-7 hours of computing time on a 4-core desktop workstation for each of 16 model locations. In addition, one cannot be certain how many additional generations are needed to achieve convergence: genetic optimisation for Albuquerque was run for additional 100 generations and while the fin control points P 2 and P 3 started to converge, the remaining fin control points did not achieve convergence.
Convergence of positions of control points can also be observed through decreasing standard deviations of positions in best candidate solutions: Table 2 shows standard deviations of positions for the best 200, the best 100 and the best 50 candidates for each model location, together with the average values of positions of control points in the best 50 candidates. Based on information in this table and diagrams in Fig 2, for each model location it was possible to select positions of control points that have converged and to select 2-3 most occurring positions for the remaining control points (except for P 12 in Houston and P 13 in Vancouver, where four alternatives were selected in each case). These selections, shown in Table 3, determine neighborhoods of the best solutions found by genetic optimisation which consist of 64-576 shading variants. An exhaustive search was additionally performed with jEPlus in these neighborhoods, which led to subtle improvements in optimal shading for all model locations other than Vancouver. Table 4 shows positions of control points of the optimal shading variants, together with heating, cooling and lighting loads in equivalent source energy terms for both the model with the optimal shading and the starting model without shading. For easier comparison, these loads are visually represented in Fig 3 as well. Fig 4 futher shows a Sketchup visualisation of the model with optimal shading found for each location.

Load changes from optimal shading
As it can be seen from Table 4, presence of shading increases heating loads in all model locations, from 26.1kWh in Miami up to 300.1kWh in Albuquerque in absolute terms. Apart from Miami, which has negligible heating load, this difference represents relative increase of between 2.0% in Fairbanks and 30.4% in Albuquerque compared to heating load of the unshaded model. Equally expected, presence of shading increases lighting loads in all model locations as well. As lighting loads are more uniformly distributed among different locations, ranging from 820.7kWh in Albuquerque to 916.0kWh in Fairbanks, relative increase of lighting load is also more uniformly distributed, ranging from 3.1% in Duluth to 11.3% in Phoenix. Increases in heating and lighting loads are, however, well compensated by decreases in cooling loads which range from 408.4kWh in Memphis to 1206.4kWh in Albuquerque in absolute terms and represent relative decrease of between 20.6% in Houston and 57.9% in Vancouver compared to cooling loads of unshaded models. When all these loads are taken together, the benefits of shading are still clear: total equivalent source energy is reduced between 202.7kWh  Optimisation of curvilinear external shading of windows H none , C none and L none are loads for a cellular office model without shading, while H opt , C opt and L opt are loads for a model with the optimal NURBS shading found after exhaustive search in the vicinity of the best solution returned by genetic optimisation. All loads are given in equivalent source energy terms, measured in kWh.
https://doi.org/10.1371/journal.pone.0203575.t004 Optimisation of curvilinear external shading of windows in San Francisco and 841.7kWh in Albuquerque in absolute terms, which represents relative decrease of between 3.4% in Burlington and 19.5% in Albuquerque. One should take into account that optimal shadings found here have overhang depths mostly between 1.5m and 2m, which renders them impractical in locations with large annual snowfall such as Chicago, Burlington, Helena, Duluth and Fairbanks. The results reported here are given for completeness of comparisons, while in practice positive effects of shading on total equivalent source energy for these locations would better be achieved by improving performance of glazing in PNNL office building models.

Some unexpected optimal positions of control points
As genetic algorithm is allowed to freely and independently changes values of genes during the process of optimisation, informed only by changes in objective function that determine fitness of candidates, it may easily end up with optimal solutions whose characteristics appear counterintuitive at first. There are a few such unexpected phenomena in optimal shadings visualized in Fig 4: a hole in the upper western part of shading for Miami, holes in the upper western and eastern parts of shading for Houston, and protrusion of the lower part of western fin for El Paso, San Francisco and Baltimore. Western hole in the optimal shading for Miami is clearly related to the position of control point P 4 , which serves as an ending point for both the western fin and the overhang. Table 5 shows heating, cooling and lighting loads of the cellular office model in Miami when the position of P 4 is varied from 0 to 8, while the positions of the remaining control points in this shading are kept intact. As expected, heating and lighting loads increase, while cooling load decreases with an increase in P 4 , since an increase in the value of P 4 also increases the shading area. However, one can see that these loads do not depend linearly on the value of P 4 . The minimal total equivalent source energy is obtained for the maximum value of P 4 , for which the shading loses its hole as shown in Fig 5. Thus, this counterintuitive hole in the optimal shading The remaining control points (P 0 , . . ., P 3 , P 5 , . . ., P 14 ) have positions (3,7,6,8,7,7,7,7,7,8,3,2,6,4). Loads are given in equivalent source energy terms, measured in kWh.
https://doi.org/10.1371/journal.pone.0203575.t005 served to help improve the results of genetic optimisation, whose best candidate solutions clearly had the value of P 4 wrongly converging to 5. Two holes in the optimal shading for Houston are related to positions of the control points P 4 and P 10 , which serve as joint points between the overhang and the fins. Table 6 shows total equivalent source energy for heating, cooling and lighting loads for cellular office model in Houston when positions of these control points are varied from 0 to 8, while positions of the remaining control points are kept intact. As in the case of Miami, when one of these control points is increased while the other one is kept fixed, the shading area is increased so that heating and lighting loads increase, while cooling loads decreases. However, total equivalent source energy ESE as their weighted sum behaves somewhat erratically with respect to individual changes in P 4 and P 10 : for example, for any fixed value of P 4 in Table 6 holds ESE P 10 ¼5 < ESE P 10 ¼6 > ESE P 10 ¼7 < ESE P 10 ¼8 ; while on the other hand, for any fixed value of P 10 holds The minimal total equivalent source energy is obtained for P 4 = 2 and P 10 = 1, which contrary to the case of Miami, shows that genetic optimisation got the value of P 4 correct and made just a minor mistake with P 10 .
Protrusions of the lower part of western fin in the optimal shadings for El Paso, San Francisco and Baltimore are consequences of a large value of the control point P 0 , which serves as the lower end of the western fin, compared to values of its remaining control points. Table 7 shows total equivalent source energy for heating, cooling and lighting loads for the cellular office model in these locations when the position of P 0 is varied from 0 to 8 while the remaining control points in optimal shadings are kept intact. Interestingly, the minimal total equivalent source energy is obtained for P 0 = 8 for El Paso which causes even larger protrusion of the western fin. The optimal value P 0 = 8 for San Francisco has already been found by genetic optimisation, while the fin protrusion is slightly reduced only for Baltimore for which the optimal value is P 0 = 4. However, one should note that differences in total equivalent source energy among shadings considered in Table 7 are less than 0.32%, which opens up an option to obtain shadings of improved design with only minor sacrifices in energy use. The remaining control points (P 0 , . . ., P 3 , P 5 , . . ., P 9 , P 11 , . . ., P 14 ) have values (4,7,5,4,6,8,7,7,8,5,4,3,4). Total equivalent source energy is measured in kWh.

Simpler near-optimal shading designs
Aims for symmetry and simplified lines are natural in building design. Although the present study defined the outer edge of the overhang using seven control points P 4 , . . ., P 10 , and thus gave genetic algorithm the freedom to explore very different shading designs, the resulting optimal shadings turn out to have very similar positions of control points P 5 , . . ., P 9 that define interior part of the overhang, suggesting that genetic algorithm may have tried to equate values of these control points altogether. Actually, the ending overhang control points P 4 and P 10 Table 8 gives the minimum total equivalent source energy for heating, cooling and lighting loads among these shading modification. When compared to optimal shadings from Table 4 (and their improvements for Miami, Houston, El Paso and Baltimore based on Tables 5, 6 and 7), it can be seen that rectification of the overhang increased total equivalent source by at most 0.23% (in Houston), and even slightly decreased it in 7 of 16 locations (Miami, Phoenix, San Francisco, Albuquerque, Vancouver, Duluth and Fairbanks). This gives some support to the expectation that truly optimal shadings have rectified overhang. While window shading has a certain effect on heating and cooling loads by casting shadows on surrounding opaque wall as well, this effect is negligible compared to the impact it has by reducing solar energy incident to the window. With respect to fins, it is important to realize Total equivalent source energy is given above in kWh. The remaining positions of control points (P 1 , . . ., P 14 ) are equal to (4,4,2,7,8,7,8,8,7,8,3, 0, 0, 0) in El Paso, (6,5,2,8,7,8,8,8,7, 0, 0, 0, 0, 0) in San Francisco and (7,3,3,7,8,8,7,8,6,8,3,1,0,0) in Baltimore. https://doi.org/10.1371/journal.pone.0203575.t007 Optimisation of curvilinear external shading of windows that they can be treated independently of each other, as they cannot both shade the window simultaneously.
The eastern fin is less pronounced in Table 4, with positions of many of its control points P 11 , . . ., P 14 close to zero. Importance of the eastern fin is tested by equating positions of the control points P 11 , . . ., P 14 with zero, which effectively removed the lower part of the eastern fin and replaced its upper part with a single arc from the middle of the window to the control point P 10 . For locations where the eastern fin is more pronounced (Miami, Houston, Phoenix, Albuquerque, Salem, Boise, Vancouver and Fairbanks), positions of the control points P 11 , . . ., P 14 are also averaged, resulting in rectification of the lower part of eastern fin. Table 9 shows heating, cooling and lighting loads in equivalent source energy terms for shadings with such modified eastern fins (and already rectified overhang). Comparison with Table 8 yields that equating control points P 11 , . . ., P 14 to zero increases total equivalent source energy by more than 0.1% in Miami, Houston, Phoenix, Albuquerque, Salem, Boise and Vancouver, so in these locations eastern fin may be deemed to have certain importance. Compared to the optimal shadings from Table 4 (and their improvements from Tables 5, 6 and 7), better shading from Table 9 for any location increases total equivalent source energy by less than 0.19%, and slightly decreases it in cases of San Francisco, Vancouver, Burlington and Duluth.
Western fin is significantly more pronounced in optimal solutions as can be witnessed from Fig 4. Unlike the eastern fin which shades from the morning sun when both the external and internal temperatures are still low, the western fin needs to block the sun rays in afternoons when the outside temperature is close to the daily maximum. Three strategies for simplifying line of the western fin are tested here: 1. equating positions of control points P 0 , . . ., P 3 to zero, which replaces the western fin with an arc located in its upper part; 2. replacing positions of control points P 0 , . . ., P 3 with their average, which results in partial rectification of the western fin; All loads are given in equivalent source energy terms, measured in kWh. Relative changes are calculated with respect to the optimal shading listed in Table 4, and their improvements from Tables 5, 6 and 7. All loads are given in equivalent source energy terms, measured in kWh. Relative changes are calculated with respect to the optimal shading listed in Table 4 and their   improvements from Tables 5, 6 and 7. https://doi.org/10.1371/journal.pone.0203575.t009 Total equivalent source energy for heating, cooling and lighting loads is given in kWh. Positions of the remaining control points P 4 , . . ., P 14 for each location are as in a better shading from Table 9. Relative changes are calculated with respect to the optimal shading listed in Table 4 and their improvements from Tables 5, 6 and 7. 3. separately replacing positions of P 0 and P 1 with their average (P 0 + P 1 )/2, and positions of P 2 and P 3 with the average (P 2 + P 3 )/2, which better keeps shape of waving fins.
For each location positions of the remaining control points P 4 , . . ., P 14 are as in better shadings from Table 9. Table 10 gives total equivalent source energy for such modifications of western fins, while Fig 6 shows a Sketchup visualisation of the model with optimal among these shadings for each location.
It is easily seen from Table 10 that the western fin is more important than the eastern fin, as equating positions of control points P 0 , . . ., P 3 to zero may increase total equivalent source energy by up to 4.82% (in Albuquerque). The western fin may be deemed irrelevant in Burlington and Duluth only, where this increase is less than 0.01%. Comparing with the optimal shadings from Table 4, and their improvements from Tables 5, 6 and 7, it can be seen the shadings, with the overhang, eastern and western fins simplified in this way, are no longer better than those found by genetic optimisation. However, increase in total equivalent source energy obtained by simplifying design is at most 0.24% in Houston, while in San Francisco, Vancouver and Burlington it is less than 0.01%. Moreover, visualisations of optimal simplified shadings in Fig 6 clearly show that their outer edges are less winding when compared to shadings from

Conclusion
An important finding of this study is that, although genetic optimisation was given wide freedom in selecting design of outer edge of shading for a cellular office in representative USA climates, the overhang control points in the resulting optimal solutions all had relatively similar positions. This at least partially rejects the initial expectation that curvilinearity of solar paths should induce curvilinearity of the optimal shading shape, as the overhang in optimal shading shapes is close to being rectangular. Another important finding is that the structure of shading shape can be significantly simplified with only a slight increase of up to 0.24% in total equivalent source energy for heating, cooling and lighting loads compared to the optimal shape found by genetic optimisation. Hence instead of 15 control points, with seven control points defining the overhang and five control points for each of the fins (with two common control points for joints of the overhang with the fins), one can achieve very similar results by identifying depths of successive control points, thus dividing the shading into six independent regions: the lower and the upper part of the western fin, joint of the western fin and the overhang, interior part of the overhang, joint of the eastern fin and the overhang, and the whole eastern fin of a shading in a cellular office, although this smaller number of regions requests higher resolution of their feasible depths. This opens up an opportunity to use depths of these regions, together with other parameters such as the window size, glazing properties, and climate characteristics, in a regression analysis in future work and possibly get explanation for the appearance of holes between overhang and fins in the optimal shading for Houston and protrusion of the lower part of western fin in optimal shadings for El Paso, San Francisco and Baltimore.