Life-cycle cost optimization of a solar combisystem for residential buildings in Nepal

ABSTRACT Current solar thermal systems used in Nepal are mostly for domestic hot water heating. Meanwhile, there is a growing need for space heating. A solar heating system can be configured to provide heat for both space heating and domestic hot water, leading to a solar combisystem. Research on residential solar combisystems in the Nepalean context barely exists in literature. This paper intends to propose, model and optimize a solar combisystem for typical single-family houses in Nepal. The optimization problem is formulated to have life-cycle cost as the objective function and a total of 11 variables, which are related to the sizing of combisystem components and the insulation thickness of building envelope. The number of hours not satisfying the thermostat heating set point is treated as the constraint. TRNSYS is used for modeling the solar combisystem. Particle Swarm optimization and the Hooke-Jeeves algorithm implemented in GenOpt are sequentially applied to solve the optimization problem. The results show that the optimization is effective to reduce the life-cycle cost by 66% for the Terai region and 77% for the Hilli region. The findings have demonstrated the importance of building envelope insulation to spread the use of solar combisystems in Nepal.


Introduction
In 2018, the total energy consumption in Nepal was about 565 Petajoule (SAARC Energy Center 2018), nearly 84% of which was contributed by the residential sector (Asian Development Bank 2017). Major residential energy end uses include cooking (60%), lighting (12%), domestic water heating (12%), and space heating (10%), where the numbers in parentheses are the percentages of the residential sector energy consumption in Nepal (Rajbhandari and Nakarmi 2015). Traditional energy sources (e.g., firewood and agricultural residues) are commonly used in many parts of the country. On the other hand, Nepal has abundant sources of renewable energy. In addition to the huge hydro power potential, the country has reliable and rich solar resources. On average, there are 300 sunny days per year and 6.8 sunshine hours per day. The northwest part of Nepal has the highest solar irradiation, reaching up to 5.5 kWh/m 2 /day in global horizontal irradiation. Other regions in the country have average daily solar irradiation in the range between 4.4 kWh/m 2 /day and 4.9 kWh/m 2 /day (ESMAP, Energy Sector Management Assistance Programe 2017). Solar systems are becoming increasingly used for service water heating in hotels, guesthouses, and even in solitary mountain lodges. It is estimated that solar thermal collectors with approximately a total of 100,000 m 2 of surface area have been installed in Nepal, which is equivalent to a collector surface area of 5-7 m 2 per 1,000 inhabitants. The installation of solar collectors was predicted to grow by 50,000 m 2 every year (Peuser, Remmers, and Schnauss 2011). Though solar thermal systems for service water heating have deep penetration in the residential sector, using solar energy for space heating is still rare. Many locations including Kathmandu have a prolonged period of cold weather in winter, with snow falling sometimes. In cold weather conditions, the occupant's thermal comfort can be satisfied with the use of a solar combisystem, which is defined as a solar heating system configured to provide heat for space heating and domestic hot water (DHW) production for a residential household (Weiss 2003).
Many studies on solar combisystems are available in literature. Best practices, standards, and guidelines for solar combisystem design and installation have been proposed through coordinated efforts via several regional or international organizations, representative examples of which include Solar Heating and Cooling Task 26 by the International Energy Agency (IEA), the ALTERNER Programme and the COMBISOL Project by the European Commission. Because the space heating load varies with different locations and building constructions, the best practices and guidelines on solar combisystem design for one location may not apply to other locations. In the remainder of this section, a review of optimization studies on residential solar thermal systems is presented with an emphasis on solar combisystem optimization. Loomans and Visser (2002) applied a genetic algorithm to find the optimal configuration of large solar hot water systems towards the minimization of payback time. Design variables included the collector type, the number of collectors, the heat storage mass and the collector heat exchanger area. Kalogirou (2004) combined the use of an artificial neural network and a genetic algorithm to maximize the life-cycle savings of a solar industrialprocess heat system. Relative to a traditionally trial-and-error method, the optimization resulted in an increase of life-cycle savings by 3.1%-4.9%, depending on whether the backup fuel was subsidized. Lima, Prado, and Montoro Taborianski (2006) used a modified simplex algorithm to optimize a thermosyphon solar water heating system in Sao Paulo, Brazil and they found that the simple payback period could lie in the range between 5.25 and 9.75 years. Fraisse et al. (2009) compared the results from the various single criterion on energetic, exergetic, environmental, and financial performance and multi-objective optimization criteria for a solar domestic hot water system. They cautioned that using a single optimization criterion instead of multiple optimization criteria (e.g. both the cost savings and energy savings) took the risk of obtaining a solution with much worse performance in other unconsidered performance criteria.
Regarding solar combisystem optimization, Bales and Persson (2003) optimized one of the generic combisystem design alternatives from IEA Task 26, which achieved the fractional energy savings of 17% to 30%, depending upon the base case design parameters. Bornatico et al. (2012) optimized a solar combisystem by combining the solar fraction, energy consumption, and cost into a single objective function. They applied two different approaches, a Particle Swarm optimization (PSO) algorithm and a genetic algorithm, and obtained similar optimal solutions. Hin and Zmeureanu (2014) combined the use of PSO and Hooke-Jeeves algorithms towards the optimization of three different objective functions separately: life-cycle cost, life-cycle energy, and life-cycle exergy. Relative to the base design presented by Leckner and Zmeureanu (2011), the optimized system reduced the life-cycle cost by 19%, the life-cycle energy by 34% and the life-cycle exergy by 34%. The work by Hin and Zmeureanu was extended further in (Rey and Zmeureanu 2016) by developing a multi-objective PSO algorithm to minimize the lifecycle cost and life-cycle energy simultaneously. The multi-objective optimization was applied in Montreal, Canada, where the two objective function values were reduced by up to 88.6% and 63.9%, respectively, relative to the base system design as presented by Leckner and Zmeureanu (2011).
Most of the reviewed work adopted the life-cycle cost (LCC) as one of the performance criteria for solar combisystem optimization. Because LCC analysis considers all costs of acquiring, owning, and disposing of a building system, it is useful to compare design alternatives that can achieve the same functionalities but differ with respect to the initial cost and operating cost (Fuller 2016). Theoretically, LCC items related to solar heating systems include capital investment costs, energy costs for operation, maintenance and repair costs, residual values, and possible financial incentives and taxes. However, not all the above cost items are considered in many previous studies.
The algorithms used for solar heating system optimization fall into two categories: traditional numerical methods and artificial intelligence (AI) based optimization methods although different classification approaches are available (Machairas, Tsangrassoulis, and Axarli 2014). Numerical methods (e.g., the Hooke-Jeeves algorithm and classical gradient-based algorithms) are typically less computationally intensive than the AI-based methods (e.g., genetic algorithms, particle swarm optimization algorithms, and ant colony algorithms) but they are limited to local search and the optimal solution depends on the starting point. With the advancement of computing capabilities, AIbased methods are predominantly used in recent literature. In particular, genetic algorithms and particle swarm optimization are the two popular methods to solve solar combisystem optimization problems (Bales and Persson 2003;Bornatico et al. 2012;Hin and Zmeureanu 2014).
Research on residential solar thermal systems performed in the Nepalean context seldomly exists in literature. Timilsina et al. (2018) compared a variable refrigerant flow (VRF) system and a solar floor heating system for a resort in Nepal. The levelized cost of energy of the solar system was found to be Nepalese Rupees (NRs) 13.18 per kWh with subsidies and NRs. 27.47 per kWh without subsidies, which were less economically favorable than the VRF system. Bhattarai et al. (2013) compared a solar air heating system with pebble heat storage to a conventional liquefied petroleum system. The results revealed that the solar system saved the operating energy cost by NRs 40,500 per year for a one-bedroom residence in Kathmandu valley. Although solar systems separately for DHW and space heating are used in Nepal, solar combisystems are new and therefore research is needed to explore the system design suitable for Nepal's climate and residential constructions.
Because of the scarcity of research on solar combisystems in the Nepalean context, this paper proposes a new solar combisystem for residential buildings in Nepal, the performance of which is modeled with TRNSYS software. Formulated to have life-cycle cost as the objective function and a total of 11 variables, the optimization problem is solved with the combined use of a Particle Swarm optimization algorithm and the Hooke-Jeeves algorithm. The variables are limited to those design parameters related to the sizing of combisystem components and the insulation thickness of building envelope. The paper highlights the importance of design optimization and building envelope insulation to spread the use of solar combisystems in Nepal. In the rest of this paper, the solar combisystem and the house prototype are described in Section 2. Then, TRNSYS modeling of the solar combisystem is presented in Section 3. The optimization problem formulation and the optimization approach are discussed in Section 4. Results are analyzed in Section 5. The paper ends with some conclusions in Section 6.

