This paper examines the effects on climate and air pollution of open biomass burning (BB) when heat and moisture fluxes, gases and aerosols (including black and brown carbon, tar balls, and reflective particles), cloud absorption effects (CAEs) I and II, and aerosol semidirect and indirect effects on clouds are treated. It also examines the climate impacts of most anthropogenic heat and moisture fluxes (AHFs and AMFs). Transient 20 year simulations indicate BB may cause a net global warming of ~0.4 K because CAE I (~32% of BB warming), CAE II, semidirect effects, AHFs (~7%), AMFs, and aerosol absorption outweigh direct aerosol cooling and indirect effects, contrary to previous BB studies that did not treat CAEs, AHFs, AMFs, or brown carbon. Some BB warming can be understood in terms of the anticorrelation between instantaneous direct radiative forcing (DRF) changes and surface temperature changes in clouds containing absorbing aerosols. BB may cause ~250,000 (73,000–435,000) premature mortalities/yr, with >90% from particles. AHFs from all sources and AMFs + AHFs from power plants and electricity use each may cause a statistically significant +0.03 K global warming. Solar plus thermal‐IR DRFs were +0.033 (+0.027) W/m 2 for all AHFs globally without (with) evaporating cooling water, +0.009 W/m 2 for AMFs globally, +0.52 W/m 2 (94.3% solar) for all‐source BC outside of clouds plus interstitially between cloud drops at the cloud relative humidity, and +0.06 W/m 2 (99.7% solar) for BC inclusions in cloud hydrometeor particles. Modeled post‐1850 biomass, biofuel, and fossil fuel burning, AHFs, AMFs, and urban surfaces accounted for most observed global warming.

1 Introduction This paper examines the effects of all (anthropogenic plus natural) contemporary open biomass burning (BB) on climate accounting for BB heat fluxes, BB moisture fluxes, black carbon (BC), tar balls, other brown carbon (BrC), other aerosol particles and gases, cloud absorption effects (CAEs), semidirect effects (SDEs), and indirect effects (IEs). Cloud Absorption Effects I and II are the effects on cloud heating of, respectively, (I) absorbing inclusions within hydrometeor particles and (II) absorbing aerosol particles interstitially between hydrometeor particles at the relative humidity (RH) of the cloud (as opposed to the RH outside of the cloud) while accounting for cloud particle scattering to the absorbing aerosol particles [Jacobson, 2012]. Absorption by aerosol particles interstitially between cloud drops increases due to increased optical focusing with increasing aerosol water content, which itself increases with increasing RH, so it is important to determine the water content of interstitial aerosol particles at the RH of the cloud for determining aerosol absorption between cloud drops. To account more fully for CAE II, it is also important to run separate radiative transfer calculations for each clear and cloudy skies in a grid cell. Several analytical and box model studies have considered CAE I [e.g., Danielson et al., 1969; Chylek et al., 1996]. Several 3‐D studies have calculated the radiative forcing [e.g., Chuang et al., 2002; Zhuang et al., 2010; Ghan et al., 2012] or climate response [Jacobson, 2006, 2010a; Erlick et al., 2006] due to BC inclusions within cloud drops and/or ice crystals. Some studies have simulated the impacts of BC interstitially between cloud drops, but treated them as if their coating were at the RH outside of clouds rather than of the cloud [e.g., Ackerman et al., 2000] or treated them at a grid‐cell‐averaged relative humidity [e.g., Ghan et al., 2012]. In the latter case, this treatment would converge to the treatment of BC at the RH of the cloud at a model horizontal resolution fine enough to resolve individual clouds (cloud‐resolving model) but would underestimate the absorption when clouds are not resolved (e.g., in a regional or global model). Thus, while some models have treated CAE II partially, no global or regional climate response study has accounted fully for CAE II, where interstitial aerosol absorption is determined at the RH of the cloud and cloudy sky radiative transfer is calculated separately from clear‐sky radiative transfer. Radiative properties of this effect were isolated with a 3‐D global model in Jacobson [2012]. CAE II with these treatments is accounted for here in a 3‐D climate response simulation and in a 1‐D model analysis of the time‐dependent formation and burn‐off of a single cloud. The semidirect effect is the change in cloudiness due to the decrease in near‐cloud RH and increase in atmospheric stability caused by absorbing aerosol particles below, within, or above a cloud. Indirect effects are the increase in cloud reflectivity (first indirect effect) and the decrease in precipitation thus increase in cloud liquid water content and lifetime (second indirect effects) due to the addition of anthropogenic aerosol particles to an evolving cloud. Previously, Jacobson [2004] calculated, through global numerical simulation, that BB particles may cause short‐term cooling, but its gases cause warming that overpowers the cooling over the long term. Intergovernmental Panel on Climate Change (IPCC) [2014] similarly reported the radiative forcing of biomass burning aerosol particles as −0.2 (−0.6 to −0.07) W/m2 in the Second and Third Assessment Reports, +0.03 (−0.09 to +0.15) W/m2 in the Fourth Assessment Report, and 0.0 (−0.2 to +0.2) W/m2 in the Fifth Assessment Report (AR5). However, all these studies and others [e.g., Naik et al., 2007; Sakaeda et al., 2011] neglected the following processes: (1) absorption by tar balls (strongly absorbing brown carbon), (2) anthropogenic heat fluxes (AHFs) and anthropogenic moisture fluxes (AMFs), and (3) cloud absorption effects I and II. In addition, IPCC [2014] considered only BC and OC as aerosol constituents, ignoring the vast number of other primary aerosol constituents and secondary aerosol precursor gases, including UV‐ and visible‐absorbing gases, such as those listed in Jacobson [1999]. Chung et al. [2012] used an observational approach to approximate the contribution of BrC as ~20% of total aerosol absorption at 550 nm, most of the rest of which was due to BC. They also found that BrC absorption offset the scattering due to all organic carbon, suggesting that radiative effects in biomass burning regions may be dominated by the effects of BC. Kirchstetter and Thatcher [2012] similarly estimated that about 14% of solar absorption in California wood smoke integrated from 300 to 2500 nm was due to organic carbon and the rest, BC. Jacobson [2012] calculated a strong contribution to the global mass absorption coefficient from tar balls. China et al. [2013] found the mass concentration of tar balls in a fire in New Mexico to be 10 times that of soot particles. Simulations here suggest that tar balls alone from biomass plus solid biofuel burning may contribute up to 0.07 K of global warming. Ten Hoeve et al. [2012] compared satellite‐derived with modeled (with the present model at high spatial resolution over Brazil) curves of aerosol optical depth (AOD) versus cloud optical depth (COD) and found that while biomass burning aerosol particles increase CODs with increasing AOD below ~0.2 due to stronger indirect effects than CAEs plus SDEs, BB aerosol particles decrease CODs at higher AODs due to stronger combined CAEs and SDEs than indirect effects. These results suggest that in regions of high BB emissions, BB particles should reduce COD, increasing solar transmission to the surface, warming the surface in a positive feedback. Anthropogenic heat fluxes (AHFs) are energy emissions resulting from fossil fuel, biofuel, and biomass burning (outdoor and in power plants), nuclear reaction, electric power use, frictional dissipation of motive energy, and human respiration. Anthropogenic moisture fluxes (AMFs) are water vapor fluxes due primarily to fossil fuel, biofuel, anthropogenic biomass combustion, evaporation of water used to cool thermal power plants and buildings, and human respiration. Another major anthropogenic water vapor source is evaporation of water pumped from underground aquifers for human consumption or agriculture. A sink is a reduced transpiration caused by deforestation without burning. However, AMF changes due to evaporated aquifer water and transpiration changes from deforestation were not treated here due to the difficulty in quantifying them spatially and globally. Several studies have examined anthropogenic water vapor emissions in urban areas [e.g., Hage, 1975; Clark et al., 1985; Lee, 1991; Moriwaki et al., 2008; Sailor, 2011]. Hage [1975] suggested that AMFs in an urban area explained observed increases in absolute humidity relative to a rural area. Moriwaki et al. [2008] estimated urban AMFs from energy consumption statistics evaluated against measurements. Whereas, all climate models account for natural moisture fluxes, global studies rarely consider the impact of AMFs. A reason for its perceived unimportance is that nonaquifer AMFs are ~8800 times smaller than natural moisture fluxes. This is the ratio of global precipitation (2.69 mm/d = 5 × 108 Tg‐H 2 O/yr), to the global AMF from Table 2 (56,680 Tg‐H 2 O/yr). Yet considering that water vapor may be responsible for ~66% (Table 1) of the natural greenhouse effect of ~33 K, and dividing by 8800 suggests that AMFs may cause a first guess of 0.0025 K (~0.31% to 0.28%) in net anthropogenic global warming (AGW) of ~0.8 to 0.9 K to date. However, this result does not consider enhanced cloudiness and reflectivity due to the water vapor, which likely offsets some or much of this additional warming. Also, it does not account for the fact that water vapor absorption bands are already largely saturated, so the addition of AMFs, which are incremental moisture fluxes, may not contribute to warming so much as if no other moisture existed in the atmosphere. A more accurate calculation of AGW due to AMFs would be with global simulations with and without AMFs alone. Such simulations were performed here; however, temperature responses were found not to be statistically significant (Table 6). Table 1. Direct Radiative Forcings for Naturally Emitted Gases From the GATOR‐GCMOM Model Over 694 Wavelength/Probability Intervals Calculated by Averaging Top‐of‐Atmosphere Net Downward Irradiance Differences When a Chemical was Present and Absent, Respectively, Globally and Over Four Times During a Single Contemporary Model Year: 12 GMT on 1 January, 1 April, 1 July, and 1 October Chemical Mixing Ratio (ppmv) Natural Direct Radiative Forcing (W/m2) Percent of Total H 2 O Variable 79.4 66 CO 2 275 30.4 25 CH 4 0.715 0.72 0.60 N 2 O 0.27 1.6 1.4 O 3 Variable 7.4 6.2 O 2 209,500 0.73 0.60 CO Variable 0.039 0.032 CH 3 Cl Variable 0.0035 0.0029 Total 120.3 100 With respect to AHF's, Ashworth [1929] hypothesized that AHFs from factories operating from Monday through half of Saturday enhanced convection, thereby explaining higher local rainfall during the week and minimum rainfall on Sundays, when the factories were closed. Washington [1972] carried out an 80 day global climate response study of the effects of uniformly distributed AHFs on near‐surface temperatures but could not distinguish the signal from noise of the model, which was run at coarse horizontal and vertical resolutions and with a fixed sea‐surface temperature. Block et al. [2004] simulated the climate response of AHFs from fossil fuels over central Europe with a regional climate model and found a strong signal. Flanner [2009] quantified realistic AHFs from recent fossil fuel and solid biofuel combustion, but not open biomass burning. In this study, the global climate responses of spatially distributed AHFs from fossil fuels, solid biofuels, and open biomass burning are calculated for an extended period with ocean temperatures and circulation responding to climate forcing. Simulations of the combined effects of AHFs and AMFs from power plants and, separately, transportation are examined. Statistical significance tests are performed.

