Improving precision in IBU calculations

[Lyrebird_Cycles]

My thanks to MHB for suggesting this and pointing me to the paper by Malowicki in my thread on CO2 calculation.

I was going to leave it be except a couple of beers made recently on a new rig using the Tinseth formula turned out to be wildly off their targets so I thought "**** it, I have to be able to come up with something better than that". Down the rabbit hole I went.

The Malowicki paper turns out to be part of his Master’s thesis and part of a recent surge in interest in the kinetics of hop utilisation, probably spurred by the surge in interest in all things hop for which smaller brewers can take much credit.

Papers open as PDFs (I hope)

Elena Hop utilisation Paper
Huang Hop Kinetics paper
Jaskula Isomerisation
Kappler degradation
Malowicki Thesis
McMurrough HPLC bitterness


Taking this information and using Glen Tinseth’s pioneering work on hop utilisation as a base, I’ve tried to develop a model of the reactions involving hops in the hot stages of brewing: kettle, whirlpool and cooling. The model leads to a calculator, I’ll post the calculator and brief explanatory notes later.

Tinseth’s calculator assumes that isomerisation is a first order effect but takes took no account of degradation kinetics. His model is also presented in a way that makes it difficult to unpack all the processes involved. The kinetic approach used by Malowicki avoids both these problems but requires a little more knowledge of how reaction kinetics work: if you are not familiar with the work of the great Svante Arrhenius, Wikipedia is your friend.

Although at first glance this is a complicated way of doing things, it actually simplifies the building of the model, because we can now use a sequence

Hop addition -> alpha acids availability -> isomerisation of available alpha acids -> degradation of available iso-alpha acids.

These are all independent first order effects so the output of each is the input to the next and the net effect of two or more processes acting on two or more components is simply the sum of those same effects acting on each component.
The models used for the isomerisation and degradation steps are minor refinements of the Malowicki models. The model used for hop addition and alpha acid availability is derived from Malowicki and from Tinseth.

Taking the hop addition first, it is evident that a major contributor to hop utilisation being less than 100% is that alpha acids are insoluble but surface active so a percentage of the calculated addition of alpha acids is unavailable for reaction and a percentage of the iso alpha acids produced are also lost (but recoverable from trub etc). Some of this effect appears to be due to surface effects in the equipment used, as an example Malowicki estimates that 20% of the alpha acid solution added to the experiments adhered to various surfaces and became unavailable. In addition some of the effect appears to be related to the presence of materials such as proteins and polyphenols which presumably scale with wort total concentration. It isn’t the sugar content of the wort as such, many experiments have shown that sugar has no direct effect. It is not known whether wort oligosaccharides and polysaccharides have any effect.

Another major effect is that of pH, literature data showing a range of changes in utilisation with pH from the oft quoted Increase by 1% for each 0.1 pH unit to data for one hop variety in the Elena paper changing utilisation by over 6% for each 0.1 pH unit. The complication here is that although isomerisation rate has not been shown to be pH dependent, degradation of iso alpha acids has been, which should mean that iso alpha recovery decreases with pH rather than increasing. The difference is likely to be due to some combination of solubility effects (alpha acids being much more soluble in higher pH), extraction effects and polyphenol interaction effects but data on these is sorely lacking. In the end I chose a simple equation to scale the initial availability with pH and chose a rate of 1% increase for each 0.052 pH units (because much of the available data is at pH 5.2 so it was easiest to use a known availability at that pH and then divide by 5.2 to scale). This combines with the increased degradation with pH (which will be incorporated below) to give roughly 2% increase in net utilisation per 0.1 pH unit.