House and solar combisystem description
The house modeled in this study has its floor area and envelope constructions based on a real design by a professional architectural firm in Kathmandu. The house has two and half stories and a total floor area of 232 m 2 . Regarding constructions, the external wall is built from 230-mm-thickness brick, plus 12-mmthickness plaster on both interior and exterior sides. The roof is constructed from 100-mm-thickness reinforced concrete slab with 12-mm cement plaster on the interior side. Following the current construction practice in Nepal, no thermal insulation is used in roof and exterior walls. The window type is singlepane glazing with wooden frame because doublepane glazing is not a common practice (Subedi 2010). Hydronic radiant floor heating is the system used for space heating. Featured with uniform warm surfaces, quiet operation, and no air drafts, radiant floor heating is more favorable for thermal comfort than forced air systems and baseboard heaters. Because of the use of radiant floor heating, the operative temperature, defined as the average of the space air temperature and the mean radiant temperature, is selected for the space thermostat control. The thermostat has its heating set point at 21°C from 7 am to 10 pm (i.e., daytime) and at 18°C from 10 pm to 7 am (i.e., nighttime), with a deadband of 4°C. Figure 1 represents the schematic of the solar combisystem in the house. Solar is the main source of energy for the system. Since the solar availability is intermittent and the magnitude of solar energy depends on weather conditions, it is not realistic to expect solar as the sole energy source in many locations. Therefore, an auxiliary heating mechanism is needed to satisfy the heating load during unfavorable weather conditions. In this work, electricity is used as the auxiliary source of heating because natural gas is not an indigenous source in Nepal and hydroelectricity is widely available. As seen from Figure 1, the solar combisystem consists of three loops: the solar collector loop, the DHW loop, and the space heating loop. These three loops are briefly described below.
The solar collector loop consists of solar collectors, a storage tank, and a circulating pump. To prevent winter freezing, propylene glycol-water solution is the liquid used in this loop. A circulating pump (P1) is used to transport glycol solution between the solar collectors and the internal heat exchanger at the bottom of the storage tank, where the heat transfer from glycol solution to tank water occurs. The on and off control of Pump P1 is based on the solar collector liquid outlet temperature (T sc;outlet , in °C) and the tank water temperature around the heat exchanger (T sc;inlet , in °C). If T sc;outlet is greater than T sc;inlet by more than 10°C, the pump is on and it keeps running until the condition T sc;outlet À T sc;inlet < 2 is met. If T sc;outlet reaches 90°C and beyond, the pump is off for the purpose of safe system operation.
The DHW loop provides domestic hot water to end users. The DHW tank has its thermostat set point at 60° C to avoid legionella growth inside the tank (Weiss 2003). A mixing valve is used to obtain the desired supply-water temperature at 45°C by mixing hot water from the DHW tank and cold water from the city mains or a local water reservoir in the house. Whenever hot water is consumed, the same amount of cold water runs through the storage tank, gets preheated, and is then supplied to the DHW tank. The DHW tank has an embedded electric heater for water heating if the cold water cannot be preheated sufficiently for the DHW needs.
The space heating loop consists of the storage tank, an instantaneous water heater, and a pump (P2). Pump P2 is on when the space operative temperature falls below the heating set point after accounting for the deadband (i.e., 19°C at daytime and 16°C at nighttime). Similarly, Pump P2 is off when the space operative temperature rises above the heating set point after accounting for the deadband (i.e., 23°C at daytime and 20°C at nighttime). Placed downstream of the storage tank, the instantaneous water heater is activated as needed to heat up the water from the storage tank to 50°C.