2 Global Gridded Emission Inventory of Anthropogenic Heat and Moisture Fluxes The realistic simulation of AHFs and AMFs requires the development of global gridded emission for each. For this work, such emissions were calculated. Table 2 summarizes the resulting global nonaquifer AHFs and AMFs for 2005. They were determined as follows. Table 2. Estimated 2005 Global Emissions of Heat, Carbon Dioxide, and Water Vapor to the Atmosphere From Anthropogenic Sources a TWh/yr Heat From Electricity Use or Frictional Dissipation of Motive Energy b TWh/yr Heat to Air Not Used to Evaporate Water From Combustion or Nuclear Reaction c Tg‐CO 2 /yr From Combustion or Bacteria or Metabolism d Tg‐H 2 O/yr From Combustion or Metabolism e Tg‐H 2 O/yr From Evaporation Public Electricity/Heat Generation Coal 6,988 4,873 7,568 3,378 12,430 Natural gas 3,520 992 2,412 1,976 6,262 Fuel oil 943.6 588 909.9 334.9 1,680 Nuclear 2,640 1,909 4,696 Hydroelectric 2,896 Geothermal 54 Wind 99.7 Solar 4.2 Biomass 220 342.2 330.1 112.6 391.5 Total 17,365 8,704 11,220 5,801 25,460 Private Electricity, Manufacturing Industries, and Construction Coal 1,634 1,140 1,770 790.2 2,908 Natural gas 3,688 1,040 2,528 2,070 6,562 Fuel oil 1,343 836 1,294 476.6 2,389 Total 6,666 3,016 5,592 3,337 11,858 Industrial processes + product use 12,570 3,910 1,745 Oil and gas production 2,660 827.8 369.5 Onroad transportation 3,260 14,400 4,633 1,754 Rail, inland ships, other nonroad 338.6 1,194 411.8 151.7 International shipping 690.5 2,434 839.7 309.4 Aircraft 498.6 1,757 594.5 233.1 Fossil fuel fires 152.6 47.49 21.20 Waste burning 96.5 30.03 10.24 Biofuel burning 25,608 8,756 2,987 Open biomass burningg 22,031 7,738 2,640 Agricultural soils 112.1 Total 28,819 94,620 44,710 19,360 37,320 Human respiration 7,150f 2,360f 965f Volcanos 562 8.96 25.15 First, worldwide 2005 public electric power generation estimates (E elec , TWh/yr) for coal, natural gas, fuel oil, nuclear, hydroelectric, geothermal, wind, solar, and biomass power were obtained from Energy Information Administration (EIA) [2012] (Table 2). Gridded global annual CO 2 emissions (E CO2 , Tg‐CO 2 /yr) for 2005 at 0.1° × 0.1° resolution were then obtained from Emission Database for Global Atmospheric Research (EDGAR) [2010] in the following sectors: fossil fuel public electricity and heat production, fossil fuel private electricity and heat production, industrial processes and product use, oil and gas production, fuel production, onroad transportation, nonroad transportation (rail, inland shipping, and other nonroad), fossil fuel fires, soil waste incineration, and agricultural soils. Emissions from the public and private electric power sectors were partitioned among coal, natural gas, and fuel oil (G. Maenhout, personal communication, 2011). CO 2 emissions for biomass electric power combustion, not included in the EGDAR terms, were estimated from biomass electric power generation as E CO2 = E elec C ef , where the emission factor for biomass was calculated in Table 3 as C ef = 1.5 g‐CO 2 /Wh. Table 3. Parameters Used or Derived in Equations to Obtain the Results in Table a Combustion Mole Ratio H 2 O:CO 2 M b Combustion Emission Factor (g‐CO 2 /Wh Electricity Generated) C ef c Carbon Content of Fuel (g‐C/kg Fuel) f C d Thermal Efficiency for Electric Power Production Only Based On H L (H H ) η th e Tank‐to‐Wheel Efficiency η tw f Lower Heating Value (MJ/kg Fuel) H L f Higher Heating Value (MJ/kg Fuel) H H Coal 1.09:1 1.08 0.343 (0.325) 22.73 23.97 Natural gas 2.0:1 0.685 0.40 (0.36) 47.10 52.23 Fuel oil 0.9:1 0.964 0.352 (0.33) 42.69 45.54 Nuclear 0.34 Biomass 0.83:1 1.5 500 0.265 (0.25) 16.33 17.34 Jet fuel 0.952:1 862 0.22 43.10 47.10 Diesel 0.9:1 870 0.22 42.79 45.77 Gasoline 0.935:1 855 0.17 43.45 46.54 Gas/diesel mix 0.925:1 859 0.18 43.20 46.30 Global gridded CO 2 emissions from anthropogenic plus natural open biomass burning at 8 day time resolution for a 5 year period were obtained from fire data and processed to the actual model horizontal resolution as described in the footnote to Table 5. This method resulted in annually averaged global biomass burning CO 2 emissions of 2.2 Pg‐C/yr (Table 2), which compares with an estimate of 1.5 Pg‐C/yr from Houghton [2005] and 2–4 Pg‐C/yr from Seiler and Crutzen [1980]. Seiler and Crutzen [1980] further estimated that wildfires in temporal plus boreal forests comprised ~35.7% of all dry matter burned but pointed out that 70–90% of such wildfires were due to human activity (e.g., campfires, debris burning, cigarettes, etc.). Thus, of total world fire emissions today, ~7.1% (3.6%–10.7%) may be natural and the rest, anthropogenic. Houghton [2005] indicates that, in 1850, CO 2 emissions from land use change may have been ~34% those in 2005. Thus, in 1850, BB emissions were lower than those today. Such emissions may have been mostly anthropogenic [e.g., Marlon et al., 2008] although partly natural as well. Today, BB emissions are much higher with only a small percent natural. This study quantifies the impacts of both anthropogenic and natural open biomass burning, but the natural component is only ~7% of the total. Gridded CO 2 emissions from solid biofuels were obtained as described in the footnote to Table 5. Aircraft emissions of CO 2 were found from flight‐by‐flight fuel burn data for 31 million flights per year at an hourly time resolution at the location of each segment of each flight [Jacobson et al., 2011]. International shipping emissions were gathered from the 0.5° × 0.5° resolution Intergovernmental Panel on Climate Change (IPCC) AR5 2005 inventory under the Representative Concentration Pathways 4.5 trajectory [Clarke et al., 2007] by back‐calculating fuel use from NO x emissions and then applying a CO 2 emission factor to the fuel use. 2 emission factor (g‐CO 2 /Wh) for fossil fuel CO 2 ‐generating electric power sources (coal, gas, and oil) was first estimated from the data in Table (1) All electric power generated is eventually converted to heat, either through transmission and distribution losses or end‐use conversion. Even mechanical motion produced by electricity dissipates by friction to heat. As such, the heat due to electric power generated in Table 2 equals the actual power generated. To obtain gridded heat release due to electric power generation, a rough global COemission factor (g‐CO/Wh) for fossil fuel CO‐generating electric power sources (coal, gas, and oil) was first estimated from the data in Table 2 with The results are given in Table 3. These factors are designed to ensure the consistency between Emission Database for Global Atmospheric Research (EDGAR) emissions and Energy Information Administration electric power generation. The value for natural gas 0.685 g‐CO 2 /Wh suggests the use of natural gas plants for both peaking power (0.65–1.17 g‐CO 2 /kWh) and base load power (0.4–0.5 g‐CO 2 /kWh). Many natural gas plants also run in spinning reserve mode, not producing any power, but emitting about 15–19% of the CO 2 that the plant would emit when generating power [Fripp, 2011]. Gridded CO 2 emissions were then divided by these emission factors (plus an emission factor due to biomass power plants for electricity) to obtain gridded electric power generation and heat release. Although electric power is not usually used where it is generated, transmission distances to end users are generally less than a few hundred kilometers, less than the width of global model grid cells used here. For heat release due to nuclear power and renewable energy sources, the global heat released from end‐use electric power was distributed proportionally to the same locations as the heat released from other fossil fuel sources, although this causes a minor spatial inconsistency since many countries did not use nuclear power or renewable energy in 2005. (2) η th is the thermal efficiency based on the lower heating value of the fuel, since combusted water is assumed to be emitted as a vapor. The thermal efficiency (Table L e = 2465 J/g is the latent heat of evaporation, V c = 1779 g‐H 2 O/kWh electricity (0.47 gal‐H 2 O/kWh) is the average cooling water evaporated per kWh of electric power generated from U.S. thermal power plants [Torcellini et al., 2003 A j = 3.6 × 106 J/kWh. Cooling water is additional water from a nearby lake or river or the ocean used to absorb some of the heat released that is not already used to evaporate water produced during combustion. Heat released during combustion or nuclear reaction was calculated accounting for the thermal efficiency of thermal power plants and for the fact that a portion of the heat goes toward evaporating cooling water. For thermal power plants (e.g., natural gas, oil, biomass, and nuclear), the net heat released during combustion or nuclear reaction (TWh/yr) waswhereis the thermal efficiency based on the lower heating value of the fuel, since combusted water is assumed to be emitted as a vapor. The thermal efficiency (Table 3 ) is thus the electric energy produced divided by the sum of electric energy plus waste heat aside from that used to evaporate water produced during combustion. In the equation,= 2465 J/g is the latent heat of evaporation,= 1779 g‐HO/kWh electricity (0.47 gal‐HO/kWh) is the average cooling water evaporated per kWh of electric power generated from U.S. thermal power plants [.,], and= 3.6 × 10J/kWh. Cooling water is additional water from a nearby lake or river or the ocean used to absorb some of the heat released that is not already used to evaporate water produced during combustion. For open biomass burning, biofuel burning, waste burning, fossil fuel fires, oil and gas production, and industrial processes and product use, the same equation was used, except that water cooling did not occur (the last term was removed) so that all heat of combustion went directly to the air. For open biomass burning, the AMFs and AHFs included both anthropogenic and natural components, of which ~7% was natural (section 2). E heatc ) or converted to mechanical energy that was eventually dissipated to heat through friction (E heatf ). The combustion and friction heat AHFs (TWh/yr) for each transportation sector (onroad vehicles, rail, nonroad vehicles, shipping, aircraft, etc.) were (3) F is the quantity of fuel used (kg/yr), H L is the lower heating value of a fuel (J/kg fuel, Table η tw is the tank‐to‐wheel efficiency of the fuel (Table f C (g‐C/kg fuel, Table (4) E CO2 here is the CO 2 emissions (Tg‐CO 2 /yr, Table m C = 12.011 g/mol and m CO2 = 44.0098 g/mol are the molecular weights of carbon and carbon dioxide, respectively. For aircraft, the fuel use here was already known in each model grid cell and each time step. For transportation, no water cooling occurred either. In that case, all energy in the fuel was dissipated either as waste heat of combustion () or converted to mechanical energy that was eventually dissipated to heat through friction (). The combustion and friction heat AHFs (TWh/yr) for each transportation sector (onroad vehicles, rail, nonroad vehicles, shipping, aircraft, etc.) wererespectively, whereis the quantity of fuel used (kg/yr),is the lower heating value of a fuel (J/kg fuel, Table 3 ), andis the tank‐to‐wheel efficiency of the fuel (Table 3 ). When the fuel use was not known spatially, it was calculated from the carbon content of the fuel(g‐C/kg fuel, Table 3 ), which was known spatially fromwherehere is the COemissions (Tg‐CO/yr, Table 2 ) from the transportation sector, and= 12.011 g/mol and= 44.0098 g/mol are the molecular weights of carbon and carbon dioxide, respectively. For aircraft, the fuel use here was already known in each model grid cell and each time step. 2 O/yr) from combustion were then calculated as (5) M is the mole ratio of H 2 O:CO 2 from combustion (Table m H2O is the molecular weight of water vapor (g/mol). For onroad transportation, the M used was for a gasoline/diesel mixture. For nonroad and shipping, it was from diesel. For waste, biofuel, and open biomass burning, it was from biomass. For other sectors, it was assumed to be from coal. Water vapor emissions (Tg‐HO/yr) from combustion were then calculated aswhereis the mole ratio of HO:COfrom combustion (Table 3 ) andis the molecular weight of water vapor (g/mol). For onroad transportation, theused was for a gasoline/diesel mixture. For nonroad and shipping, it was from diesel. For waste, biofuel, and open biomass burning, it was from biomass. For other sectors, it was assumed to be from coal. 2 O/yr) from evaporation of cooling water in thermal power plants (coal, gas, oil, biomass, and nuclear) were calculated as (6) Finally, water vapor emissions (Tg‐HO/yr) from evaporation of cooling water in thermal power plants (coal, gas, oil, biomass, and nuclear) were calculated as Another anthropogenic source of carbon dioxide, water vapor, and heat is human respiration. However, because these gases and heat are removed from the air to produce organic material during photosynthesis, and people eat food containing organic material in order to respire, these emissions do not represent a long‐term net addition to the atmosphere, only a displacement from agricultural fields to populated regions. As such, the emissions were not included in the model calculations, but their magnitudes are provided in Table 2 for comparison and completeness. Human respiration occurs by the generalized metabolic reaction C 6 H 12 O 6 + O 2 →6CO 2 + 6H 2 O + heat (2880 kJ/mol‐C 6 H 12 O 6 ). A typical person breathes out an average of 1 kg‐CO 2 /day. To ensure mass and energy consistency, such emissions are accompanied by 0.409 kg‐H 2 O/day and 126.2 W of heat release. Average emissions account for various activities during the day. For example, sleeping releases ~70 W; resting releases ~100 W; walking releases ~200 W; and more strenuous activity releases ~300 W [Sailor, 2011]. Table 2 indicates that the global emissions of CO 2 , H 2 O, and heat from human metabolism in 2005 were ~2360 Tg‐CO 2 /yr, 965 Tg‐H 2 O/yr, and 7150 TWh/yr (~0.0016 W m−2), respectively. These numbers represent 5.0%, 1.7%, and 5.6% of the total anthropogenic emissions of CO 2 , H 2 O, and heat, respectively, from the table. Flanner [2009] calculated the global AHF in 2005 excluding open biomass burning and metabolism and assuming no water cooling of power plants, as 0.028 W m−2 (125,000 TWh/yr). The equivalent number here for 2005 is also 0.028 W m−2 (125,000 TWh/yr). Of this, 0.023 W m−2 (99,800 TWh/yr) is from heat due to combustion, electric power use, and frictional dissipation of motive force, and 0.0057 W m−2 (25,600 TWh/yr) is from combustion heat used to evaporate power plant cooling water (not water from combustion, which is accounted for in water vapor direct emissions here) (Table 2). Open biomass burning (anthropogenic plus natural) results in an additional AHF here of 0.005 W m−2 (21,000 TWh/yr), giving a total AHF before evaporating cooling water of 0.033 W m−2 (146,000 TWh/yr) and after evaporating of 0.027 W m−2 (121,000 TWh/yr). Of the total AHF after evaporating cooling water (121,000 TWh/yr), 76% was due to direct combustion and the rest to electric power use or frictional dissipation of motive energy (Table 2). Of the same total, 17.2% was due to open biomass burning combustion, 20.0% was due to biofuel combustion, 20.4% was due to transportation (onroad, nonroad, shipping, rail, and air), 29.5% was due to public plus private electric power generation and use of the electric power, and 10.4% was due to industrial processes and product use. Of the public plus private electric power generation combustion heat, about 41.2% was due to coal, 26.0% was due to natural gas, 12.8% was due to nuclear power, and 10.4% was due to fuel oil. The major anthropogenic emission source of water vapor considered was evaporation during cooling of coal, nuclear, natural gas, and biofuel power plants and industrial facilities. Globally, about 25,600 TWh/yr (E elec L e V c /A j ) of combustion heat was used to evaporate cooling water in 2005. The second largest anthropogenic source of water vapor considered was fossil fuel, biofuel, and anthropogenic plus natural biomass burning combustion. Such combustion releases carbon dioxide, water vapor, and air pollutants. Of the evaporation plus combustion sources of water vapor, which totaled ~56,680 Tg‐H 2 O/yr worldwide in 2005, about 65.8% was from evaporation. However, to put this in perspective, the total anthropogenic emission rate of water vapor was only about 1/8800th the worldwide natural emission rate of ~500 million Tg‐H 2 O/yr. An additional source of anthropogenic water vapor is evaporation of water pumped from underground aquifers for irrigation and human consumption. When such water evaporates, it cools the ground; however, the heat is regained in the atmosphere when the water vapor recondenses. Because these two effects cancel exactly in the global average and since water vapor is a greenhouse gas, the net impact of aquifer water vapor should be to cause interim climate warming. This AMF was not treated here as discussed previously. The direct radiative forcing at the top of the atmosphere (TOA) due to nonpumped‐aquifer anthropogenic H 2 O can be estimated simplistically by calculating the direct forcing from all natural plus anthropogenic water vapor in a global model (+79.4 W m−2, Table 1) and multiplying this by the ratio of the anthropogenic to natural emission rate of H 2 O, 1/8800, giving +0.009 W m−2. Possibly by coincidence, this TOA forcing is similar to the surface energy flux of water vapor, 0.0087 W m−2, calculated by multiplying the total anthropogenic H 2 O emissions in Table 2 by L e and dividing by the number of seconds per year and the area of the Earth. Although these numbers are similar, one is calculated at the TOA and the other at the surface. However, most water vapor, thus its radiative forcing, is confined to the boundary layer. Adding the direct heat emitted to the air after evaporating water of 0.027 W m−2 gives a total all‐source AHF plus AMF direct forcing of 0.036 W m−2. For comparison, the direct forcing due to anthropogenic CO 2 is 1.82 W/m2 and that due to CH 4 is 0.48 W/m2 [IPCC, 2014, Table 8.2]. Table 2 indicates that, of the combustion plus evaporation water vapor, 34.4% was from coal power generation, 29.8% was from natural gas power generation, 8.6% was from fuel oil power generation, 8.3% was from nuclear power, 5.3% was from solid biofuel burning (not for electric power), 4.7% was from open biomass burning, and 4.3% was for all transportation. Based on its global atmospheric loading and the rainfall rate, water vapor has an e‐folding lifetime in the atmosphere against removal of about 9 days. It heats near‐surface air in three ways: by trapping thermal‐IR radiation as a greenhouse gas, by releasing latent heat upon condensation, and by raising geopotential heights upon increasing surface air pressure.

3 Description of the 3‐D Simulations The 3‐D global model used for this study was GATOR‐GCMOM [Jacobson, 2002a, 2003, 2010a, 2012; Jacobson et al., 2007]. The model simulates time‐dependent feedbacks among meteorology, gases, aerosol particles, hydrometeor particles, the surface ocean and deep ocean, vegetated and bare soil, freshwater surfaces, and radiation transfer. The version of the model used is very close to that used and evaluated in Jacobson [2012], so most of the description is not repeated here. Three discrete (multiple size bin) aerosol size distributions and three discrete hydrometeor distributions were treated (Table 4). The three aerosol distributions were an emitted fossil fuel soot (EFFS) distribution, an emitted combined biofuel soot and biomass burning soot (BFBB) distribution, and an ultimate internally mixed (UIM) distribution, which is the most internally mixed of all distributions. Each aerosol distribution contained 14 size bins ranging from 2 nm to 50 µm in diameter. The three hydrometeor distributions included liquid, ice, and graupel, each with 30 size bins (0.5 µm to 8 mm in diameter). Table 4. Aerosol and Hydrometeor Discrete Size Distributions Treated in the Model and the Parameters (Number Concentration and Chemical Mole Concentrations) Present in Each Size Bin of Each Size Distribution Species Aerosol Emitted Fossil Fuel Soot (EFFS) Aerosol Emitted Biofuel and Biomass Burning Soot (BFBB) Aerosol Ultimate Internally Mixed (UIM) Cloud/Precipitation Liquid Cloud/Precipitation Ice Cloud/Precipitation Graupel Number X X X X X X BC X X X X X X POM X X X X X X SOM X X X X X X Tar balls X X X X X H 2 O(aq)‐h X X X X X X H 2 SO 4 (aq) X X X X X X HSO 4 − X X X X X X SO 4 2− X X X X X X NO 3 − X X X X X X Cl− X X X X X X H+ X X X X X X NH 4 + X X X X X X Na+(K+, Mg2+, Ca2+) X X X X X NH 4 NO 3 (s) X X X X X (NH4) 2 SO 4 (s) X X X X X Soil dust X X X X Poll/spores/bacteria X X X X NaCl(s) X X X X H 2 O(aq)‐c X H 2 O(s) X X Each size bin of the EFFS distribution contained black carbon (BC), weakly to moderately absorbing primary organic matter (POM), secondary organic matter (SOM), hydrated liquid water, H 2 SO 4 (aq), HSO 4 −, SO 4 2−, NO 3 −, Cl−, H+, NH 4 +, NH 4 NO 3 (s), and (NH 4 ) 2 SO 4 (s). Each size bin of the BFBB distribution contained these same components plus tar balls, Na+, K+, Ca2+, and Mg2+. Each size bin of the UIM distribution contained the same components as the EFFS and BFBB distributions plus soil dust, pollen, spores, and bacteria. Each size bin of each hydrometeor distribution contained the same components as in all three aerosol distributions plus condensed as well as hydrated liquid water and/or ice. Gases, such as HNO 3 , HCl, NH 3 , H 2 SO 4 , and organics could condense onto or dissolve into EFFS, BFBB, and UIM aerosol particles and dissolve within liquid hydrometeor particles or chemically react on ice and graupel hydrometeor particle surfaces. Global simulations were run at 4°SN × 5°WE horizontal resolution and with 68 vertical sigma‐pressure atmospheric layers from the ground to 0.219 hPa (≈60 km), including 15 layers from 0 to 1 km and 500 m resolution from 1 to 21 km. The center of the lowest atmospheric layer was 15 m above ground. The model also treated 10 ocean layers. A 2‐D time‐dependent mixed‐layer ocean dynamics module driven by surface wind stress and conserving potential enstrophy, vorticity, energy, and mass was used to solve for mixed‐layer velocities, heights, and horizontal energy transport [Ketefian and Jacobson, 2009]. The model also solved air‐ocean exchange of CO 2 and all gases, vertical chemical and energy diffusion through all ocean layers, and 3‐D ocean equilibrium chemistry and pH. The model was initialized on 1 January 2006 with 1° × 1° reanalysis meteorological fields [Global Forecast System, 2013] and run forward for 20 years with no data assimilation. Since this period is too short for long‐lived greenhouse gas mixing ratios and deep ocean temperatures to equilibrate with atmospheric emissions, the simulation results here are transient in nature. However, the period is sufficiently long for aerosol emissions to be in rough equilibrium with their atmospheric burdens. Atmospheric composition of all gases, aerosol particles, and hydrometeor particles changed over time. Natural emissions of soil dust, sea spray, pollen, spores, bacteria, biogenic gases, lightning gases, and volcanic gases and particles varied with time and space primarily as a function of meteorology, as described and quantified in Jacobson and Streets [2009]. Table 2 provides heat, moisture, and carbon dioxide fluxes from open biomass burning, biofuel burning, and fossil fuel burning. Table 5 provides contemporary (circa 2005) gas and particle emission rates from these sources. References for and methods of calculating emission rates are described mostly in the footnote to Table 5. BB emission rates were determined from 5 years of gridded 8 day cycle fuel‐burn data (Table 5) repeated every 5 years. Emissions from fossil fuels and solid biofuels were repeated each year. Although ~7% of BB is natural, the emission rates in Table 5 and used in the simulations are for all BB. In addition, tar ball (TB) emission rates from biomass and biofuel burning were assumed to be 30% of total POM (including tar ball) emission rates by mass. Adachi and Buseck [2011] found that 1–14% of biomass burning particles by number were tar balls within 4 min of burns. The fact that most tar balls are in the accumulation mode (mean diameter 210 ± 100 nm), whereas other organics and BC are generally smaller, justifies a larger mass fraction than number fraction of tar balls. Table 5. World Fine‐Particle (<1 µm Diameter) Global Emission Rates of BC, POC, and Other Particle Components in Soot From Fossil Fuels, Solid Biofuels, and Open Biomass Burning Fossil Fuel Biofuel Biomass Burning Total (a + b + c) Loading (mg/m2) Lifetime (days) BC 3.29 3.33 2.70 9.32 0.197 3.93 POC‐C 4.64 9.40 12.1 26.1 0.874 Tar ball‐C 4.03 5.20 9.23 0.192 3.87 S(VI) 0.056 0.51 0.69 1.25 3.63 Na+ ‐‐ 0.031 0.014 0.045 ‐‐ ‐‐ K+ as Na+ ‐‐ 0.22 0.87 1.09 ‐‐ ‐‐ Ca2+ as Na+ ‐‐ 0.22 0.16 0.38 ‐‐ ‐‐ Mg2+ as Na+ ‐‐ 0.096 0.09 0.19 ‐‐ ‐‐ NH 4 + ‐‐ 0.006 0.027 0.032 ‐‐ ‐‐ NO 3 ‐ ‐‐ 0.27 0.23 0.50 ‐‐ ‐‐ Cl− ‐‐ 0.51 0.38 0.89 ‐‐ ‐‐ One simulation discussed required 1850 emissions. Natural emissions were as just described. The 1850 biomass burning emission rates were assumed to be ~34% those of 2005 emissions and were mostly natural (see section 2). The 1850 anthropogenic emissions were obtained from Clarke et al. [2007]. A baseline simulation (Simulation A), several sensitivity simulations, and four random variation simulations were run to isolate several individual contributions to global warming. The baseline simulation included emissions of all circa 2005 anthropogenic and natural gas emissions, particle emissions, heat fluxes, and moisture fluxes and 2005 urban surfaces (USs). The sensitivity simulations were the same as (Simulation A) except as follows: (Simulation B) 1850 emissions of the following together: (i) all BB (including its heat and moisture fluxes), (ii) anthropogenic aerosol particles and their precursor gases (AAPPGs), (iii) anthropogenic greenhouse gases (AGHGs), (iv) AHFs, and (v) AMFs, but no USs; (Simulation C) no present‐day anthropogenic or natural biomass burning (BB) (or associated gases, particles, AHFs, or AMFs); (Simulation D) no AHF from any anthropogenic source; (Simulation E) no AMF from any anthropogenic source; (Simulation F) no AHF or AMF from power plants or electricity use; (Simulation G) no AHF or AMF from transportation; (Simulation H) no AHF from present‐day biomass burning; (Simulation I) no CAE I (absorption by black and brown carbon inclusions within cloud drops or ice crystals); and (Simulation J) no CAE I or present‐day BB gases, particles, AHFs, or AMFs. Four additional 20 year baseline simulations were run, the same as Simulation A except that each had a different random perturbation in the sixth digit of one variable in the initial conditions. The purpose of these simulations was to determine whether differences found in the sensitivity test simulations here were statistically significant, based on a two‐sided t test, relative to deterministic chaotic variation that occurs in climate models [e.g., Chervin and Schneider, 1976]. AAPPGs consisted of emitted fossil fuel soot (EFFS), biofuel soot and gases (BSG), biomass burning (BB) particles and gases, all other anthropogenic NO x , SO x , NH 3 , CO, H 2 , and speciated volatile organic carbon compounds. EFFS consisted of fossil fuel particulate black carbon (BC), primary organic matter (POM), H 2 SO 4 (aq), HSO 4 −, and SO 4 2. Other chemicals could enter EFFS particles, as described in Table 4. Biofuel soot and gases (BSG) and biomass burning (BB) were assumed to consist of the same plus tar balls, H 2 O(aq), Na+, K+, Mg2+, Ca2+, Cl−, NO 3 −, and many gases (listed in the footnote to Table 5), including CO 2 . AGHGs include anthropogenic CO 2 , H 2 O, CH 4 , N 2 O, CO, and chlorofluorocarbons. Urban surfaces (US) were assumed to have the same characteristics as in Jacobson and Ten Hoeve [2012].

4 3‐D Climate Response Simulation Results Figure 1 compares modeled globally distributed aerosol optical depth and cloud fraction with Moderate Resolution Imaging Spectroradiometer (MODIS) satellite data. Detailed comparisons of modeled aerosol absorption optical depth with satellite and other data are given in Jacobson [2012] so are not repeated here. Modeled total 550 nm aerosol optical depth is slightly lower than the observed global value but not significantly. The spatial distribution and magnitude of modeled cloud fraction are roughly similar to the data although the model is biased high in some biomass burning regions of Brazil and Africa and low in some biofuel burning regions of central Asia. These biases may reflect uncertainties in emissions in these regions. Figure 1 Open in figure viewer PowerPoint (a) MODIS versus (b) modeled 550 nm clear‐sky annually averaged aerosol optical depth. (c) MODIS versus (d) modeled annually averaged cloud fraction. 4.1 Effects of Biomass Burning Table 6 and Figure 2 summarize the globally, land‐ and ocean‐averaged surface temperature changes between the baseline simulation and several of the sensitivity simulations, averaged over 20 years. This section focuses on results from the 20 year simulations that isolated the transient effects of BB particles, gases, heat, and moisture on global climate and atmospheric composition (Simulation A versus C). Averaged over this period, BB was modeled to warm the surface in the global average by ~0.4 K (Figure 2a). The warming would likely increase over a longer simulation period, as suggested by the slight slope of temperature change in Figure 2a, since CO 2 mixing ratios and deep ocean temperatures are not in equilibrium with CO 2 emissions after 20 years. The simulation‐averaged top‐of‐the‐atmosphere (TOA) net downward minus upward solar plus thermal‐infrared irradiance change resulting in this surface temperature change averaged over this period was +0.9 W/m2. This irradiance change is transient and would go to zero if the climate system had reached equilibrium. Note also that this irradiance change is not a direct or indirect forcing but a climate response computed simultaneously with the temperature change and accounting for all feedbacks, including to the temperature of the surface and deep oceans. Table 6. The 20 Year and Globally Averaged Differences in Near‐Surface (at 15 m) Air Temperature (T), Cloud Optical Depth (COD), Cloud Fraction (CF), Precipitation (Pr), and Near‐surface Wind Speed (WS) Between the Baseline Simulation and Each of Several Sensitivity Simulations Simulation (Effect on Climate From) ΔT (K) ΔCOD (‐‐) ΔCF (‐‐) ΔPr (mm/d) ΔWS (m/s) Stat. sig. ? AGHG, AAPPG, AHF, AMF, US since 1850 +1.00 +1.55 −0.018 −0.052 −0.072 Yes All BB (and AHF, AMF from BB) +0.409 −0.020 −0.007 −0.024 −0.039 Yes All AHF from all anthropogenic sources +0.033 −0.026 −0.0007 −0.0009 −0.008 Yes All AMF from all anthropogenic sources −0.020 +0.009 −0.00002 +0.001 +0.0007 No AHF + AMF from power plants and electricity use +0.034 −0.027 −0.0004 +0.001 −0.008 Yes AHF + AMF from transportation −0.010 −0.029 −0.00002 +0.004 −0.002 No Figure 2 Open in figure viewer PowerPoint National Aeronautics and Space Administration‐Goddard Institute for Space Studies (NASA‐GISS), 2014 (a) Cumulative‐average (thick line) and current (in 5 day increments, thin line) near‐surface air temperature change over 20 years between Simulation A (with BB) and Simulation C (no BB). Any point on the cumulative line is the average value between the start and that time. (b) Modeled temperature differences due to all present‐day biomass burning (BB), (c) all present‐day anthropogenic heat fluxes (AHFs), and (d) present‐day AHFs and AMFs from power plants plus electricity use (PP + E AHF + AMF). (e) Modeled 15 m air temperature differences from 20 year equilibrium simulations using present‐day anthropogenic greenhouse gas (AGHG) emissions, anthropogenic aerosol particle and their precursor gas (AAPPG) emissions, AHFs, AMFs, and urban surfaces (US) minus those from 1850 emissions and surfaces. (f) Net observed temperature difference 2010–2013 relative to 1910–1930 mean [),]. Numbers in parentheses are globally averaged changes. In most cases land (l) and ocean (o) changes are given as well. Two sets of global, 3‐D sensitivity simulations were run to isolate the physical processes causing this global warming. One test (Simulation A versus H) indicated that about 7% of the warming from biomass burning was due to the heat emitted during burning (AHF) and the resulting impacts of the heat on cloud reduction and local temperatures. A second test (Simulation J versus I compared with Simulation A versus C) indicated that about 32% of the warming due to biomass burning was due to absorption by black and brown carbon inclusions within cloud drops and ice crystals (CAE I), which causes clouds to burn off. Both factors were missing from previous assessments of biomass burning impacts on climate, as were AMFs and CAE II (section 1). In the base BB case (Simulation A versus C), all biomass burning caused ~0.4 K warming in the global average, ~0.6 K warming over land, and up to ~1 K warming over the Arctic (Figure 2b). Figure 3 shows differences in other parameters due to BB. Figure 4 shows that temperature changes between Simulations A and C were strongly statistically significant in the Arctic and over most of the globe. Arctic warming resulted in a decrease in the average global sea ice depth by ~3.1% and snow depth by ~0.6%, decreasing global surface albedo by ~1.3%. Some of the Arctic warming was due to an increase in absorbing aerosol particles over the Arctic and within Arctic snow and sea ice. Some warming was also due to an increase in the poleward transport of heat and air caused by biomass burning. The increase in air mass, as reflected by higher surface pressures over the Arctic (Figure 3g), increased subsidence and compressional heating of the Arctic surface. Figure 3 Open in figure viewer PowerPoint (a) Modeled differences in 550 nm aerosol optical depth (AOD), (b) 550 nm aerosol absorption optical depth (AAOD), (c) 650 nm cloud optical depth (COD), (d) 650 nm cloud absorption optical depth (CAOD), (e) cloud fraction, (f) precipitation, and (g) sea‐level pressure between Simulations A and C, with minus without all present‐day biomass burning. Numbers in parentheses are globally averaged changes, land (l) changes and ocean (o) changes, respectively. Figure 4 Open in figure viewer PowerPoint t test relative to four random perturbation simulations. The numbers in parentheses are domain‐averaged confidence levels. (a–c) Spatially‐distributed confidence levels for several of the temperature difference plots in Figures 2 b, 2 c, and 2 d, respectively, as determined from a two‐sidedtest relative to four random perturbation simulations. The numbers in parentheses are domain‐averaged confidence levels. BB increased gas and particle emissions, thereby increasing 550 nm aerosol optical depth (AOD) over land to 4 times the AOD over the ocean (Figure 3a). The increase in 550 nm aerosol absorption optical depth (AAOD) (Figure 3b) was about one quarter the increase in total 550 nm AOD due primarily to both BC and tar ball absorption. Percentage wise, BB increased globally averaged AAOD by ~55% and AOD by ~15%. BB reduced globally averaged cloud optical depth by ~0.2%, with 10 times the reduction over land than over the ocean (Figure 3c). Part of this increase was due to an increase in cloud absorption optical depth (CAOD) by ~16% globally (Figure 3d), although the locations of the highest increase in CAOD were not the same as the locations of the greatest reduction in COD. BB also decreased globally averaged cloud fraction by ~1% (Figure 3e), with a ~2% decrease in cumulus cloud fraction and a ~0.13% decrease in stratiform cloud fraction. BB also decreased globally averaged precipitation by ~0.8% (Figure 3f). The greater number of aerosol particles due to biomass burning increased the globally averaged number concentration of cloud condensation nuclei worldwide by ~7.3% and the number concentration of ice deposition nuclei by ~4%, increasing the number concentration of liquid drops and ice crystals in clouds that were not destroyed by absorptive heating. In addition, the water vapor feedback due to the rise in surface temperatures cause by BB and the water vapor released during BB combustion both enhanced atmospheric column water vapor (kg/m2) by 3.6%, increasing column cloud liquid water and ice by ~2.2% and ~0.6%, respectively, increasing the size of liquid drops and ice crystals. Both factors increased the number of collisions and bounceoffs among crystals and between crystals and liquid drops, increasing the global lightning flash rate by ~6%, contributing to enhanced nitric oxide in the troposphere. The methodology for calculating flash rate is described in Jacobson and Streets [2009], who also compare modeled baseline flash rate with satellite data globally. Increases in lightning flashes due to biomass burning have also been modeled at high resolution [Ten Hoeve et al., 2012]. BB increased globally averaged near‐surface ozone by ~0.6 ppbv globally and ~1.1 ppbv over land. It increased near‐surface dry particulate matter smaller than 2.5 microns in diameter (PM 2 .5 ) by ~1.9 µg/m3 globally and ~2.4 µg/m3 over land. The population‐weighted increases in ozone and dry PM 2.5 were ~0.9 ppbv and ~2.9 µg/m3, respectively. A population‐weighted value is the product of a grid cell value for a parameter (e.g., ozone) and the population of the grid cell, summed over all grid cells and divided by the world population. Excess premature mortalities from ozone and PM 2.5 were calculated as in Jacobson [2010b] by combining modeled concentrations in each time step and in each surface grid cell with population in the cell and relative risks of premature mortality for short‐term exposure to ozone (compiled by Ostro et al. [2006]) and long‐term exposure to PM 2.5 [Pope et al., 2002] above a threshold (35 ppbv for ozone and 0 µg/m3 for PM 2.5 , but with one fourth the relative risk below 8 µg/m3 than above). With this method, BB was calculated to cause ~20,000 (10,000–30,000) premature mortalities/yr from excess ozone and ~230,000 (63,000–405,000) premature mortalities/yr from excess PM 2.5 for a total of ~250,000 (73,000–435,000) excess premature mortalities/yr, representing 5–10% of global mortalities due to outdoor plus indoor air pollution. This estimate is very rough since the model did not resolve cities spatially due to large model grid cells, and relative risk data are also uncertain. Table 2 indicates that ~17.8% of the AHF from all sources, ~17.3% of all anthropogenic CO 2 , and ~4.7% of all AMF were from anthropogenic plus natural (~7% of the total) open biomass burning. The AHF causes direct warming of the air. While nonpermanently deforested CO 2 helps to regrow foliage, it does so with a time lag, resulting in a net accumulation of CO 2 and heating of the atmosphere [Jacobson, 2004]. Absorption due to BC and BrC inclusions alone within cloud drops and ice crystals (CAE I) caused ~32% of BB warming here. 4.2 Global Effects of All AHF and AMF Globally, the near‐surface temperature changes due to all AHFs (+0.03 K) (Figure 2c, Table 6) and power plant plus electricity use AHFs + AMFs (+0.03 K) (Figure 2d, Table 6) were statistically significant over many land areas (Figures 4b and 4c, respectively) although less significant than in the BB case. However, modeled temperature changes due to all AMFs and AHFs + AMFs from transportation were not statistically significant (Table 6). As such, conclusions are not drawn here about those simulations. The reason for the lower statistical significance of temperature changes due to AMF is the offsetting impacts of the warming due to the thermal‐IR absorption and the cooling arising from the high solar albedo associated with low‐cloud formation due to water vapor. These offsetting factors caused net changes in temperature due to AMFs to be small through much of the world. AHFs, on the other hand, increase saturation vapor pressures in the global average, slightly reducing cloudiness (Table 6). In addition, transportation emits 22% of the CO 2 , 4% of the H 2 O, and 33% of the heat that the power sector emits (Table 2). Thus, the transportation sector resulted in lesser temperature changes than did the power sector. 4.3 Effects of All Anthropogenic Gases and Particles on Global Climate The difference between Simulations A and B provides a model estimate of net anthropogenic global warming due to emission and land‐use differences between 1850 and the present. Differences between those simulations include emissions from anthropogenic greenhouse gases and anthropogenic aerosol particles (BB, EFFS, BSG, and their precursor gases, and particle sulfate, nitrate, ammonium, and organic carbon from other anthropogenic sources). They also include AHFs, AMFs, and urban surfaces. The globally averaged near‐surface air temperature difference was ~1 K (Table 6, Figure 2e). This is slightly higher than the data analysis estimate of 0.92 K for 2010–2013 relative to 1900–1920 in Figure 2f and somewhat higher than the data analysis estimate of 0.78 ± 0.06 K for 2003–2012 relative to 1850–1900 [IPCC, 2014]. Figures 2e and 2f compare the spatial distribution of modeled global warming since 1850 due to all anthropogenic gas, aerosol, heat, and moisture fluxes and urban surfaces with that of observed warming since 1900–1920. The model was able to replicate the general major locations of observed warming although differences can be seen, such as a high bias in model predictions over the northeastern Pacific Ocean and a slight high bias over the Arctic. The net observed warming consists of warming due to anthropogenic greenhouse gases and absorbing aerosol particles and the urban heat island effect less cooling due to nonabsorbing aerosol particles. Major contributors to the warming were greenhouse gases, biomass burning (gases, particles, AHFs, and AMFs) (~0.4 K, Figure 2e, Table 6), emitted fossil fuel soot, and biofuel soot and gases (including AHF + AMF from it). Cooling aerosol particles and their precursor gases, including SO 2 , NO x , and weakly absorbing organic gases, offset a substantial portion of the warming. However, it should be pointed out that individual warming and cooling contributions to global warming are not usually linearly additive. In other words, if chemical A causes Y K of warming and chemical B causes Z K of warming, warming due to A plus B does not equal Y + Z K since some of warming due to A and B overlap (is the same warming), and nonlinear feedbacks affect A + B differently from A or B individually. Finally, IPCC [2014] includes a compilation of results from models of different levels of complexity attempting to simulate anthropogenic climate change. Since climate response models used in those calculations did not include full treatments of CAEs I or II, AHFs, or AMFs, their accuracy in predicting the observed mean global temperature change without such warming may be because they simultaneously underpredicted cooling due to aerosol particles and/or their cloud‐induced changes. Given the uncertainty in modeling indirect effects, this is a distinct possibility. Similarly, results here are subject to uncertainties in emissions, grid resolution, indirect effect treatment, and other factors. 4.4 Direct Radiative Forcing of Black Carbon In Jacobson [2010a], global climate response calculations were run accounting for indirect effects, semidirect effects, CAE I (but not CAE II), and AHFs and AMFs of biofuel soot and fossil fuel soot, but not of BB. EFFS alone was found to increase global temperatures by 0.3–0.5 K over 15 years. Bond et al. [2013] report +0.28 K and +0.44 K temperature increases from two other studies that examined the climate impacts of fossil fuel BC. In the latter case, BC was quadrupled and the temperature response divided by four. Here the modeled solar plus thermal‐IR direct radiative forcing (DRF) of BC alone (not BC + POM) from all sources (preindustrial plus current fossil fuels, solid biofuels, and biomass burning) in the clear skies and interstitially between cloud drops at the RH of the cloud was calculated as +0.52 W/m2 (94.3% solar). DRF was calculated from one simulation with all anthropogenic and natural emissions present, but with two radiation calls each time step during the simulation. The first radiation called included all BC from all sources; the second included no BC from any source. The net downward minus upward solar plus thermal‐infrared irradiance change at the TOA between the two radiation calls was averaged over the simulation to obtain the DRF. This DRF, which was calculated with three discrete aerosol size distributions and over a time‐dependent simulation, is similar to the +0.55 W/m2 determined from 18 size distributions calculated for a short‐time period [Jacobson, 2001]. It also compares with the estimated DRF of +0.71 (+0.08 to +1.27) W/m2 due to industrial era BC + BrC and +0.88 (+0.17 to +1.48) W/m2 due to all BC + BrC from Bond et al. [2013]. The solar plus thermal‐IR DRF of BC inclusions in cloud drops and ice crystals (CAE I) was calculated here as an additional 0.062 W/m2 (99.7% solar), which is in the range of the 0.05–0.07 W/m2 estimate from Bond et al. [2013]. It is also close to the estimate of +0.07 W/m2 by Chuang et al. [2002] but higher than that (+0.01 W/m2) of Ghan et al. [2012], possibly due to the much lower BC loading and different method of treating absorption in the latter study. As shown next, though, the change in climate response of absorbing aerosols within a cloud is inversely proportional to the change in radiative forcing in the cloud drops, so the cloud drop radiative forcing is not an indicator of the warming due to BC in a cloud.

5 Description of 1‐D Simulations Vertical 1‐D simulations were run with GATOR‐GCMOM, collapsed from 3‐D (references in section 3) to 1‐D, to examine the effects of CAE I and II on the burn‐off of individual clouds and the feedback of such burn‐off to ground surface temperature. One purpose of these simulations was to find example conditions leading to cloud burn‐off. Another was to determine whether surface temperature changes are proportional to direct radiative forcing changes in clouds containing absorbing aerosol particles. Time‐dependent processes treated in the 1‐D version of the model included aerosol particle microphysics and chemistry, cloud hydrometeor particle microphysics and chemistry, gas chemistry, radiative transfer through aerosol particles, clouds, gases, and subsoil processes. The model treated 89 layers between the surface and 60 km. Each layer in the vertical column was initialized with gases, meteorological parameters, and a 60 bin discrete size‐ and composition‐resolved aerosol distribution. Only the UIM aerosol distribution was treated along with the three hydrometeor distributions in Table 4, so all aerosol components were initially in the UIM distribution. Aerosol processes treated included homogeneous nucleation, condensation, dissolution, internal‐particle chemical equilibrium, aerosol hydration of liquid water as a function of RH (up to 100%), temperature, composition, ion crystallization within aerosol particles, aerosol‐aerosol coagulation, aerosol‐hydrometeor coagulation, sedimentation, and dry deposition [Jacobson, 2002b, 2003; Jacobson et al., 2007]. Figure 5a shows the initial vertical profile of RH and ambient black carbon. No emissions were considered during the short simulation period, but particle and gas sedimentation, dry deposition, rainout, and washout were. The 1‐D model did not consider horizontal advection or large‐scale pressure changes and did not treat convection (thus vertical velocity was zero but it did treat vertical diffusion due to buoyancy). It assumed a relatively stratiform cloud. Initially, no cloud existed in the column. Because some layers were slightly supersaturated initially, a cloud subsequently began to form. During cloud formation, the RH decreased not only due to the reduction in water vapor and latent heat release during condensation but also due to the radiative heating of air by absorbing inclusions and by water vapor in the cloud. Temperatures (thus RHs) in the cloud were also affected by radiative heating and/or cooling due to hydrometeor, gas, and aerosol scattering and absorption within the cloud. Figure 5 Open in figure viewer PowerPoint (a) Initial RH and aerosol BC mass concentration (summed over all sizes) for 1‐D simulations of the effects of absorbing aerosol on cloud burn‐off. (b) Time‐dependent BC mass concentration in aerosol and inclusions in hydrometeor particles during Sim1, where absorption by BC aerosol (interstitially between cloud drops and in cloud‐free air) and inclusions in cloud drops were accounted for. The increase in BC concentration in aerosol after cloud evaporation was due to the release of BC cores from evaporated cloud drops as aerosol. BC was also lost over time due to wet and dry removal. At 4500 s, only one cloud layer remained. (c) Time‐dependent change in the vertical profile of cloud liquid content (LWC) during Sim1. (d) Same as Figure 5 c but for Sim2, where BC absorption within cloud drop inclusions was not calculated, but BC inclusions were still treated and absorption by BC aerosol interstitially between cloud drops and in cloud‐free air was still treated. (e) Time‐dependent difference in the vertical temperature profile between Sim1 and Sim2. This is the temperature change due to CAE I. (f) Same as Figure 5 e but for RH. (g) Cloud absorption optical depths due to BC inclusions alone during Sim1. (h) Total cloud optical depths during Sim1 (“with absorption by inclusion, interstitial, and outside‐cloud BC”), Sim2 (“with absorption by interstitial and outside‐cloud BC only”), and Sim3 (“without BC at all”). (i) Time‐dependent net down minus up TOA irradiance change between Sim1 and Sim2 (“BC inclusions”) and between Sim1 and Sim3 (“inclusion, interstitial, and outside‐cloud BC”). (j) Same as Figure 5 i but for time‐dependent surface ground temperature change. This latter case accounts for indirect, semidirect, and CAE I and II. (k) Relationship between the total (inclusion, interstitial, and cloud‐free BC) direct radiative forcing in Figure 5 i and the ground temperature change in Figure 5 j that is due to all BC in cloud inclusions + aerosol. Hydrometeor particles were size and composition resolved and formed on top of and coagulated with aerosol particles. Three hydrometeor size distributions could exist in each layer—liquid, ice, and graupel. Since hydrometeor particles formed on top of aerosol particles and the aerosol inclusions were tracked in hydrometeor particles of each size, each hydrometeor particle contained multiple chemical components (Table 4). Aerosol inclusions within each hydrometeor particle and aerosol particles interstitially between hydrometeor particles fed back to radiative heating rates and photolysis. As such, the model treated the first and second indirect aerosol effects, the semidirect effect and both cloud absorption effects. Cloud and precipitation microphysics were solved explicitly, as described in Jacobson et al. [2007], who relied on microphysical algorithms derived significantly in Jacobson [2003]. The primary size‐ and composition‐resolved hydrometeor processes included (1) condensation/evaporation, sublimation/deposition: (2) coagulation among liquid, ice, and graupel and their aerosol components, (3) contact freezing of drops by size‐resolved interstitial aerosol particles, (4) heterogeneous plus homogeneous freezing of liquid drops, (5) coagulation of hydrometeor particle and their aerosol inclusions with interstitial aerosol particles, (6) breakup of temporary‐coalesced liquid drops, (7) sedimentation of liquid, ice, and graupel hydrometeor particles and their aerosol inclusions, (8) coagulation of precipitation hydrometeor particles with interstitial and below‐cloud aerosol particles (washout), (9) removal of precipitation and their aerosol inclusions (rainout), (10) below‐ and in‐cloud evaporation/sublimation to smaller hydrometeor particles and aerosol cores, (11) gas washout, (12) aqueous chemistry within liquid cloud and aerosol particles, and (13) heterogeneous chemistry on ice crystals. Since the model tracked inclusions within hydrometeor particles, the radiative transfer calculation accounted for the optics of all absorbing aerosol constituents within aerosol particles and within and between cloud and precipitation particles. Thus, the model treated both the microphysical and radiative effects of aerosols on clouds and precipitation. The radiative transfer code [Toon et al., 1989] was used to solve for atmospheric heating rates and actinic fluxes over each of 694 wavelengths/probability intervals with gas absorption coefficients from Jacobson [2005]. Aerosol and cloud optical properties were calculated by integrating spectral optical properties over each size bin of each aerosol and hydrometeor particle size distributions. In aerosol particles, black carbon for optical calculations was treated as a core surrounded by a mixed shell [Ackerman and Toon, 1981]. UV and short‐visible absorption by organic carbon were accounted for, allowing for treatment of brown carbon radiative effects [Jacobson, 2010a]. Since all aerosol component concentrations were tracked within each size of each hydrometeor particle type (liquid, ice, and graupel) throughout the evolution of clouds and precipitation, it was possible to calculate cloud absorption due to BC inclusions in clouds. The most physical method of treating absorbing inclusions is assuming they are randomly dispersed throughout a hydrometeor particle. This method was applied in the present model with the iterative dynamic effective medium approximation [Chylek et al., 1984], as described and applied in Jacobson [2006, 2010a, 2012].

6 1‐D Simulation Results Three time‐dependent simulations were run. In the first (Sim1), all processes were accounted for. In the second (Sim2), all processes except absorption by BC inclusions in hydrometeor particles were accounted for (absorption by BC interstitially between cloud drops and outside of the cloud was still accounted for, and BC was still tracked through hydrometeor inclusions to ensure all other feedbacks operated in exactly the same way as in Sim 1). The third simulation (Sim3) was the same as the first, except that the BC concentration was set to zero in all layers (thus, no aerosol BC outside the cloud, BC interstitially between cloud drops, or BC inclusions within drops). During Sim1, direct radiative forcing of BC inclusions and of BC aerosol was also calculated with a separate radiative transfer call during each time step of the simulation. The solar zenith angle during the simulations ranged from 22.86° at the start to 40.92° after 5000 s. The time step used for all aerosol, hydrometeor, and radiative processes was 5 s. Figure 5 shows the most relevant elements of the simulations. The air was initially supersaturated from 0.4 to 0.75 km (Figure 5a). In Sim1 and Sim2, BC was initially moderately concentrated near the surface, and its concentration declined with height to 200 ng/m3 at 4 km (Figure 5a), an elevated concentration similar to in Beijing [Ramana et al., 2010]. The initial particle number (down to 2 nm diameter) at the surface was 97,700 cm−3 and declined with increasing height. During Sim1 and Sim2, the BC concentration evolved in the cloud and aerosol over time as nucleation scavenging and aerosol‐hydrometeor coagulation moved BC from interstitial aerosol to hydrometeor particle inclusions. As the cloud dissipated, though, aerosol inclusions aside from those that had been lost by precipitation (about 18% of all BC during the simulation) were released back to the aerosol size distribution, as illustrated in Figure 5b. Figures 5c and 5d show the time‐dependent change in the vertical profile of cloud liquid content (LWC) during Sim1 and Sim2, respectively. The lesser growth of the cloud during Sim1 was due solely to solar absorption by BC inclusions in hydrometeor particles (CAE I) and the subsequent increase in temperature (Figure 5e) and decrease in RH (Figure 5f) in the cloud. The cloud absorption optical depth due to BC inclusions in Sim1 (Figure 5g) decreased as the cloud dissipated. At 5000 s, some cloud inclusion absorption remained, as one cloud layer persisted because thermal‐IR cooling slightly exceeded solar warming in the layer. The cloud top cooling also enhanced turbulent mixing of moisture from below. Both factors increased the RH in the layer, enhancing its growth. Such a feature has been modeled in 3‐D with respect to anvils [Ackerman et al., 2000]. In some other simulations run, the cloud dissipated entirely in the Sim1 case but not in the Sim2 or Sim3 cases. In all cases, the same overall conclusion was found. The more rapid dissipation of the cloud in Sim1 is reflected by the significantly faster and sustained decrease in cloud optical depth (COD) in Sim1 than in Sim2 or Sim3 (Figure 5h). In all cases, the COD began to rebound starting after 1000 s due to the persistence and growth of one cloud layer for the reason described above, but the COD during Sim1 remained much lower than during Sim2 or Sim3. All other cloud layers are dissipated. Figure 5i shows the time‐dependent top‐of‐the‐atmosphere (TOA) direct radiative forcing (DRF) due to BC inclusions and, separately, aerosol during Sim1. Figure 5j shows the corresponding surface temperature change. The DRF due to BC inclusions in hydrometeor particles was calculated by running the radiative transfer code a second time through the column model each 5 s time step during Sim1 with absorption and scattering due to BC inclusions removed from the hydrometeor particle optical calculations. The DRF of BC aerosol was calculated by running the radiative transfer code a third time through the column model each time step with the optical properties of BC aerosol above and below the cloud and interstitially between hydrometeor particles within the cloud shut down. Figure 5i shows that the DRF due to BC inclusions decreased over time then steadied out due to the persistence of the one cloud layer. The DRF due to the sum of BC inclusions and BC aerosol (interstitially between cloud drops plus outside of clouds) was larger than but followed that of BC inclusions because the DRF of BC aerosol alone did not change much during the simulation. The reason is that most of the column aerosol BC was above or below the cloud during the entire simulation (Figure 5b). In addition, although BC aerosol within the cloud was converted to hydrometeor particle inclusions during the growth of the cloud, much of that BC was returned to aerosol as hydrometeor particles evaporated and released their cores back to aerosol. The decrease in BC DRF upon dissipation of the cloud and release of aerosol cores to the air was due to a decrease in optical focusing together with the removal of about 18% of column BC by precipitation. By optically focusing more light into a cloud drop that eventually hits an inclusion, BC inclusions in drops amplify the BC mass absorption coefficient compared with uncoated BC, by up to a theoretical peak of 2.5–4 [Fuller et al., 1999]. Mikhailov et al. [2006] measured an amplification factor of up to 3.5 for a soot‐water drop mixture at 100% RH. Chylek et al. [1996] modeled an amplification factor of BC inclusions within a single drop for a particular case of 2–2.5. The maximum enhancement factors for the coated core Mie code used here similarly range from 2 to 4, depending on the BC versus total particle diameters [e.g., Jacobson, 2012, Figure 13]. These data and theoretical calculations suggest that, as cloud drops shrink back to aerosol cores, their amplification factor decreases significantly. This is the reason for the decrease in DRF in the cloud‐inclusion‐only case in Figure 5i. Whereas total TOA BC radiative forcing (forcing due to BC in hydrometeor particle inclusions plus forcing due to BC in aerosol particles outside of the cloud or within the cloud at the RH of the cloud) decreased during the simulation, particularly during the first 1000 s, the ground surface temperature due to both BC inclusions alone and to all BC (accounting for indirect effects, semidirect effects, and CAEs) increased continuously since, as BC inclusions burned off the cloud, sunlight poured through the remainder of the cloud, warming the remaining cloud layers and surface rapidly. Figure 5k shows the anticorrelation between instantaneous BC DRF and surface temperature. The fact that the slope is negative indicates that instantaneous (rather than average) radiative forcing is not a proxy for temperature change in the presence of absorbing inclusions in clouds. Figure 5i also shows the TOA irradiance change between Sim1 and Sim2 and also that between Sim1 and Sim3. Irradiance changes are differences between simulations, whereas radiative forcings are differences calculated during one simulation when the radiative transfer code is called a second time with a parameter removed. Although irradiances changes in Figure 5i are anticorrelated with temperature change during some stages of cloud evolution, they are correlated proportionately with the temperature changes in Figure 5j over much of the simulation. Thus, irradiance changes may be a better surrogate for surface temperature changes in the case of absorbing inclusions in cloud drops than direct forcing. One result of these simulations is that instantaneous DRF does not correlate positively with surface temperature change for cases when clouds with absorbing inclusions are present. The reason is that, when absorbing inclusions burn off cloud drops, the DRF of cloud‐drop inclusions decreases because drop heating decreases drop liquid water content and optical focusing, causing drop solar absorption itself to decrease, yet the ground temperature escalates rapidly because the cloud burn‐off now allows sunlight to pour to the surface. Clouds with absorbing inclusions are present over much of the world's land and coastal areas, near shipping lanes and near aircraft flight tracks. Results from these 1‐D simulations help to explain the warming effect of biomass burning found in this study since BB emits significant quantities of absorbing BC, tar balls, and other BrC.

7 Conclusions This paper examined the effects on climate with 20 year transient simulations of open biomass burning (BB) when heat and moisture fluxes, all major gas and aerosol constituents (including black and brown carbon, tar balls, among others), cloud absorption effects, (CAEs), semidirect effects, and indirect effects of BB were accounted for. Biomass burning was calculated to cause 20 year global warming of ~0.4 K because CAE I (~32% of BB warming), CAE II, semidirect effects, AHFs (~7% of BB warming), AMFs, and tar balls from biomass burning together outweighed indirect effects, contrary to previous biomass burning studies that did not treat CAEs, AHFs, AMFs, or brown carbon. AHFs from all sources and AMF + AHF from power plants and electricity use each accounted for statistically significant +0.03 K global surface air temperature warming. Modeled temperature changes from all AMFs alone and AMF + AHF from transportation were not statistically significant. BB was estimated to cause ~250,000 (73,000–435,000) premature mortalities/yr, with over 90% due to particulate matter and the rest due to ozone. The combination of fossil fuel soot, biofuel soot and gases, BB, AHFs, AMFs, anthropogenic greenhouse gases, other anthropogenic aerosol particles, and urban surfaces together was found to account for nearly all net observed global warming to date. Solar plus thermal‐IR direct radiative forcings were calculated as follows: all‐source BC outside of clouds and interstitially between cloud drops at the relative humidity of the cloud: +0.52 W/m2 (94.3% solar); BC inclusions in cloud drops and ice crystals: +0.06 W/m2 (99.7% solar); all AHF globally without evaporating cooling water: +0.033 W/m2 (thermal‐IR); AHF with evaporating cooling water: +0.027 W/m2 (thermal‐IR); and all AMF globally: +0.009 W/m2 (thermal‐IR). In addition, global direct radiative forcings were calculated as in Table 1 for natural H 2 O, CO 2 , O 3 , N 2 O, O 2 , CH 4 , CO, and CH 3 Cl. One‐dimensional simulations helped to explain the strong effect on cloud burn‐off of BB and BC from other sources. They showed that when absorbing aerosols exist in clouds, instantaneous direct radiative forcing (DRF) and surface temperature change are anticorrelated because when absorbing aerosol burns off a cloud, the aerosol DRF decreases due to a decrease in optical focusing, yet surface temperature escalates rapidly due to the pouring in of sunlight to the surface. As such, particle burn‐off of clouds may be a major underrecognized source of global warming.

Acknowledgments This work was supported by U.S. Environmental Protection Agency grant RD‐83337101‐O, the NASA High‐End Computing Program, and the National Science Foundation. We thank Greet Maenhout for providing additional detailed information from the EDGAR emissions inventory. Data and results from this study may be obtained upon request from the authors.