Determining the optimal shaking rate of a reciprocal agitation sterilization system for liquid foods: A computational approach with experimental validation

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor ContentslistsavailableatScienceDirect

Food and Bioproducts Processing

jo u r n al h om ep a g e :w w w . e l s e v i e r . c o m / l o c a t e / f b p

Determining the optimal shaking rate of a

reciprocal agitation sterilization system for liquid foods: A computational approach with

experimental validation

Ferruh Erdogdu

, Mustafa Tutar

, Sigurd Oines

, Igor Barreno

, Dagbjorn Skipnes

aDepartmentofFoodEngineering,AnkaraUniversity,Ankara,Turkey

bMechanicalandManufacturingDepartment,MGEPMondragonGoiEskolaPoliteknikoa,Spain

cIKERBASQUE,BasqueFoundationforScience,Spain

dDepartmentofProcessTechnology,Nofima,AS,Norway

eCSCentroStirlingS.Coop,Aretxabaleta,Spain

a r t i c l e i n f o

Articlehistory:

Received22February2016

Receivedinrevisedform19July2016 Accepted22July2016

Availableonlinexxx

Keywords:

Cannedfoods Optimization Modelling

Reciprocalagitation-shaking

a bs t r a c t

Anewcanningprocesswhereareciprocatingagitationiscarriedoutinhorizontallyoriented containershasbeenrecentlydemonstratedtoreduceprocessingtimesandenableenergy savingswithlessdegradationinthequalityofprocessedfoodproducts.Reciprocalagita- tionbyimposingadditionalforcesenhancesconvectivemixingwithincreasedproduction efficiency.Thereciprocalagitationusesthehorizontalaccelerationinadditiontogravity andsumoftheseforcesleadtoaconsiderableincreaseintheheattransferrates.Inthelit- erature,therehavebeenexperimentalapproachestoevaluateheattransferenhancement.

However,duetothebalanceamongtheseforces,theremightbeanoptimumreciprocal agitationratefortheincreasedheattransferdependinguponthephysicalpropertiesofthe liquidprocessed.Therefore,theobjectivesofthisstudyweretodeterminetheoptimum agitationratesbydevelopingacomputationalmodelforheattransfer.Forthispurpose,a multi-phasemodelsimulationwasperformedusingafinitevolumemethodbasedondis- cretizationofgoverningflowequationsforliquidandgasphaseinanon-inertialreference frameofmovingmesh.Experimentalstudiesformodelvalidationwerecarriedoutina reciprocallyagitatedretortusing98.2mm×115mmcanscontainingdistilledwaterwith 2%headspaceasamodelcase.Themodelresultswereinagreementwiththeexperimen- tal data,andtheoptimumreciprocalagitationratewasdetermined.Theresultsofthis studyaretobeusedtooptimizetheprocesswithrespecttoimprovethehealth-promoting compoundsofprocessedfoods.

©2016InstitutionofChemicalEngineers.PublishedbyElsevierB.V.Allrightsreserved.

1. Introduction

Traditional canning has been a convenient way and pro- videdageneralistandeconomicmethodforprocessingand preservationoffoodproducts.Consumerdemands forhigh qualityfoods,however,forcethefoodprocessorstoimprove

∗ Correspondingauthor.Tel.:+905338120686;fax:+903123178711.

E-mailaddresses:[email protected],[email protected](F.Erdogdu).

and innovate their processing.It isawell-known fact that the shorter the process time ata givenprocess condition, while still achieving the required safety for consumption, theless thedamagetothesensoryand nutritivequalityof thefoodproducts.Basedonthisconcept,followingtheuse ofretortsforcanning,theagitationretortswereintroduced

http://dx.doi.org/10.1016/j.fbp.2016.07.012

0960-3085/©2016InstitutionofChemicalEngineers.PublishedbyElsevierB.V.Allrightsreserved.

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor in1920swiththeagitationmechanismbasedonhorizontal

axialrotation(Atesetal.,2014).Verticalrotationofthecans waslaterintroducedwiththeend-over-endrotationprinciple (Clifcorn et al., 1950). Introductionofagitation mechanism incanningforliquidorliquid–solidparticlescontainingfood productswastheresultofacertaindisadvantageofthestatic retortsystems(Rosnesetal.,2011).Theprimarychallengein thestatic processingistheslowheatpenetrationresulting inalackofconsistencyinsensoryand nutritiveproperties (Ohlsson,1980).

Consideringtheeffectiveheattransferratesobtainedwith agitation,areciprocatinghorizontalagitationwithrapidback andforthmotionofthehorizontalcontainersinanoscillat- ingwayhasbeenproposedtoincreasetheheattransferrate further,andanagitatingretortwithhighfrequencylongitu- dinalmechanismwasdevelopedin2006(Atesetal., 2014).

Bothhorizontalandend-over-endbasedagitationretortssuf- ferfromthelimitationofthattheappliedforcestoenablethe motionwithinthe containerwere abalance betweengrav- ityandcentrifugalforces(Waldenand Emanuel,2010).Due tothisbalance,theagitationincreasestheheattransferrate uptoanoptimumwhilefurtheragitationmightnotaffector mightinfluencetheprocessinanegativewaydependingespe- ciallyupontheviscosityoffoodproduct.Adetailedanalysis andcomparisonamongthegravityandcentrifugalforcesfor thecaseofaxialrotationeffectsinhorizontalaxialrotation ofcanswerereportedbyErdogduandTutar(2012)andTutar andErdogdu(2012).Thereciprocalagitation,however, used thehorizontalaccelerationinadditiontogravity,andthesum oftheseforcesenabledaconsiderable increaseintheheat transferrateswithreductionsintheprocesstime(Waldenand Emanuel,2010).

Thefirststudiesinthefoodengineeringliteratureusing thereciprocatingagitationsystemswereexperimentalbased to demonstrate the possible process time reductions and improvementsintheheattransferrates.Bermudez-Aguirre etal.(2013a)demonstratedtheimprovementinheattransfer coefficientunderstaticandhorizontalgentle-rockingmodes.