TRNSYS simulation
There are many building energy performance simulation programs (Crawley et al. 2008;DOE 2020), some of which (e.g., TRNSYS, EnergyPlus, and ESP-r) support the modeling of solar thermal systems. TRNSYS software (Klein et al. 2017) is selected in this study because of the following: 1) it has a rich and extensively validated library of components for solar thermal systems; 2) it supports the dynamic building thermal simulation; and 3) it has the feature of coupling external optimization software for simulationbased optimization. TRNSYS supports modeling of many diversified components, commonly known as TYPEs. A TRNSYS TYPE may represent a physical component found in thermal and electrical energy systems, or a utility routine to handle inputs of weather data, outputs of simulation results, or other timedependent forcing functions. A system model can include many TRNSYS components, the connections between which are represented by their input-output relationships. TRNSYS employs a sequential modelling technique. At each simulation time step, the simulation operates on the system components one by one according to a prescribed path. Table 1 lists the major components of the solar combisystem and their corresponding TRNSYS TYPEs. Certainly, the TRNSYS model includes several other accessory components (e.g., weather file, valves, and controls) not listed in Table 1 for the purpose of brevity. A brief description of these components is presented below while details can be found in (Thapa 2019).
In TRNSYS, TYPE 1b models the thermal performance of flat-plate solar collectors using a standard quadratic efficiency equation. The quadratic equation expresses the collector's thermal efficiency as (Klein et al. 2017): where a 0 , a 1 , and a 2 are correlation coefficients that take the value of 0.769, 13.01 and 0.0489 respectively (Banister and Collins 2015), I T and ΔT refers to the solar irradiance (kJ/hr/m 2 ) on the collector surface and the difference (°C) between the fluid inlet temperature and the ambient air temperature. TRNSYS TYPE 534 models the heat transfer in the immersed heat exchangers and the parasitic heat losses to the surroundings through the tank surface. Vertical temperature stratification of tank water is modeled by a number of nodes, representing the number of equally divided tank water volume. Each node interacts thermally with its neighboring nodes through thermal conduction and fluid movement (either forced movement from inlet flow streams or natural destratification mixing due to temperature inversions in the tank). Increasing the number of nodes leads to a high resolution of temperature stratification but increases the computation time. Five nodes are used in the simulation model in this work. TRNSYS TYPE 158 is used to model the DHW tank. The thermostat is located at the three-fourth of the tank height. The auxiliary electric heater is activated whenever the sensed water temperature is below 60°C. Similar to the storage tank, the DHW tank is modeled with a total of 5 nodes. The DHW profile ( Figure 2) from ASHRAE 90.2 is employed in this work with a daily hot water consumption of 300 L (Fairy and Parker 2004;Weiss 2003).
TRNSYS TYPE 138 is combined with TYPE 14k to simulate the operation of the instantaneous water heater, which heats up the water temperature to 50° C for space heating if the inlet water temperature from the storage tank is lower. TYPE 14k stores the instantaneous water heater availability schedule. If the schedule is indicated as unavailable (i.e., during the non-heating season), the instantaneous water heater is not turned on, irrespective of the storage tank water temperature.
TRNSYS Type 56 is used to model the house. The house is modeled as a single zone building. The floor area, house volume, and window area are based on the real house design described in Section 2 but with the following simplification made during the model development: all windows on each orientation are combined and represented by a big window with the total window area. This simplification does not affect the thermal load but reduces the modeling efforts. The radiant floor with embedded tubes is modeled in TYPE 56 as an active layer above the concrete slab. Figure 3 shows the pictorial view of the simulation studio modeling of the solar combisystem using TRNSYS 18. Table 2 lists the 11 design variables to be optimized, together with their types (i.e., discrete vs. continuous), boundary values (for continuous variables) and discrete sets of values (for discrete variables), and units. Except for the three building-envelope-related variables (i.e., the thicknesses of roof insulation, floor insulation, and wall insulation), all other variables are related to the solar combisystem design and  operation, most of which are component sizing (e.g., collector area, DHW tank volume, storage tank volume, and auxiliary heater capacities). The broad set of variables are thus selected because they all potentially affect the energy performance of the solar combisystem. Considering that most residential buildings in Nepal are not insulated in the field, the buildingenvelope-related variables are included to explore and demonstrate their impact on system sizing and thereby the cost. In Table 2, the lowest value of insulation thickness for roof and external walls is set at 1 mm, which actually represents the case of no insulation because TRNSYS does not allow the construction layer thickness to be set at zero.