Tinseth attempts to account for the known effect of lower utilisation in higher strength worts with his “Bigness factor”, given as
1.65 x 0.000125^(1-SG).
Unpacking this, it appears to be based on a wort of 13.75 Plato having a bigness factor of 1. The loss due to wort concentration scales either side of this so that “wort” of 0 oP would have a bigness factor of 1.65. The only other publication I can find which gives empirical data for the effect of wort concentration is the Kappler paper which paints a completely different picture: the data given show a linear reduction in recovery of added iso alpha from 90% at 10oP to 52% at 18oP, thus implying that worts below 8oP would have no further effect and worts above 26oP would have no recovery at all, both of which seem unlikely.
I decided to redraft the Tinseth curves, primarily because I don’t work in SG but also to allow me to separate out the effect of higher gravity worts boiling at higher temperatures (which is accounted for elsewhere). The re-drafted version deviates slightly from Tinseth’s curve due to the fact that SG is not a simple function of oP.

The overall equation for alpha acida availability is Avnet = Av0 * 1.038^-P * pH / 5.2 where Av0 is a factor to account for initial loss and P is the wort concentration in Plato.

I chose the factor 1.038 to minimise the deviation at the original concentration of 13.75oP, it gives a deviation less than 1% at all wort concentrations less than 18oP. Interestingly the section of the new curve between 14oP and 18oP fits the data in the Kappler paper between these values as well as the Tinseth curve.

Each of these factors can be modified for different brewing conditions: the initial availability factor first, since it exists purely to account for differences in hop utilisation. It could, for instance, be used to account for hops being added in different forms: addition as dried flowers would be reflected in an Av0 which is about 10% lower than that for addition of pellets, etc., although it is probably simpler to discount the label alpha level by the appropriate factor when adding the hops. Other effects which will affect this factor are such things as brewery size (to a first approximation, surface effects scale as the 2/3 power of brewhouse volume) and surface composition. Some other brewhouse effects such as the use of a calandria will be covered below under temperature.

The base of the concentration scale factor (1.038) can be also be altered, for instance to account for different wort compositions. At present there is little data to support this but part of the reason for making this explicit in the model is that it will allow such data to be collected, for example by calculating the effect of boiling two worts which are similar except for total protein concentration and quantifying the effect this has on hop utilisation. One of the reasons I have included it is that it at least conceptually accounts for the observed deviation of alpha acid utilisation in the first ten minutes of boiling: assuming that at least part of the concentration effect is due to alpha acids adsorbing onto the surface of the hot break, the effect will not stabilise until the hot break is fully formed.

I will post again on the kinetics of the next steps in the chain (the isomerisation and degradation reactions) in the next couple of days, meanwhile you can brush up on Arrhenius.

 

[MHB]

Malowicki's paper started me looking at hops/bitterness in a whole new way.
The one big (and perhaps insurmountable) problem shouldn't be getting a good idea what is in the wort at the end of the boil - but what is left at the end of fermentation, by definition what an IBU is.

There hasn't been a lot of research done on trans Iso-alpha breakdown products, but quite a lot of anecdotal evidence that long boils or early hop additions (mash hopping, first wort hopping...) lead to a very smooth round bitterness, I suspect that the longer time hot leads to more of the trans iso products and that they are playing a role.

This should be fun, looking forward to where it goes
Mark

 

[Lyrebird_Cycles]

You realise this is all your fault, don't you?

I'm not even going to try to address the issues post hot side, too many variables.

The Kappler paper contains some stuff about differential rates of trans - and co - isohumulone recovery at different boil times. IMO the theory that hop bitterness quality is primarily due to low coisohumulone is a bit suspicious. Yes, hops with good bitterness texture are often low cohumulone but that might be a case of confusing correlation with causation, my favourite analogy is that rangers have high skin cancer rates but it isn't because red hair doesn't shade the skin as well.

It may be that longer boil times and especially adding hops when there are lots of other surfaces and polyphenols present (mash hops) primarily influence the polyphenol interactions. The AWRI has been working on polyphenol chemistry in wine for a decade, they're learning a lot. The chemistry in the boil is even more dynamic and complex, we'll work it out one day.

Partly this is personal: I don't have quite a good enough palate to be a wine show judge*, but I'm very good on tannin structure. I have a much better palate for beer and am especially good on bitterness texture. I think the two are related.

* I did the AWRI advanced wine assessment course and scored 2 SDs above the average (of trained winemakers who want to be show judges) on the reds but that's still not good enough for the show circuit, those guys are really good.

 