Atesetal.(2014),forexample,comparedthenovelagitating retort and static retortprocesses forbacterial inactivation, anditwasconcludedthatagitatingretortprocesssignificantly loweredtherequiredprocesstime.ThestudybySinghetal.

(2015a)focusedonevaluatingtheheattransferenhancement underreciprocal agitation whileSingh et al. (2015b) devel- oped anexperimental methodology todetermine the heat transfercoefficientincannedparticulatefluidsunderrecipro- catingfrequenciesupto3Hz.SinghandRamaswamy(2015a) focusedontheeffectofproductrelatedparametersonheat transferwhileSinghandRamaswamy(2015b)determinedthe effectoftheorientationofcans duringreciprocating agita- tion thermal processing. Singh et al.(2016) introduced the conceptofreciprocalagitationprocesstoimprovethequal- ityofcannedgreenbeansduringthermalprocessing.Singh andRamaswamy(2016)carriedoutanoptimizationstudyfor theheattransferrateandreciprocationintensityforthermal processingofliquidparticulatemixtures.Thesestudieswere basedonexperimentalapproacheswhileasimilarsituation wasexploredbyLiffmanetal.(1997)andPeschetal.(2008)ina computational–theoreticalwayforconvectionduetohorizon- talshakingandheatedfluidlayerssubjectedtotime-periodic horizontalaccelerations,respectively.