Objective function
LCC is the objective function of the optimization problem. As Equation (2) shows, the LCC consists of the initial cost (IC) of the combisystem and building insulation and the energy cost (OC) of the combisystem operation (i.e., the energy for running the two pumps, the electric heater in the DHW tank, and the instantaneous water heater) over the life period of 25 years. The formulation of LCC in Equation (2) has several simplifications. First, replacement cost is not considered. Solar collectors can last 25 years (Hin and Zmeureanu 2014) and hence they do not need to be replaced in the middle of the project life period. Components that possibly need to be replaced include water tanks and pumps. All design alternatives have the same number of water tanks and pumps but may have different sizes. The replacement cost difference due to sizes is negligible because 1) in practice, the labor cost (which is a major part of replacement cost) of pump and tank replacement does not vary with the equipment sizes within a certain range; 2) after being discounted to the present value, the difference of displacement equipment costs in 10-15 years (Hin and Zmeureanu 2014) becomes trivial. Second, maintenance cost is assumed to be the same for all design alternatives and therefore is not considered. Lastly, salvage values at the end of the project life are minimal and thus not considered either.
where X represents the vector of design variables to be optimized (see Table 2). Table 3 provides the unit costs of the major elements contributing to the initial cost in Equation (2). The data are from different sources. The unit costs of solar collector, water tanks, and pumps are provided by a solar thermal system installer in Nepal and they are the installed costs, including material/equipment and labor costs. The cost of an instantaneous water heater depends on its rated power. A limited set of data based on the Rheem products are used to develop the following regression equation that correlates the cost (in $) and the rated power (in kW) of instantaneous water heaters: where the symbol PW IWH refers to the rated power of instantaneous water heater (IWH) in kW. The IWH cost in Equation (3) refers to the equipment cost only. It is reasonable to assume that the installation cost does not change with the IWH sizes, hence ignoring the installation cost has no impact on the optimal solution. The operating cost in Equation (2) is calculated as where AC X ð Þ is the annual energy cost of operating the combisystem, a is the effective interest rate, and n is the project life (i.e., 25 years).
The annual energy cost is calculated by multiplying the electricity price and the annual electricity consumption of the combisystem obtained from TRNSYS simulation. The average electricity price for residential buildings is 0.096 USD per kWh in Nepal (Nepal Energy Forum 2020).
The effective interest rate considers the time value of money during the analysis period. It is calculated as: where e is the escalation rate of electricity price and i d is the discount rate after accounting for the price inflation. It has been observed that for the past halfdecade, the escalation rate of electricity price is 3.44%  per year (Nepal Energy Forum 2020). The discount rate (i d ) is set at 9.71% (Nepal Rastra Bank 2019).