[MastersBrewery]

Thanks for melting my brain before my second coffee. I'll have to read all this atleast twice more, and with caffeine levels in the green. Most of this will inevitably go over my head, however for me IBU doesn't quatify all that hop additions bring to a brew. Flavour and aroma are left to ones imagination; colour, bitterness, and gravity, go only part way in describing a beer, but are the only scales we have.

 

[Lyrebird_Cycles]

Agreed, hence the comments on polyphenols above. This was brought home to me many years ago when a brewery for which I was working changed from using primarily POR as bittering hop to the then new NZ Super Alpha (AKA Dr Rudi) to save money. In my opinion the beers went to **** overnight.

There's a small difference in cohumulone level and the NZ hops have that nasty green pinecone character that's so fashionable right now but not much of that survived the 90 min standard boil. What did survive is what I called "triploid bitterness", a short, hard, coarse bitterness that was totally inimical to the styles of beer we had established.

 

[klangers]

Good read. This is bringing me back to the last time I played with reaction kinetics, but that was on supercritical transesterification of vegetable oils to produce biodiesel.

I'll have to look at the individual papers later.

 

[GibboQLD]

This is right up my alley, thanks for the papers & post -- looking forward to some "light" reading in my lunchbreak.

 

[Lyrebird_Cycles]

Back again, now we proceed to the heart of the matter, the isomerisation of the alpha acids to produce iso alpha acids. This is fortunately the simplest part of the process because as far as is known the kinetics equation in the Malowicki paper is a full description of the process:

d [alpha] / dt = -k1[alpha]

k1 = A1 e-Ea1/R.T

The activation energy (Ea) given by Malowicki is 98.6 kJ / mole. He gives the pre exponential factor (A1) as 7.9 x 10^11 but that is taken from a line of best fit through multiple temperature data points. For my model I have recalculated the pre exponential factors to fit the line through the data point closest to the normal conditions for brewing (T = 373 K) so the new pre exponential factor is A1 =7.2 x 10^11.

Don’t worry about the seemingly outlandishly large number, the pre exponential factor is the exponent of the y intercept on a plot of Ea/R vs 1/T so it represents the activity that would be present at infinite temperature.

Now all we need to do is work out the temperature profile to which the alpha acids are subjected and we can calculate the net isomerisation. The proper way to do this is to integrate the expression given in the first equation over the temperature curve but this would mean repeating this exercise every time you make a change and Excel doesn’t do calculus. If you have access to Matlab or similar be my guest but for the rest, you’ll have to trust my calculus.

Fortunately there’s an easy approximation at hand: since our isomerisation curves are first order logarithmic curves, as long as the rate of change of temperature is more or less constant we can integrate over a time period where temperature is changing by taking the logarithmic mean of the rate factor at each of the start and finish temperatures and applying that to the known integral of the rate equation. Excel can handle this manipulation.

The concept of the logarithmic mean will be familiar to anyone who has done heat exchanger sizing calculations. The log mean of two numbers is given by

Mean (x,y) = y – x / ( ln(y) – ln (x)) or y – x / ln (y/x) The two expressions are equivalent.

All that remains is to work out the temperatures to use. I have included factors in the calculator to account for the effect of boil vigour and kettle design, atmospheric pressure due to altitude and changes in Mean Sea Level Pressure and boiling point elevation with wort concentration.

Boil vigour and kettle design is the largest of these so we’ll tackle that first. Wort boils are normally conducted quite vigorously, principally in order to drive off unwanted volatiles like DMS. Normally this is stated as percentage volume loss per hour and figures of up to 10% per hour are common. The boiling point normally given is the equilibrium vapour point but a vigorous boil is by definition not in equilibrium so the boiling temperature will be higher. Another way of thinking of this is that the boiling temperature is the average of the boiling point of the liquid and the temperature of the steam which is in the process of travelling through the wort on its way to the surface. Since the steam temperature is higher than boiling point, the weighted average of the two temperatures is also higher than boiling point, how much higher will depend on such things as the depth of the kettle. In commercial scale breweries it is common for the boil to occur at about 102 oC.