Eventhoughtherewerecertainfindingsreportedforthe effect of reciprocal agitation on the temperature increase and enhanced heat transfer rate (Bermudez-Aguirre et al.,

2013a;Atesetal.,2014;SinghandRamaswamy,2015a,b,2016;

Singh et al., 2016),development ofacomputational model (withoneexceptionwheretheheattransfercoefficientbased lumpedmodelwithoutconsideringthetemperaturedistribu- tion wasintroducedbyBermudez-Aguirreetal.(2013b)and determiningtheoptimalagitationrateswerenotfocusedin detail. For determiningthe optimal conditions, oneexcep- tion was reported bySingh and Ramaswamy(2016) where the optimal conditionsofreciprocation intensity forliquid particulate mixtures were experimentally determined.The optimizationstudiesbasedonacomputationalmodelaresig- nificantsincethecomputationalmodel mightalsobeused also forprocess developmentpurposes. Therefore, the pri- maryobjective ofthis study wastodetermine theoptimal agitationrateinareciprocalagitationprocessusinganexper- imentally validated computational model. The secondary objectives were first todevelopa computationalnumerical modelforheatandmomentumtransferinsidetherecipro- callyagitatedcanstodeterminethetemperaturedistribution and velocitychanges and thenexperimentally validatethe model.

2. Materials and methods

Forthegivenobjectives,thestudyconsistedofexperimental and computationalparts.Intheexperimentalpart, awater filledcanwasprocessedinboilingwaterandagitatingcondi- tions.Inbothcases,thehorizontallyorientedcancontained waterasatestliquidtorepresentalowviscosityNewtonian liquid. Thetime–temperaturedataobtainedatthegeomet- riccenterinthefirstexperimentswereusedtodevelopand validatethecomputationalmodel,todecideuponthecom- putationalparameterswiththemeshindependencystudies.

Followingthemeshindependencystudy,thecomputational model was validated with the temperature data obtained under horizontal agitating conditions, and the model was appliedto horizontalagitation ratesfrom 20 to140rpmto obtain the agitationrateinthe directionofthe axisofthe horizontal canresulting in maximumheat transfer. In the reciprocalagitation systems,the crankshaft,usedtoderive the horizontalmotion,angularvelocityisrelatedtoengine revolutionspermin(rpm).

2.1. Experimentalmethodology

Thefirststepinthisstudywastodecideuponthecomputa- tionalparametersandtesttheaccuracyofthecomputational method. For this objective, an experimental study with a cannedwatersample(98.2mm×115mmcansfilledwithdis- tilledwaterwith2%headspace)wascarriedoutinaMicroflow 911 EAT Shakaretort (Steriflow, Roanne, France)in boiling water under stationaryconditions. The retort system was heatedbydirectsteaminjection,equippedwithapreheating tanktoprocesswater.Thisprocessedwaterwasthencircu- latedthroughaheatexchanger(onlyusedforcooling)tothe retortandspreadbyaperforatedplatetoobtainwaterraining overthecans.Thecanwasfixedinahorizontalpositioninthe boilingwater.Type-Tthermocoupleconnectedtoadatalog- gerE-ValFlez(Ellab,Copenhagen,Denmark)waslocatedatthe geometricalcenterusingringgasketsandlocking-receptacles.

Theexperimentalset-upwasshowninFig.1.

The canmaterial wasa steelsheet witha thicknessof 0.19mmandthermalconductivityvalueof15–16W/m²K.This enabledtheassumptionofthenegligibleconductioneffectof

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor Fig.1–Experimentalset-upwherethecanwasfixedina

horizontalpositionwithinstalledtype-Tthermocouples.

canwallontheheattransfer,andthemediumtemperature wasacceptedtobethecansurfacetemperature.

In the second group of the experiments, the same can wasusedunderhorizontallyagitationconditionswherethe shaking rate changed from 0rpm (in the first 26s of the process) to 80rpm (at the 35s of the process). Since the horizontal-acceleratedagitation rateswere includedinthis

modelvalidationpartofthestudy,theagitationrateswerefirst convertedtothetangentialvelocityvalues.Forthispurpose, sinceahorizontalagitatedsystemusedaslider–crankderived mechanism(Fig.2a–modifiedfromReaderandHooper,1982), displacement of the can during the agitation process was definedwiththefollowingequation:

x_p=

r−r·(1−sin(ω·t))+n·

1−

1−cos² (ω·t) n²

(1)

wheren=L/r,xisthedisplacement(m),ωistheshakingrate (rpm), t isthetime (s),risthe crankradiusofthe system (0.075m),andListhelengthoftheconnectivityrod(0.5m).

Since(n²»cos²(ω·t)),thisequationwasthensimplifiedwith:

xp=[r·sin(ω·t)] (2)

ThecomparisonofEqs.(1)and(2)forthechangeofdis- placementwithrespecttothe(ω·t)valuesdidnotshowany significantdifferenceresultinginverysamedisplacementval- ues.Therefore,thesecondequationwithitssimplifiedform

Fig.2–(a)Ahorizontalagitatedsystemwithaslider–crankderivedmechanism(modifiedfromReaderandHooper,1982);

(b)reciprocalagitationrateusedinthesecondsetofmodelvalidationexperiments;(c)tangentialvelocitychangeinthe transitionperidfrom0ro80rpmreciprocalagitationrate.

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor waschosentoderivethetangentialvelocityequationforthe

reciprocalagitationsystem:

v=r·ω·cos(ω·t) (3)

However,thisequationbroughttheissueofnon-zeroveloc- ityatthebeginningofthehorizontalagitationprocess.This resultedintheinstabilitiesandnon-convergenceproblemsin thenumericalcomputations.Topreventthisandtomakesure thattheprocessstartwitha‘0’velocityinitially,themovement was displaced with oneperiod, and the following velocity equationwaspreferredtostartwith:

v=r·ω·sin(ω·t)

limt→0v=0 (4)

Asexplainedintheexperimentalset-upformodelvali- dation,followingtheinitialsteady26s,reciprocalagitation transitionofthesystemwasfrom0toahigheragitationrate of80rpm[0.623m/s].Thistransitionperiodtook9sandwas assumedtofollowanexponentialincrease.Theexponential increasefrom0velocityto80rpmreciprocalagitationratewas preferredtoenablethesmoothtransitionbetweentheserates.

Theexponentialincreasewastheonlywayforthissmooth transitionandtopreventtheovershootafter9softhetran- sition period. Amongvarious trials, the linear increase for exampleresultedinasharptransitiontothe80rpm,which might be rather difficult to control physically. To conform thetransitionperiodwith80rpmofreciprocalagitationrate after9ssmoothly, thegivenvaluesbelowforthefollowing transitionstageequation,wherethe velocitychangeinthe transitionperiodwasshown,enabledthis:

v= 80[rpm]

60[s/min]·[1−exp(−k·t)] (5)

where(k=1)isthe constant(s⁻¹)and t=(t−26)(s).Fig. 2b showsthehorizontalagitationratethroughtheexperiments whileFig.2cdemonstratethevelocityprofilefrom0to80rpm inthetransitionperiodof9s(from26to35s).Thechangein thehorizontalagitationrate,asreportedinFig.2b,c,andthe variable–experimentallyrecordedmediumtemperaturewere usedinthemodelvalidationcasetocomparethenumerical resultswiththeexperimentaloneobtainedatthecentreof thecan.

Forbothcases,3-experimentswerecarriedout,theaver- agevalueswiththestandarddeviationwereusedinthemodel validation. Sincethe standard deviationsofthe average of thetemperaturechangebasedon thesethreeexperiments, additionalexperimentswereavoided.

2.2. Governingequationsandthecomputational model

The numerical methodology and full scale model experi- mentaltestingverificationsproposedausefulcomputational algorithmfordynamicmonitoringofheadspace(air)andliq- uid(water)interactionsthroughtheagitationandsolvedthe fluid-thermal energy interactions in order to optimize the reciprocalagitationprocess.Thetwo-phasevolumeoffluid (VOF)approachaccompaniedwiththefinitevolumemethod (FVM)basednumericaldiscretizationschemewasutilizedin thesimulationoftwo-phaseflowundervaryingphysicalcon- ditionsthrough unsteady, three-dimensionaland turbulent

flowsimulationsforthegivenRayleighnumberrangeover1E9 attheinitialphaseoftheheatingasexplainedbelow.

The basic mathematical model for the discretization processincludedthesolutionoffundamentalgoverningequa- tionsoffluidflowmotion,knownascontinuityequationand momentumconservationequations,i.e.,Navier-Stokes(N-S) equationsforincompressiblefluidinanon-inertialframe:

2.3. Continuityequation

∂

∂t +∇v_r = 0 (6)

2.4. Momentumequation

∂

∂t(vr)+∇(vrvr)=−∇P+∇¯¯r+ F (7)

wherewasthedensity(kgm⁻³),twasthetime(s),vrwas the relativevelocityvector ofafluidparticle (ms⁻¹),Pwas the staticpressure (Pa), ¯¯_r wasthe stresstensor (described below),F wastheexternalbody force(N)includinggravita- tionaleffectsandaccelerationduetothenon-inertialframe motion.Thestresstensor, ¯¯_rwas:

¯¯

=[(∇v_r+∇v^T_r)−2

3∇v៭_rI] (8)

wherewasthedynamicviscosity(Pas).Itwasdefinedtobe atemperaturedependentpolynomialfunction.Iwastheunit tensor,andthesecondtermontherighthandsidewasthe effectofvolumedilation.Thevolumedilationwasneglected inthesolutionssincetherewasnoeffectintheprocess.Energy conservationequation,alsosolvedforthepresentflow,was writtenintermsofrelativeinternalenergy(Er)andrelative totalenthalpy(Hr):

2.5. Energyequation

∂

∂tE_r+∇(v_rH_r)=∇(k∇T+¯¯_r)+S_h (9) where

E_r=h−P +1

2(v²_r−u²_r) (10)

Hr = Er+ P

(11)

Velocityevolutionswerethentransformedfromstationary torotatingframeusing:

vr = v− ur (12)

where vr wasthe relative velocity(ms⁻¹)arising from the meshmotion(velocityviewedfromthemovingmeshofthe oscillatory reciprocatingmotion),v wasthe absoluteveloc- ity (ms⁻¹)(velocity viewedfrom thestationaryframe),and

ur wasthelongitudinalvelocity(ms⁻¹)(velocityduetothe movingmesh).Theabovegoverningequationsweredirectly discretizedwithafinitevolumemethod(FVM)inconjunction withaninterfacetrackingmodel(asdescribedbelow)forthe air–liquid system. Reynolds-averaged Navier-Stokes (RANS) basedformoftheseequationswerediscretizedtogetherwith thetransportequationsofturbulencekineticenergyandits dissipationwithinthefinitevolumeschemebyusingtheRANS

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor basedk-εturbulenceclosuremodel(YakhotandOrszag,1986)

accompaniedwith the utilized interface tracking model at higheragitationrates.

Julienetal.(1996)reportedthe‘soft’and‘hard’turbulence conditionsledbythehigherRayleighnumbers,andthelower rangeswheretheRayleighnumber(Ra)wassmallerthan1E7, wascharacterizedbysoftturbulenceconditionsbyKooijetal.

(2015).TheinitialRayleighnumber forthefirstexperimen- talcasewas 2.66E9(thisvaluewas obtainedatthe 0.5sof themodelvalidationsimulationwheretherewasnomove- mentofthecan)whereitwasbeyondthesoftturbulencecase.

Eventhough theRayleighnumber, asaproductofGrashof andPrandtlnumber,determinedtheturbulenceinducement inthenaturalconvectionflowinthecylindricalcavity,local cellReynolds(ReL)numberchangewasalsocontrolled.Itwas upto30alongthesurfaceofthecanattheinitialphaseofthe stationarymodelvalidationcase(0.5s)whileitincreasedupto 100towardstheend.TheturbulenceReynoldsnumber(Re_T) was,ontheotherhand,highenoughtoresolvethepresent flowwitha turbulencemodel, withits instantaneous local valueupto18,400initiallyatthecentreofthehorizontalcylin- dricalcavity.TheRayleighandlocalReynoldsandturbulence ReynoldsnumberweredeterminedusingEq.(13):

Ra= g·ˇ·T·D³ (/)² ·Pr Re_L=V_c^1/3·v·

ReT=k²·

·ε

(13)

where g was the gravitational acceleration (m/s²), ˇ was thermalexpansioncoefficientforwater(1/K), Twasthemax- imumtemperaturedifferencebetweentheheatingmedium andtheinitialtemperatureofthesystem(K),Dwasthechar- acteristicdimension(diameterofthecylindricalcavity,m),Pr wasPrandtlnumber,andwasthedynamicviscosity(Pas), was the density (kg/m³), was the velocityencountered inagivencell(m/s),andV^1/3_c wasthecharacteristiclength ofthelocalcell,andkεwereturbulencekineticenergyand dissipationrate,respectively.BesideshighRayleighnumber encounteredattheinitialphaseoftheprocess,thelaminar flowconditionwasstilltestedforconvergenceduringtheini- tialtestsimulations,butthesetrialsresultedinconvergence problems.Therefore,basedontheRayleighnumberinforma- tionforturbulenceconditionsandconsideringtheresultsof theinitialsimulationsformodelvalidationpurposes,thetur- bulencemodelwasactivatedinthesimulations.Inaddition, forthesimulationstudycarriedoutundersteadyconditions formodelvalidation–meshindependencystudy,theturbu- lenceReynoldsnumberwasaround80towardtheendofthe simulation.

Fortheturbulencemodel,thefollowingturbulenceparam- eterswereapplied:

- Initial turbulence intensity (I) was assumed to be 5%, Basedonthemaximumtangentialvelocityvalueof1.1m/s (obtainedbyEq.(4)),theturbulencekineticenergyvalue(k) was0.00453m²/s²:

k= 3

2·(v_max·I)²=0.00453 (14)

- Turbulencedissipationrate(ε)was0.025m²/a:

ε=C^3/4 ·k^3/2

L =0.025 (15)

whereLwastheturbulentlengthscale(0.002m),andCwas turbulencemodelconstant(0.09).

Thetrackingofinterface betweenair–water phaseswas accomplishedthroughthevolumeoffluid(VOF)methodpro- posedbyHirtandNichols(1981).Inthismodel,asinglesetof momentumequationswassharedbythefluidsandthevol- umefractionsofeachofthefluidsineachcomputationalcell weretrackedthroughthedomain.Thefieldsforallvariables andpropertiesaresharedbyphasesandrepresentsvolume- averagedvaluesaslongasthevolumefractionofeachofthe phases isknown ateach location. Thus, thevariables and propertiesinanygivencell areeitherpurely representative ofoneofthephasesorrepresentativeofmixtureofphases dependingonthevolumefractionvalues.Thevolumefrac- tionsofwaterandairinthecomputationalcellsumtounity.

Interfacetrackingwascarriedoutbysolvingcontinuityequa- tionforvolumefractionofoneofthephaseswhereairwas specifiedasprimaryphaseandthusthevolumefractionofthe liquidphasewassolved.InadditiontoVOFmethod,Ubbink’s compressive interfacecapturingscheme (Ubbink andIssaa, 1999)forarbitrarymeshes(CICSAM)wasalsoapplied.

For the numerical solution procedure, a finite volume method (FVM)basedsolver (Ansys Fluent V15, Ansys,Inc., Canonsburg,PA,USA)wasusedtosolvetheprecedingpartial differentialgoverningequationsofthepresenttwo-phaseflow problem.Intheproposedcomputationalmodel,thecollocated FVMwasemployedtodiscretizethegoverning3Dflow-energy equations.Alltherequiredthermaland physicalproperties forairand waterphaseswere temperaturedependent and reportedinErdogduandTutar(2012).Initially,waterinthecan inbothexperimentalconditionwasatrestandhadtheinitial temperatureof300.92and301.38K inthesteadyandagita- tioncases,respectively.Whiletheboilingwatertemperature andvariablemediumtemperatureswereusedinthemodel validationsimulation,auniformconstantwalltemperatureof (Tw=373.15K)wasusedtodeterminetheeffectofreciprocal agitationratesonthetemperatureevolution.Forthecaseof agitationprocess,theheatingmediumtemperaturewasvari- able,butitwasstillusedtobeasaconstantwalltemperature overthecansurfaceduetotherapidmovementofthecandur- ingtheprocess.Overtheinitialperiodoftheagitationprocess wheretherewasnoreciprocalmovementofthecan(Fig.2c), theheatingmediumtemperatureandtheinitialtemperature ofthecanweresimilar.Therefore,thegivenassumptionwas assumedtoholdtrueduringtheinitialperiod.Thesurfaceten- sionvaluealongtheinterfaceofairandwaterwasassigned tobe0.72N/m,andthetimestepsizeusedinallsimulations was1E−4s.

3. Results and discussion

Usingtheresultsofthefirstexperimentaldataset,thesimu- lationschemesweredecided:

- apressurebasedsolverwiththeabsolutevelocityformula- tion,

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor Fig.3–Initialphasecontoursofthegeometrywithair(atthetop)andwaterphases(a)andthemeshstructuresadoptedfor (b)theminimum–199,728and(c)themaximum–337,800sizedmeshes.

- the pressure-velocitycouplingwas carriedout withPISO (pressureimplicit withsplittingofoperator)schemewith skewness-neighbourcoupling,

- transientformulationwasfirstorderimplicit,and

- spatial discretization forgradient was Green-Gaussnode based; for pressure PRESTO; for momentum first order upwind; forvolume fractionCICSAM;and forturbulence kineticenergyfirstorderupwindschemeswereused.Even though the first order upwind scheme is toodissipative to stabilizethe computation, the initialsimulation stud- iesconfirmedthatthegivensolutionschemerssuitbetter for straight convergence and stabilized computation for the chosen timestep sizeand meshresolution. Besides, regardingtheorderofthediscretizationscheme,thesys- temuncertaintyaswellastheturbulencemodeluncertainty mightbelargerthantheerrorcausedbynumericaldissipa- tion.

Usingtheseschemes,thecomputationalmodelwasfirst appliedtostudymeshindependencyandhencetodetermine thefinalmeshconfigurationbasedonthefirstsetofexperi- mentalresults.Then,themodelvalidationstudywascarried outunderareciprocalagitationcondition,andthereciprocally agitatingspeedsfrom20to140rpmwerethentestedfortem- peraturechangeduringtheagitationtodeterminetheeffect ofagitationandoptimumagitationrate.

Fig.3showstheinitialphasecontoursofthegeometrywith headspace–air(atthetop)andwaterphasesandthemesh structuresadoptedfortheminimum(199,728cells)andmax- imum(337,800cells)numberedmeshconfigurations.Fig.4a showstheresultsofmeshindependencystudywithrespect toexperimentaldataobtainedatthecentreofthehorizontal canlocatedinboilingwater.Therewasnotasignificantdiffer- encebetweenthe199,728and286,080cellswhilethe337,800 cellstructuredmeshover-predictedthetemperature(Fig.4a).

Thisdifferencemightbeduetothecontextofmesh-density –round-off errorrelation.Though round-off errorsmay be accumulatedmorewithhighernumberofmeshcells,further

round-offanalysismightberequiredtoidentifytheireffect onthemeshstructureandresolutionfordifferenttime-step sizes.However,itshouldbeemphasizedthattheflowsystem uncertaintyaswellastheturbulencemodeluncertaintyfor

280 290 300 310 320 330 340 350 360

0 30 60 90 120

Temperature (K)

Time (s)

Experimental 199728 cells 286080 337800

(b) (a)

290 300 310 320 330 340 350 360

0 15 30 45 60 75

Temperature(K)

Processme (s)

Experiment Model Medium

Fig.4–Comparisonofexperimentaldatawiththe simulationresults(a)experimentaldataobtainedatthe centreofahorizontallyplacedstaticcaninboilingwater with;(b)experimentaldataobtainedatthecentreofa horizontallyplacedstaticcanunderreciprocallyagitating conditions(themeshstructureusedinbothcomputational modelshad199,728cells).

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor -1.5

-1 -0.5 0 0.5 1 1.5

0 1 2 3 4

Tangenal Velocity (m/s)

TIme (s)

20 rpm 80 rpm 140 rpm

Fig.5–Comparisonoftangentialvelocitychangeversus reciprocalagitationrates.

thechoiceofmeshresolutionandtimestepsizeaccompa- niedcould besignificant onthe flowresultsinaddition to theorderofthespatialdiscretizationschemeandrelaxation parameters.Therecould/maybenostraightforwardsolution forminimumnumericaldiffusionandhigheraccuracywith useofveryhighmeshresolutionaccompaniedwithsmaller timestepsizeandhigherorderspatialdiscretizationscheme.

Themeshindependencysimulationsandtheinitialsimula- tions havedemonstratethat the selectedtimestep sizeof 1E−4sforthegeneratedmeshresolutionof199,728meshcells were accurateenough toobtain resultswhich wouldbe in goodcorrespondencewiththeexperimentaldata.Therefore, basedonthemeshindependencyresults,themeshstructure with199,728cellswasusedinthesecondpartofthemodel

290 300 310 320 330 340 350 360

0 20 40 60 80 100

Tavg (K)

Time (s)

0 rpm 20 40 60 80 140

330 335 340 345 350

40 50 60 70 80

Tavg (K)

Time (s)

60 rpm 80 140

Fig.6–Effectofreciprocalagitationrateonthevolume averageincreaseoftemperature.

Fig.7–Phasecontoursinthevariousx-andz-planesofthe computationalgeometryatthe(a)beginning(1s);(b)30s;

and(c)90softhe20rpmreciprocalagitationcase.

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor Fig.8–Phasecontoursinthevariousx-andz-planesofthe

computationalgeometryatthe(a)beginning(1s);(b)30s;

and(c)90softhe80rpmreciprocalagitationcase.

validationandsimulationstodeterminetheeffectofrecipro- calagitationrate.

Followingthis,themodelwasvalidatedcomparedtothe experimentaldata obtained under different reciprocal agi- tationconditions(summarized inFig.2).Fig. 4b showsthe comparisonofthecancentretemperaturedatawithrespectto

Fig.9–Temperature(K)contoursinthecentralx-and z-planeofthecomputationalgeometryatthe(a)beginning (1s);(b)30s;and(c)90softhe20rpmreciprocalagitation case.

themodelresultsforthefirst75softheprocess.Asobserved inthisfigure,themodelresultsdemonstratedthevalidityof the developedcomputationalmodel. Even thoughthesim- ulation resultscompared well withthe experimentaldata, there wasa differencebetween the simulationresults and experimentaldata.However,consideringthecomplexnature

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor Fig.10–Temperature(K)contourswithinthecentralx-

andz-planeofthecomputationalgeometryatthe(a) beginning(1s);(b)30s;and(c)90softhe80rpmreciprocal agitationcase.

oftheprocessandexperimentalconditions,themodelpre- dictionscaught thetrend oftheexperimentaltemperature data.After75softheprocessing,themodelvalidationcase studydidnotcontinuetorunduetoveryhighrequirementof computationaltime. Ittook 39.3h tocomplete a 1s ofthe

simulationforvalidationpurposeunder80rpmreciprocalagi- tationinanIntelZeon4-Core,3.7GHz–32GBRAMsystem.

3.2. Effectofreciprocalagitationrate

Aftervalidatingthemodel,effectofthereciprocalagitation ratesfrom20to140rpmonthetemperatureevolutioninthe cansweretested.Fig.5showsthecomparisonoftangential velocity values for 20, 80 and 140rpm reciprocal agitation ratesversustimewhileFig.6showstheeffectofreciprocal agitation rates on the volume average temperature (Tavg) increaseofthecanforthefirst100softheprocess.Inthese simulations,theboundarytemperaturewassettobeboiling conditions as reported above in the first model validation case. As demonstrated in Fig. 6, the effect of reciprocal agitation rate was noticeable, and the increased agitation rates increased the temperature evolution especially until 80rpmreciprocalagitation.Theeffectoffurtherincreasing the reciprocal agitation over 80rpm did not result in any significantdifference.Thiseffectofreciprocalagitationuntil amaximumagitationvaluemightbeexplainedbythebalance between the agitation, gravitational buoyancy and viscous forcesgoverningthereciprocalagitationprocess.

Fig.11–Temperature(K)contoursinthecentralx-and z-planeofthecomputationalgeometryatthe(a)beginning (1s);(b)90softhe0rpm(naturalconvectionheatingcase).

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor AsindicatedbyBoonpongmaeeandMakotani(2009)and

TutarandErdogdu(2012),theflowfield duringanagitation processinvolves complexinteractionsofcentrifugal–rota- tional,gravitational buoyancyand viscousforces.Basedon thisconcept,the followingforceanalysiswasperformedto betterunderstandthetemperatureevolutionduringtherecip- rocalagitationprocessasafunctionofgravitational,agitation andviscousforces:

f_cbf f_gbf =ω²r

g (15)

f_inertial f_viscous= ω²r⁴

² (16)

whererwasthecrankradiusofthesystem(0.075m),gwas thegravitationalacceleration(9.81m/s²),andwaskinematic viscosity(m²/s), f_cbf,f_gbf,f_inertial andf_viscouswere centrifugal buoyancy,gravitationalbuoyancy,agitationrelatedforcesand viscousforces,respectively,andωwasthedimensionalspeed ofrotation(1/s).ωwasshownby(ω=2f/60)wherefwasthe horizontalagitationrate(rpm).Eqs.(15)and(16)represented FroudeandTaylornumbers,respectively.

The(ω²r) value,infact, showed the effect ofhorizontal agitationrateoverthegravitationalacceleration,anditwas definedtobethereciprocationintensity(g0):

g₀=ω²·r·

1+ r

L (17)

whereLwasthelengthoftheconnectivityrod(0.5m)of theslidercranktypesystemusedintheexperimentalstudies.

Considering that [(r/L=0.15)<1],the reciprocation intensity wasapproximatedby:

g₀=ω²·r (18)

This was used in the definition of Froude number to demonstratethereciprocationeffectovergravitationalforce.

Froudenumberincreasedfrom0.03to1.64withtheincrease ofthehorizontalagitationratefrom20to140rpm.Simulta- neously,Taylornumberincreasedfrom2.26E8to1.18E10,from 3.65E9to2.29E10andfrom1.13E10to6.87E10atthe100sof 20,80and140rpmshaking,respectively.Theincreaseratio ofTaylornumberwas5.22,6.27 and6.08(theseratios were betweentheTaylor numbervaluesobtainedatthe100and 1softheprocess)forthese3-reciprocalagitationrates.This change in the Taylor number indicated that the reciprocal agitationforcesstartedshowingtheirimpactevenatlowagi- tationrates,butthiseffectincreasedtoacertainhighestvalue at80rpm.Theeffectofinertialforcesobtainedbythehori- zontalagitationoverviscousforcesasaninternalresistance oftheprocessedliquidwasnotsignificantbeyond80rpm.The furtherincreaseafter80rpmagitationratedidnotmakeany significanteffectontheTaylornumberincreaseinthegiven process,anddidnotleadtoanyfurthertemperatureincrease.

In fact, the volumeaverage temperature increaseobtained with80rpmagitationwasoverthecaseof140rpmtowards theendoftheprocess(Fig.6)asalsoindicatedbytheincrease

Fig.12–Instantaneousvelocityvectors(m/s)onx–zandy–zmid-planesoftheflowdomainat(a)1s(verybeginningofthe agitation)and(b)90s(towardstheendoftheagitation)at20rpm.

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor ratioinTaylornumber.Froudenumberwas0.54at80rpm,and

thisvaluewas50%ofthegravitationalaccelerationeffect.

AsindicatedbyWaldenandEmanuel(2010),thereciprocal agitationapparentlyusedthehorizontalaccelerationinaddi- tiontogravityforces,andthesumoftheseforcesenableda considerabletemperatureincreasewhilethefurtherincrease after80rpmwas preventedviatheeffect ofviscousforces.

Eventhoughthereciprocationintensitywas1.64timeshigher thantheregularnaturalconvectioncasegovernedbythegrav- itationalaccelerationat140rpm,thetemperatureincreasedid notgain anymorebenefitfrom this highrate ofagitation.

Thelowerviscosityoftheliquid(water)alsoplayedasignifi- cantroleinthisconcept,anditseffectisexpectedtobemore pronouncedforthecaseofnon-Newtonianhigherviscosity liquids.Trevino(2009)alsonotedthattheoscillatingtechnol- ogywithrackingmovement(reciprocallyagitation)mightnot bebeneficialforprocessinglow(1%starch–watermixture)and high (5% starch–watermixture)viscosity productswhere a comparisonoftheeffectsofoscillatingandstaticretortther- malprocessingwascarriedout.However,itwasalsostated thataneffectiveheattransferratemightbepossibletoobtain formediumviscosity(3%starch–watermixture)products.

Tutar and Erdogdu (2012) explained that, in agitation relatedprocesses,headspace–airbubblemovesthroughvia thegiveneffectofagitationandviscosityforces.Thisleadsto themixingtoincreasetheheattransferrate.Forthecaseof lowviscosityNewtonianliquids,however,thismixingeffectis

hardlyseen,andtheairbubblegenerallymightmovethrough thetopovertheprocessunderacertaineffectofagitationrate.

Figs.7and8showthephase(Phase2–headspace)contours (headspaceshownwithredcontouratthebeginningofthe process)atthebeginning(1s),30and90softhe20and80rpm cases inthe various x- and z-planesof the computational geometry,respectively.AsobservedinFig.7,headspacemoved atthetopofthegeometrycontinuouslywiththegivenmove- mentofthecanat20rpmreciprocalagitationrate.However, the80rpmagitationledtoanabruptchangeoftheheadspace distributionduetothesuddenstartofthehigheragitation.

Theconsiderabledifferencebetween20and80rpmagitation rateswasalsoshowninFig.5.Fig.8ashowsthesuddendisrup- tionoftheheadspaceatthebeginningcomparedtothecase of20rpmandacertaininhomogeneousdistributionthrough theprocess(Figs.8b,c).Atthe90softheprocess,basedon thephasecontours,itmightbeassumedthatalowamount ofheadspacewasmixedinthewaterat20rpmagitationrate whilethiswas higherat80rpmagitation rate.Thismixing alsobroughtaconsiderableandhomogeneoustemperature increasecomparedtothecaseof20rpm.

Figs.9and10showthetemperaturedistributionthrough the canfor20and 80rpm,respectively.Themorehomoge- neous distribution natureof the temperature contoursare observedat80rpmagitationratewhilethe20rpmagitation rate resembled more like a natural convection case with the distinctly stratified temperature contours. The natural

Fig.13–Instantaneousvelocityvectors(m/s)onx–zandy–zmid-planesoftheflowdomainat(a)1s(verybeginningofthe agitation)and(b)90s(towardstheendoftheagitation)at80rpm.

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor convectioncasewassummarizedinFig.11.Figs.12and13,on

theotherhand,comparativelyrepresentedtheinstantaneous velocityvectorsobtainedonthexzandyzmid-planesofthe flow domain at1 (very beginning of the agitation motion) and90s(towardstheendoftheagitationmotion)at20and 80rpm,respectivelytofurthercomparetheagitationeffects forthefieldofflowmotion.Whentheagitationstarted,the headspacewasimmediatelyaffected,anditcouldnotmain- tainitsshapeandorientationwithrespecttothewaterflow domain,whichwasacceleratedlikeasolidbodymotiondue tothehorizontalaccelerationasseeninFigs.12aand13a.As theagitationmotioncontinued,dynamicpressuredeveloped intheflowdomain,anddifferentspatialandtemporalevolu- tionofthefluidflow,dependinguponthesimulationtimeand agitationrate,wereclearlyidentifiedfromtheinstantaneous velocity vectors in the vicinity of the air–water interface whoseshapeandpositionbecamehighlyunstablewiththe increasingagitation rate (Figs. 12b and 13b). Superposition ofinertial,agitation andnaturalforcesduetogravitational forces(mixingconvectiveforces)became moreevidentand effectiveontheair–waterflowdomainathigheragitationrate of80rpm,leadingtohigherinertialandconvectiveinstabili- tiesofthepresenttwo-phaseflowsystemwithmorecomplex, three-dimensionalflowbehaviourinthewholedomainsys- tem.Thisbehavioureventuallywouldmakeapositiveeffect onthetemperatureevolutionfor80rpmcomparedto20rpm aspreviouslyobservedinthetemperaturecontoursinFig.10.

With this positive effect, forced convective heat transfer mechanism on the temperature evolution became more dominantastheagitationrateincreasedupto80rpmanddid notchangeathigheragitationratesasalsoexplainedabove.

4. Conclusions

This study introduced determining the optimal reciprocal agitation rate based on the experimentally validated com- putational model. The computational model was used to determinethetemperaturechangeinacanfilledwithwater with 2% headspace undergoing a reciprocal agitation pro- cess,andthetemperatureevolutionduringtheprocesswas determinedtobeundercontrolofreciprocationintensityand theratioofagitationandviscousforces.Thecomputational resultsindicateda certainlimitofreciprocal agitationrate intheviewoftemperatureincreaseindicatingthesignificant effectofviscosityand inertialforcesobtainedbytheagita- tion.ForaNewtonianlowviscosityliquidcase,represented bywater,the80rpmreciprocalagitationratewasdetermined tobeanoptimumrate.

Itwouldbevaluabletodeterminetheoptimumagitation conditionsforhighviscositynon-Newtonianliquidsconsid- eringthatasignificantportionofthefoodproductsprocessed incanslieinthiscategory.Besides,afurtherstudytodemon- stratetheeffectofheadspacevolumetoincreasetheagitation ratesforliquidandparticulatefoodproductswouldalsobe required.

Acknowledgement

This study was developed within the framework of the ERA-net SUSFOOD Sunniva project, “Sustainable food pro- ductionthroughqualityoptimizedrawmaterialproduction andprocessingtechnologies forpremiumqualityvegetable productsandgeneratedby-products”.Thefollowingacknowl- edgementsarerecognizedbytheauthors:

FerruhErdogduacknowledgestheMinistryofFood,Agri- cultureandLivestock(GDAR)ofTurkeyforthetravelsupport toattendtheprojectmeetings.

MustafaTutaracknowledgesthesupportfromtheDaniel andNinaCarassoFoundation,France(www.fondationcarasso.

org)andtheDepartmentofEconomicDevelopmentandCom- petitiveness(ELIKA).

SigurdOinesandDagbjornSkipnesacknowledgethesup- portfromtheResearchCouncilofNorwaythroughgrantno.

NO10829.

References

Ates,M.B.,Skipnes,D.,Rode,T.M.,Lekang,O-I.,2014.Comparison ofbacterialinactivationwithnovelagitatingretortandstatic retortaftermildtreatments.FoodControl43,150–154.

Bermudez-Aguirre,D.,Lima,F.,Reitzel,J.,Garcia-Prez,M., Barbosa-Canovas,G.V.,2013a.Evaluationoftotalheattransfer coefficient(hT)duringinnovativeretortprocessing:static, gentlemotionandrockingmode.In:IFTAnnualMeeting, Abstractnumber:031-02.

Bermudez-Aguirre,D.,Lima,F.,Reitzel,J.,Garcia-Prez,M., Barbosa-Canovas,G.V.,2013b.Developmentofa

mathematicalmodeltodescribeheattransferusingthree differentprocessingmodes.In:IFTAnnualMeeting,Abstract Number:031-138.

Boonpongmaee,T.,Makotani,Y.,2009.Heattransferinrotating cylindricalcellswithpartitions.Int.J.HeatFluidFlow30, 211–217.

Clifcorn,I.E.,Peterson,G.T.,BoydJ.m.O‘Neil,J.H.,1950.Anew principleforagitatinginprocessedcannedfoods.Food Technol.4,450–460.

Erdogdu,F.,Tutar,M.,2012.Acomputationalstudyforaxial rotationeffectsonheattransferinrotatingcanscontaining liquidwater,semi-fluidsystemandheadspace.Int.J.Heat MassTransf.55,3774–3788.

Hirt,C.W.,Nichols,B.D.,1981.Acomputationalmethodfor pressure0surfacehydrodynamics.J.Press.VesselTechnol.– Trans.ASME103,136–141.

Julien,K.,Legg,S.,McWilliams,J.,Werne,J.,1996.Hardturbulence inrotationRayleigh-Benardconvection.Phys.Rev.E53,R5557.

Kooij,G.L.,Botchev,M.A.,Geurts,B.J.,2015.Directnumerical simulationofNusseltnumberscalinginrotating

Rayligh-Benardconvection.Int.J.HeatFluidFlow55,26–33.

Liffman,K.,Metcalfe,G.,Cleary,P.,1997.Convectiondueto horizontalshaking.In:CSIRO-1997:InternationalConference onCFDinMineralandMetalProcessingandPower

Generation,pp.165–168.

Ohlsson,T.,1980.Optimalsterilizationtemperaturesforsensory qualityincylindricalcontainers.J.FoodSci.45,1517–1521.

Pesch,W.,Palaniappan,D.,Tao,J.,Busse,F.H.,2008.Convectionin heatedfluidlayerssubjectedtotime-periodichorizontal accelerations.J.FluidMech.596,313–332.

Reader,G.T.,Hooper,L.C.,1982.StirlingEngines.SponPress, Oxfordshire,UK.

Rosnes,J.T.,Skara,T.,Skipnes,D.,2011.Recentadvancesin minimalheatprocessingoffish:effectsonmicrobialactivity andsafety.FoodBioprocessTechnol.4,833–848.

Singh,A.P.,Singh,A.,Ramaswamy,H.S.,2015a.Modificationofa staticsteamretortforevaluationheattransferunder reciprocationagitationthermalprocessing.J.FoodEng.153, 63–72.

Singh,A.,Singh,A.P.,Ramaswamy,H.S.,2015b.Arefined methodologyforevaluationofheattransfercoefficientsin cannedparticulatefluidsunderrapidheatingconditions.

FoodBioprod.Process.94,169–179.

Singh,A.,Ramaswamy,H.S.,2015a.Effectofproduct-related parametersonheat-transferratestocannedparticulate non-Newtonianfluids(CMC)duringreciprocatingagitation thermalprocessing.J.FoodEng.165,1–12.

Pleasecitethisarticleinpressas:Erdogdu,F.,etal., Determiningtheoptimalshakingrateofareciprocalagitationsterilizationsystemfor Singh,A.P.,Ramaswamy,H.S.,2015b.EffectofCanorientationon

heattransfercoefficientsassociatedwithliquidparticulate mixturesduringreciprocationagitationthermalprocessing.

FoodBioprocessTechnol.8,1405–1418.

Singh,A.,Singh,A.P.,Ramaswamy,H.S.,2016.Acontrolled agitationprocessforimprovingqualityofcannedbeans duringagitationthermalprocessing.J.FoodSci.81, E1399–E1411.

Singh,A.P.,Ramaswamy,H.S.,2016.Simultaneousoptimization ofheattransferandreciprocationintensityforthermal processingofliquidparticulatemixturesundergoing reciprocatingagitation.Innov.FoodSci.Emerg.Technol.33, 405–415.

Trevino,J.,(masterofsciencethesis)2009.EffectofOscillating andStaticRetortThermalProcessingTechnologyUsingand InstitutionalSizePouch.ClemsonUniversity,Clemson,SC, USA.

Tutar,M.,Erdogdu,F.,2012.Numericalsimulationforheat transferandvelocityfieldcharacteristicsoftwo-phaseflow systemsinaxiallyrotatinghorizontalcans.J.FoodEng.111, 366–385.

Ubbink,O.,Issaa,R.J.,1999.Amethodforcapturingsharpfluid interfacesonarbitrarymeshes.J.Comput.Phys.153, 26–50.

Walden,R.,Emanuel,J.,2010.Developmentsinin-container retorttechnology:theZinetecShakaprocess.In:Doona,C.J., Kustin,K.,Feery,F.E.(Eds.),CaseStudiesinNovelFood ProcessingTechnologies:InnovationinProcessing,Packaging andPredictiveModelling.WoodheadPublishingLtd., Cambridge,UK.

Yakhot,H.Q.,Orszag,S.A.,1986.Renormalizationgroupanalysis ofturbulence.I.Basictheory.J.Sci.Comput.1,1–51.