Constraint
The optimization is essentially a process to search for the combination of variable values that minimize the life-cycle cost. If no constraint is applied, it is possible that a combisystem design with the smallest collector area, the smallest storage size and the smallest instantaneous water heater is selected because that combination leads to the lowest initial cost. However, the design is very unlikely satisfactory to meet the space heating load and may cause many hours not meeting the thermostat setpoint. To address the above issue, the number of hours under the heating set point (HUSP) needs to be specified as a constraint of the optimization problem (Hin and Zmeureanu 2014): In the above equation, the HUSP defines the total number of hours that the space operative temperature falls below the heating set point over the entire heating season (from November to March). An upper limit of 300 hours is used in this paper by referring to the appendix G of ASHRAE Standard 90.1. Although ASHRAE Standard 90.1 is for commercial buildings, the number specified there can be reasonably used as a reference for residential simulation.

Optimization approach
The explicit constraint (Equation 6) presents a challenge to solve the formulated optimization problem. The penalty function approach is employed to convert the original constrained problem into an unconstrained one. The new objective function with the added penalty term is then formed as: where the expression GT HUSP X ð Þ; 300 ð ÞÞ takes the value of 1 if HUSP X ð Þ > 300 and the value of 0 otherwise, and 200,000 is an arbitrarily large number to penalize any design alternatives violating the constraint.
GenOpt, an open source generic optimization program specially designed for optimization problems that rely on building simulations to evaluate the objective function (Wetter 2016), is used to solve the optimization problem. GenOpt enables the coupling between TRNSYS and the optimization algorithms via the use of readable text files. GenOpt has a library of local and global optimization algorithms. Because the formulated optimization problem has both continuous and discrete variables and the performance space shape is unknown, the hybrid use of Particle Swarm Optimization (PSO) and Hooke-Jeeves optimization (HJ) is an appropriate choice, as demonstrated by several previous studies (e.g., Lee and Cheng 2012;Hin and Zmeureanu 2014). PSO is applied first to perform a global search of the optimal solution. Using the optimum found from the PSO as the starting point, the HJ optimization is then applied to perform a local search towards the LCC minimization.
PSO is a stochastic optimization algorithm that uses a set of potential solutions to perform a global search of optimal one for non-linear optimization problems. Each such solution is called a particle and the set of particles is called a population. At the end of each generation, each particle is updated based on its own (local) best objective function value ever found in the evolution process and the global best value of all particles in its neighborhood. A continuous variable of the i-th particle is updated with the following equations (Wetter 2016): ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi In Equations (8)-(11), x indicates a continuous variable withx max andx max being its upper boundary and lower boundary, respectively, t is the generation index, υ is the particle's velocity, ρ is a uniformly distributed random number between 0 and 1, c 1 is the cognitive acceleration constant, c 2 is the social acceleration constant, � is the constriction coefficient, κ is the constriction gain, λ is the maximum velocity gain for continuous variables, φ ¼ c 1 þ c 2 , p l and p g is the variable corresponding respectively to the local and global best objective function.
In the PSO algorithm used in this study, a discrete variable is encoded with a binary string. Let ψ i ¼ 0; 1 f g m be the binary representation of the discrete variable y i , where m indicates the binary string length. The binary string is updated every generation according to the following equations (Wetter 2016): In Equations (12)-(15), the superscript j is the position index of the binary string, i.e, j 2 1; 2; . . . ; m ð Þ, v max is the maximum velocity gain for discrete variables, q l and q g is the binary string of the variable corresponding respectively to the local and global best objective function. All other symbols have the same meanings as those defined for continuous variables in Equations (8)-(11).
The HJ optimization consists of a sequence of exploratory moves about a base point which, if successful, are followed by pattern moves. The exploratory moves are performed by perturbing the current point (X k ) by a small amount (positive or negative) for each variable and observing whether the objective function value LCC X k ð Þ improves or worsens. After all the variables have been considered, a new point (X kþ1 ) will be reached. If X kþ1 ¼ X k , the step length is reduced by dividing a parameter value called mesh size divider and the exploratory procedure is repeated. If X kþ1 �X k , a pattern move is made. The pattern move attempts to speed up the search by using the information on the best search direction. The next pattern point is given by where a is a positive acceleration factor with a common choice of value being 2. The objective function is evaluated at the pattern solution X k . If LCC X kþ2 ð Þ < LCC X kþ1 ð Þ, continue the pattern move after updating the points X k and X kþ1 ; otherwise, if LCC X kþ2 ð Þ � LCC X kþ1 ð Þ, repeat the process by starting the exploratory search from the new point X kþ1 . The algorithm terminates when the step length is reduced a predefined number of times (i.e., number of step reductions).
The GenOpt software allows the hybrid application of PSO and HJ in a single setup by specifying the number of simulations as one of the parameters to manipulate simulation-based optimization. The simulation number, however, includes the number of simulations called by both PSO and HJ. It is thus difficult to control how much efforts are taken by each of the optimization algorithms. To address this problem, the application of PSO and HJ are applied sequentially in two separate steps. The first step is to run the PSO which uses the number of generations as the stopping criterion. After the PSO is complete, the obtained optimum is manually specified as the starting point for the HJ, which uses the number of simulations as the stopping criterion. Table 4 shows the optimization parameters used in this paper. The coupling of GenOpt and TRNSYS needs a number of text input files, namely the simulation template file, the command file, the configuration file, and the initialization file. Based on the input files, GenOpt launches the TRNSYS simulation program, reads the function value being minimized from the simulation result file, checks possible simulation errors and determines a new set of input parameters for the next run. The output file contains the variable values and the corresponding LCC results for all iterations. Details of these GenOpt input and output files can be found in the user manual (Wetter 2016).