The reaction rate is thus also higher than the reaction rate that would obtain at boiling point. For greatest accuracy we would need to do a weighted mean of the reaction rates for each isothermal mass in the kettle. The use of calandria boiling exaggerates this effect as the temperature of the wort in the calandria is much higher than the average temperature. The trouble with this approach is that the calculator would have to be changed for each kettle design. A simplification is to roll all these factors into one “boil vigour” factor and add it to the calculated boiling point to give the reactive boiling temperature. This figure is explicit in the calculator, by default it is set to 2 oC which corresponds to an average commercial boil but it can be moved up or down to suit the kettle in use.

Altitude and MSLP also affect the boiling point, for the simple reason that the boiling point is the equilibrium vapour point so when the pressure goes down so does the boiling point. Altitude is the larger effect in most cases but also the easiest to deal with since you aren’t likely to move your brewery often. The easiest way to handle this is to add the effects of pressure change with altitude from the expression

Palt = P0 * (1 - Alt * 2.25577 x 10-5)5.25588

using the local MSLP (which by definition is adjusted to sea level) as P0 and apply the net pressure to the Antoine equation:

TboilA = 1730.63 / (7.1962 - log10(Palt)) - 233.426 with Palt in kPa and Tboil in oC.

So that’s included in the calculator. To use this fully you’ll need to look up the BOM for local MSLP on brew day and use that as P0.

That being said neither of them is large enough to lose sleep over, they can be safely set to zero altitude and 1013 hPa unless you are on top of Kosciusko or in the eye of a cyclone and I would suggest they’re both bad places for brewing anyway.

Lastly comes the effect of wort solids on boiling point. Boiling point elevation is a colligative property, meaning it depends on the number of non-water molecules present. In wort this is principally the sugar so we need to know the average molecular weight of wort solids. I have simply worked backwards from known boiling points of worts at moderate Plato levels to give an average MW for wort solids as 320, which combined with the boiling point elevation constant for water of 0.51 C / mole gives a boiling point expression.

TboilC = TboilA + 0.016*oP.

In some commercial situations where a lot of glucose syrup is used the average MW could be well below the figure given so the factor0.16 would need to be adjusted upwards.

For the sake of completeness the calculator includes a factor for boiloff (in % per hour) and calculates the boiling temp at the beginning and end of the boil. If makeup water is added half way through the boil this complicates things, you can either split into two part boils or take the average of the concentrations before make up addition and at end of boil.

The rate constant used for boiling is the log mean of the rate constant for boil start temp and boil end temp. Boil length is user input.

Whirlpool and cast out / cooling calculations are handled a similar way with a few caveats. It is assumed that the temperature falls to boiling point as soon as the boil stops, then declines linearly to the whirlpool end temperature. The rate constant used for this section is again the log mean of the two rate constants.

For cast out / cooling, things are a little more complex as all the wort sees the start temperature and only the last runnings sees the end temperature. If the flow rate is constant the integral is computable but again excel doesn’t do calculus. The compromise position is a double log mean calculation.

At this stage I intend handling first running adds by calculating on three points: first running (T = mashout, oP = first running), end running (T = mashout, oP= boil start) and boil start (T = boil, oP = boil start). This requires an extra 4 user inputs and it only gets used if you do first running hops so I’m not sure if it’s worth it.

Mash adds are a whole new problem, as I see it there’s too much interference from grain solids for it to be tractable.

 

[GABBA110360]

**** thats bitter

 

[Lionman]

I wish this was written in english :blink:

If you can tell me at the end the perfect way to use hops I'll buy you a beer. :lol:

Seriously though, nice work. I think.

 

[MHB]

Ya the limitations of excel, I stuck the Malowicki's integrated equation (12) into a spread sheet with 400, 1 minute intervals.
Just fine for 1 hop addition, just copy and paste the whole slab for a second, third... hop, then sum the lot. Pushed out to quite a few rows but got sensible answers.