Results and analysis
The solar combisystem optimization has been performed for two regions in Nepal: Terai region and Hilli region. The Terai region lies in the south of the country, with a width of about 26 km to 32 km and an altitude of about 60 m to 305 m from the mean sea level. This region has a warm winter, with temperature ranging from 8°C to 23°C. Biratnagar, Janakpur, and Birganj are representative cities of this region. However, no weather file corresponding to any of the cities in the Terai region is available in TRNSYS software. Hence, the weather file for the city of New Delhi in India, which has similar climate as the Terai region, is selected for TRNSYS simulation. The Hilli region lies to the north of the Terai region and has an altitude ranging from 305 m to 3000 m above the sea level. Kathmandu, the capital of Nepal, belongs to this region and it will be used as the representative location in TRNSYS simulation. Kathmandu has an average daily temperature from 2°C to 12°C in winter. Figure 4 shows the evolution of the life-cycle cost over the 50 generations from the PSO. For each generation, the objective function values for all 15 solutions are plotted except for one outlier solution generated at the 18 th generation. The outlier solution has a very high LCC (~$219,000) because it has triggered the penalty function. In Figure 4, solutions in the  ,710 USD and 16,390. USD Because of the stochastic nature of PSO, the best solution may worsen from one generation to another, which can be apparently observed during the first 15 generations. After the 25 th generation, the optimization is almost converged and the best solution has changed only 50 USD in the last 25 generations, which implies that the maximum number of generations can be reduced to save computation time.
Using the best solution found from the PSO as the starting point, the HJ optimization algorithm is applied. Because the HJ algorithm is not able to deal with discrete variables, the variables on building envelope insulation, solar collector area, DHW tank auxiliary power that have been considered discrete during PSO, are specified with their optimal values found from PSO.
The HJ algorithm runs for a total of 15 iterations. Figure 5 shows that the life-cycle cost is reduced from 14,707 USD to 9,644 USD at the end of the HJ optimization. Because the HJ optimization is a local search algorithm that maintains only one solution, the LCC does not worsen from one iteration to the next iteration.  of the optimization algorithm. During the evolution, several individuals clearly deviate away from the converged region, which is normal for PSO because of its global exploration of the solution space. Using the best solution found from the PSO as the starting point, the HJ algorithm further reduces the LCC value from 17,812 USD to 12,515. USD Table 5 lists the initial and optimal values of the 11 variables obtained at the end of the optimization for both regions. The initial values are defined based on design guidelines, conventional construction practices, and engineering judgment. The initial designs are technically feasible and do not violate the constraint on thermal comfort (Equation 6). Note that the water flow rate of the space heating loop is the only variable that takes different initial values to address the space heating load in the two regions. Table 5 leads to the following major observations: • The solar collector area is reduced from the initial value of 20 m 2 to 8 m 2 and 14 m 2 , respectively for the Terai region and the Hilli region.
The system is sized a smaller optimal collector area in the Terai region because it has a smaller heating load than the Hilli region. The solar collector slope for both regions is found to be about double of the latitude (28°). Hin and Zmeureanu (2014) pointed out that any slope in the range between the location latitude and 15° above the latitude can be regarded as a reasonable design. To better understand the impact of collector slope on LCC, a sensitivity analysis has been made by perturbing the collector slope from 20º to 65º in a step of 5º while keeping all other variables the same value as the obtained optimal solution for the Hilli region.
The results of sensitivity analysis show that the LCC changes slightly (~$14), which indicates that any value of collector slope in the range between 20º and 65º is reasonable. • The sizes of both the storage tank and the DHW tank are decreased from the initial design. This is mainly for two reasons: 1) both tanks have high capital cost while the electricity price is very low; and 2) the solar collector area is small, and it is not necessary to have a large storage tank. • The solar collector fluid flow rate is increased from 200 kg/hr for the initial design to the uppermost boundary for both regions. The most possible reason is that the increase of the flow rate slows down the increase of the storage water temperature, which is favorable for the solar collector efficiency. • The flow rate of the space heating loop is optimized to near its lower boundary for the Terai region but to a much higher value for the Hilli region. This flow rate interacts with the space heating load and the instantaneous water heater size. Using insulation in building envelope significantly reduces the space heating load, which drives the decrease of the water flow rate from the initial design. • Building-envelope-related variables have a large impact on the optimal design. Adding insulation to external walls, roof and floor decreases the space heating load. The reduction of heating load triggers the use of smaller instantaneous heater and solar collectors, which have high initial cost. Meanwhile, adding insulation costs much lower than the use of a large combisystem. Therefore, it is reasonable to expect that relative to the initial design, a higher level of insulation is used in the envelope, as demonstrated from the optimal insulation values in Table 5. Both roof insulation and wall insulation are increased from the minimum of 1 mm (a proxy of no insulation because TRNSYS does not accept a zero value for the layer's thickness) for the initial design to the maximum of 50 mm for the optimal design. The optimal thickness of ground insulation remains at its initial value of 30 mm for the Terai region but doubles the initial value for the Hilli region.

Conclusions
The work presented in this study includes the design, modeling, and optimization of a solar combisystem for typical single-family houses in Nepal. Life-cycle cost, consisting of the initial cost of the solar combisystem and the energy cost of operating the system over 25 years, is the objective function to be minimized. Major optimization variables include the thermal collector area, the DHW tank size, the storage tank size, the flow rate of water circulating between the thermal collectors and the storage tank, the auxiliary electric heater sizes, and the thickness of insulation used in external walls, roof and floor. The number of hours not satisfying the thermostat set points is treated as the constraint of the optimization problem. The Particle Swarm optimization algorithm and the Hooke-Jeeves algorithm are combined to solve the optimization problem.
The solar combisystem has been modeled and optimized in two different climatic regions in Nepal: the Terai region and the Hilli region. For the Terai region, it is favorable to have a system with a small solar collector area of 6 m 2 and the DHW and storage tanks with a volume of about 130 L. For the Hilli region, the optimal design has a collector area of 14 m 2 and larger tanks (150 L ~ 170 L). Major findings from this work and their implications include the following: • Optimization is an effective approach to explore the design space for LCC minimization. The optimal solution has reduced the LCC of solar combisystem by 66% in the Terai region and by 77% in the Hilli region, relative to the corresponding initial design based on guidelines, conventional construction practices, and engineering judgment. • The optimization has demonstrated that building envelope insulation is highly important towards LCC minimization of solar combisystems. Therefore, it is recommended to include buildingenvelope-related variables (e.g., thermal insulation thicknesses) in solar combisystem optimization. Given the practice of not using thermal insulation in residential buildings, energy-efficient building codes are needed to catalyze the deployment of solar thermal systems for space heating in Nepal.
This paper considers LCC optimization only. It is worthwhile to pursue optimal designs based on other performance criteria such as life-cycle energy and life-cycle environmental impact assessment and compare the optimal solutions from different criteria. Multi-objective optimization is a promising approach to be considered in future when the trade-off between conflicting performance criteria needs to be addressed.

Disclosure statement
No potential conflict of interest was reported by the author(s).