Plotting the Alpha, Iso-Alpha and Trans Iso all on one graph then playing around with the temperature can be quite startling, especially when you get to higher temperatures (say 115+ oC). Get too hot and the Iso is denatured as quickly as it forms.

Beware of the solubility limits, you can stick as much Alpha as you like into the kettle, a fair fraction of it might even isomerise but as soon as we cool to fermenting temperatures there are going to be very real limits on solubility (say about 100 IBU), then a fair wack of that will get lost during the ferment, even more on cooling to packaging and serving temperatures. Personally I would be very skeptical about any claims much over 85 IBU, no matter what we do.

My biggest educational regret is that I didn't do a lot more maths when I had the chance, calculus is about my upper limit, even then its a plod and a lot of dredging of the memory banks.
Mark

 

[Lyrebird_Cycles]

Yes, I thought about trying to incorporate a factor for solubility but there are several problems.

Alpha acids are very poorly soluble, iso alphas are somewhat better. There is good analytical evidence that beers with iso alpha levels above the alpha solubility limits exist. This implies that there is a dynamic interchange between the alpha acids in solution and the alpha acids stuck elsewhere in the system (presumably some still in the hops and some on the surfaces of the kettle and the trub flocs).

I know of no study which has made quantitaive determinations of the effect of this.

 

[yankinoz]

Thanks for your work.

Something that's indirectly related: where did the common wisdom come from that Rager is best for partial boils, Tinseth for full boils? I can't see it in the respective formulas. Arguably the different curves they've assumed for utilisation could affect their respective accuracies in high and low gravity boils, but if that's the case the recommendation would best be phrased as gravity-dependent.

 

[MHB]

This is the experiment that got me started thinking about bitterness calculations, its not very technical but well worth reading.
http://www.morebeer.com/brewingtechniques/library/backissues/issue7.1/bonham.html

I have always thought that "Boil Vigour" could be replaced with rate of evaporation, at least that gives some quantifiable indication of the amount of energy being pumped into the wort.
Bigness, could be viewed as boiling point, as the BP goes up proportionally to the amount dissolved into the wort.

Mark

 

[Lyrebird_Cycles]

I have always thought that "Boil Vigour" could be replaced with rate of evaporation, at least that gives some quantifiable indication of the amount of energy being pumped into the wort.
Bigness, could be viewed as boiling point, as the BP goes up proportionally to the amount dissolved into the wort.

Mark

I thought about using raw evaporation rate but then I'd need another factor to compensate for other kettle design effects. It seemed to me that rolling them into one variable would suit the majority of users.

It would be a trivial matter to separate them out if that's what people want.

If you read part 2 in detail you'll see that I have already separated out the part of the wort dissolved solids effect (aka "bigness") that is ascribable to increased boil temperature, it is

A: only a small part of the overall effect because the average molecular weight of wort dissolved solids is quite high (wort is generally < 1M).

B: Usually works in the opposite direction because increased temperature increases isomerisation rates. As you will see when I finish part 3, pH also plays a role here.

 

[Lyrebird_Cycles]

Part 3 won't be up today, I realised that I'd made no provision for cooling in place. Not sure how to do both types of cooling in one calculator.

 

[Lyrebird_Cycles]

This last section is largely about the kinetics of degradation and as such is the major reason for doing this in the first place given that the Malowicki paper showed that degradation has a measureable effect on final IBU. No previous calculator has tried to incorporate this.

The degradation kinetics in the Malowicki paper results in simple equations analogous to the eqs given in the post above:

d[iso] / dt = -k2 [alpha]

k2 = A2 e -Ea2/RT

with Ea2 and A2 given as 108.0 kJ /mole and 4.1 x 10 12 respectively. Again I’ve moved the curve to pass through the point nearest typical wort temperatures which changes the pre exponential factor A2 to 3.45x 10^12

These results were taken in buffer at constant pH. The Huang paper shows that this is a simplification because the activation energy changes strongly with pH, the difference being almost 20 kJ / mole for a change of 1 pH unit. The physical interpretation of this is that changing the pH alters the temperature dependence of the reaction rate: at higher pH a given change in temperature will have a larger effect on reaction rate.

If the changes in temperature dependence were the only effect of pH, all would be well and we could just change the energy term and be on our way. For this to be the case, a change in activation energy would have to have no effect on the pre exponential factor. Unfortunately a check of the curves in fig 8 of the Huang paper blows this out of the water: if the activation energy were the only thing to change, the lines would intersect at the Y axis
(1/T = 0 -> ln k2 = ln A2). The lines do intersect but the intersection point appears to be at about 1/T = .0024k -1 ( T = 417K)

We therefore have to take into account both a change in reaction rate and a change in the temperature dependence of that reaction rate. The easiest way to do this is to assume that the line at pH 5.2 in Malowicki’s results passes through the same point as the lines in Huangs results, and since we know both A2 and Ea2 for the Malowicki result we can use a simple transformation to find the Y axis intersection points for the other lines which will give the A values for each Ea value: Calculate the kinetic constant for the known equation at 417K ( = .09969 ) then substitute this into a transformed version second equation above:

k2 = A2 e -Ea2/RT so A2 = k2 e Ea2/RT (note the change of sign in the exponent)

And out pops the new A2. Again Excel can handle this manipulation.

The net result of this is that kinetic “constant” in the rate equation above will alter for each different temperature and pH but the form of the equations will stay the same. We can thus proceed as before, calculating a new k2 for each combination of pH and temperature over the process life of the alpha acids and using the log mean to integrate over any time period where either or both change, using the equation

[iso]t = [alpha]0 * k1/(k2-k1) * (e –k1t – e –k2t)

to obtain the iso alpha level at time t with [alpha]0 being the output of the alpha acid availability calculation outlined in the first post.

The beauty of this approach is that the isomerisation and degradation calculations can be continued as the wort is further processed meaning for instance the whirlpool time will affect the isomerisation and degradation of hops added during the boil.

This requires three rate equations: the alpha eq as previously, the iso equation, also as previously but with a different starting point, and a separate equation

[iso]t,deg = [iso]0 * e – k2t

To account for degradation of iso alpha acids formed in previous steps.

Lastly, I have made two changes to the procedures outlines in the two posts above.

One: at MHB’s suggestion, I have separated out the two factors influencing boil vigour so that the vigour automatically increases with boiloff rate.

Two: I have split the cooling calculation to allow cooling via heat exchanger or cooling in place.
The heat exchanger calculation assumes a constant flow rate so average residence time is one half the total cooling time. It uses a double log mean of the cooling temps at beginning and end because the temperature profile seen by each part of the cooling wort changes with its time through the heat exchanger. The temperatures involved are the hot side temperatures, it is assumed that the cold side is cool enough for the reaction rate to be effectively zero.
The in place calculation assumes that all of the wort is exposed to the same temperature profile and that the temperature linearly ramps from the start temp to a temperature at which the reaction rates are effectively zero. I have set this at 50 degrees C for convenience.

An excel file with the calculator is attached below, have a play and get back to me with any suggestions.

View attachment Bitterness Calculator_V2.xlsx

 

[Lyrebird_Cycles]

I was worried that the hop utilisation figure only used the starting OG so it understated the effect of concentration. Accordingly, I have modified the calculator, it now uses the end boil figure rather than the start for the effect of OG on hop utilisation.

New version attached:


View attachment Bitterness Calculator_V3.xlsx

 

[manticle]

I've become super lax about my additions of late - just going off guestimates and using my palate as a yardstick but understanding the mechanisms at work at this level of detail is fascinating.
Cheers for your work so far.

 

[MHB]

This is the old Excel spreadsheet I put together (about 10 years ago) it shows the interplay of the Alpha > Iso Alpha > Degradation.
Its basically equation 12 (with the sign reversed), it illustrates the conditions of the experiment, without any of the additions/modifications/improvements you have made.
Interesting work to date much appreciated.
Mark

View attachment IBU's.xls