- Email: [email protected]
A Stochastic Spatio-Temporal Model of the Flow-Front Dynamics in a Vacuum Assisted Resin Transfer Moulding Process
Proceedings,18th IFAC Symposium on System Identification Proceedings,18th IFAC Symposium on System Identification July 9-11, 2018. Stockholm, Sweden on Proceedings,18th IFAC Symposium Symposium on System System Identification Identification Proceedings,18th IFAC online at www.sciencedirect.com July 9-11, 2018. Stockholm, Sweden Available Proceedings,18th IFAC Symposium July Sweden July 9-11, 9-11, 2018. 2018. Stockholm, Stockholm, Sweden on System Identification July 9-11, 2018. Stockholm, Sweden
ScienceDirect
IFAC PapersOnLine 51-15 (2018) 383–388
A Stochastic Spatio-Temporal Model of the A Stochastic Spatio-Temporal Model of the A Stochastic Spatio-Temporal Model of the Flow-Front Dynamics in a Vacuum Assisted A Stochastic Spatio-Temporal Model of the Flow-Front Dynamics in a Vacuum Assisted Flow-Front Dynamics in a Vacuum Assisted Resin Transfer Moulding Process Flow-Front Dynamics in a Vacuum Assisted Resin Transfer Moulding Process Resin Transfer Moulding Process Resin Transfer Moulding Process ∗∗ Michael Nauheimer ∗∗ Rishi Relan ∗∗ ∗∗ Uffe Høgsbro Thygesen ∗∗
Michael Nauheimer Relan Thygesen ∗ ∗∗ ∗∗ ∗∗ ∗ Høgsbro ∗ ∗ Rishi ∗∗ Uffe ∗∗ Michael Nauheimer Rishi Relan Uffe Høgsbro Thygesen Henrik Madsen Olesen Kirkeby Michael Nauheimer Relan Høgsbro Thygesen ∗∗ Bendt ∗ Klaus ∗ ∗ Rishi ∗∗ Uffe ∗∗ Henrik Madsen Olesen Kirkeby ∗∗ Bendt ∗ Klaus ∗ Michael Nauheimer Rishi Relan Uffe Høgsbro Thygesen ∗∗ ∗ ∗ Henrik Madsen Madsen Bendt Bendt Olesen Olesen Klaus Klaus Kirkeby Kirkeby Henrik ∗∗ ∗ ∗ Henrik Madsen Bendt Olesen Klaus Kirkeby ∗ ∗ Siemens Gamesa Renewable Energy, Denmark ∗ Siemens Gamesa Renewable Energy, Denmark ∗ Siemens Gamesa Denmark (e-mail: [email protected], [email protected], Gamesa Renewable Renewable Energy, Energy, Denmark ∗ Siemens (e-mail: [email protected], [email protected], Siemens Gamesa Renewable Energy, Denmark (e-mail: [email protected], [email protected], [email protected]) (e-mail: [email protected], [email protected], [email protected]) ∗∗ (e-mail: [email protected], [email protected], [email protected]) Technical University of Denmark, Denmark [email protected]) ∗∗ DTU Compute, Technical University Denmark, ∗∗ DTU Compute, [email protected]) ∗∗ DTU Compute, Technical [email protected], University of of [email protected]) Denmark, Denmark Denmark (e-mail: [email protected], DTU Compute, Technical University of Denmark, Denmark ∗∗ (e-mail: [email protected], DTU Compute, Technical [email protected], University of [email protected]) Denmark, Denmark (e-mail: [email protected], [email protected], [email protected]) (e-mail: [email protected], [email protected], [email protected]) (e-mail: [email protected], [email protected], [email protected]) Abstract: With With more more emphasis emphasis on on the the use use of of green green energy, energy, the the size size of of the the turbines turbines and and blades blades Abstract: Abstract: With more emphasis on the use of green energy, the size of the turbines and blades in the wind turbines is continuously increasing. With increasing blade size, the casting process Abstract: With more emphasis on the use of green energy, the size of the turbines and blades in the wind With turbines is continuously increasing. With increasing blade the casting process Abstract: more on risk the use of green energy, the size of size, the turbines andmade blades in the turbines is continuously increasing. With increasing blade size, the casting process becomes more complicated and the the of faults faults increases. Production of such such blades, of in the wind wind turbines is emphasis continuously increasing. With increasing blade size, the blades, casting process becomes more complicated and risk of increases. Production of made of in the wind turbines is continuously increasing. With increasing blade size, the casting process becomes more complicated complicated and the the risk risk of faults faults increases. Production of such such blades, blades, made of fibre reinforced polymer composites, without the possibility of visual inspection of the infusion becomes more and of increases. Production of made of fibre reinforced polymer composites, without the possibility of visual inspection of the infusion becomes morefor complicated and the(possibly risk of faults Production ofprocess. such blades, of fibre reinforced polymer without the possibility of inspection of the infusion process calls a sensor sensorcomposites, system virtual) for monitoring monitoring the Therefore, in fibre reinforced polymer composites, without the increases. possibility of visual visualthe inspection ofTherefore, the made infusion process calls for a system (possibly virtual) for process. in fibre reinforced polymer composites, without the possibility of visual inspection of the infusion process calls for a sensor system (possibly virtual) for monitoring the process. Therefore, in this paper, we propose a stochastic spatio-temporal modelling approach to model the flow-front process calls for a sensor system (possibly virtual)modelling for monitoring the to process. Therefore, in this paper, we propose a spatio-temporal approach model flow-front process calls a sensor system (possibly for monitoring the to process. Therefore, this paper, wefor propose a stochastic stochastic spatio-temporal modelling approach to model the the flow-front dynamics inside the vacuum assisted resin virtual) transfermodelling moulding process, knowledge offlow-front which in is this paper, we propose a stochastic spatio-temporal approach model the dynamics inside the vacuum assisted resin transfer moulding process, knowledge of which this paper, we propose a stochastic spatio-temporal modelling approach model theof dynamics inside the vacuum vacuum assisted resinoftransfer transfer moulding process, to knowledge offlow-front which is is essential for determining the current state the infusion process. dynamics inside the assisted resin moulding process, knowledge which is essential for determining the current the infusion process. dynamics inside the vacuum assistedstate resinof moulding process, knowledge of which is essential determining the state of the process. essential for for determining the current current state oftransfer the infusion infusion process. © 2018, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. essential for determining the current state of the infusion process. Keywords: Wind Wind turbines, turbines, Partial Partial differential differential equations, equations, Stochastic Stochastic differential differential equations, equations, Keywords: Keywords: Likelihood, Wind turbines, turbines, PartialKalman differential equations, Stochastic differential equations, equations, Maximum Extended filter, Grey-box models Keywords: Wind Partial differential equations, Stochastic differential filter, Grey-box models Maximum Extended Keywords: Wind turbines, PartialKalman differential Maximum Likelihood, Likelihood, Extended Kalman filter, Grey-box models Maximum Likelihood, Extended Kalman filter,equations, Grey-boxStochastic models differential equations, Maximum Likelihood, Extended Kalman filter, Grey-box models 1. INTRODUCTION INTRODUCTION (2003)). In general these sensors are accurate but limited 1. (2003)). In these sensors are accurate but limited 1. INTRODUCTION INTRODUCTION (2003)). In general general these sensors are accurate butmoulded limited to measuring on or close to the surface of the 1. (2003)). In general these sensors are accurate but limited to measuring on or close to the surface of the moulded 1. INTRODUCTION (2003)). In general these sensors are accurate but limited to measuring measuring on or or closeinto tothethe the surface of the the based moulded parts. It has been shown past that sensors on Large scale composite shell structures like wind turbine to on close surface of moulded parts. It has been shown in the past that sensors based on Large scale composite shell structures like wind turbine to measuring on or close to the surface of the moulded parts. It has been shown in the past that sensors based on Large scale composite shell structures like wind turbine optical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer parts. It has been shown in the past that sensors based on Large scale composite shell structures like wind turbine optical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer parts. It has been shown in the past that sensors based on Large scale composite shell structures like wind turbine optical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer making it possible to measure through thickness of the structure. These structures are usually manufactured usoptical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer making it possible to measure through thickness of the structure. These structures are usually manufactured usoptical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer making it possible to measure through thickness of the structure. These structures are usually manufactured usshell. Recently, (Matsuzaki and Shiota (2016)) proposed ing the vacuum assisted resin transfer moulding (VARTM) making it possible to measure through thickness of the structure. Theseassisted structures are usually manufactured us- shell. Recently, (Matsuzaki and Shiota (2016)) proposed ing the vacuum resin transfer moulding (VARTM) it using possible to measure thickness ofdata the structure. These structures usually manufactured us- making Recently, (Matsuzaki and through Shiota (2016)) (2016)) proposed ing the vacuum vacuum assisted resin transfer moulding (VARTM) ashell. method two-sided visual observations for 4D process. The VARTM process involves infusing a therthershell. Recently, (Matsuzaki and Shiota proposed ing the assisted resinare transfer moulding (VARTM) ashell. method using two-sided visual observations for 4D data process. The VARTM process involves infusing a Recently, (Matsuzaki and Shiota (2016)) proposed ing the vacuum assisted resin transfer moulding (VARTM) a method using two-sided visual observations for 4D data process. The VARTM process involves infusing a therto accurately reconstruct the 3D flow moset polymer polymer into a fibre fibre reinforcement preform.aEpoxy Epoxy method using two-sided visual observations forresin 4D data process. The VARTM process involves infusing ther- aassimilation to accurately reconstruct the 3D flow moset into reinforcement preform. aand method using two-sided observations forresin 4D data process. The VARTM process involves infusing aEpoxy ther- assimilation assimilation to accurately reconstruct the 3D resin flow moset polymer into aaatofibre fibre reinforcement preform. Epoxy permeability filed of aavisual fibre preform. resin is transferred the mould through several inlets. assimilation to accurately reconstruct the 3D resin flow moset polymer into reinforcement preform. and permeability filed of fibre preform. resin is transferred to the mould through several inlets. assimilation to accurately reconstruct the 3D resin flow moset polymer into a fibre reinforcement preform. Epoxy and permeability filed of a fibre preform. resin is transferred to the mould through several inlets. Progression of the epoxy inside the mould is only driven by and permeability filed of a fibre preform. resin is transferred to the mould through several inlets. However, some blades produced at the Siemens Gamesa Progression of the epoxy inside the mould is only driven by and permeability filed of a fibre preform. resin is transferred to the mould through several inlets. blades produced at Gamesa Progression of the epoxy epoxy inside the the mould mould is is only driven by However, the pressure pressure of difference between and thedriven pressure Progression the inside only by However, some some blades produced at the theinSiemens Siemens Gamesa Renewable Energy factories are casted one piece using the difference between the mould mould is and the pressure some blades produced at the Siemens Gamesa Progression of the epoxy inside the only driven by However, Renewable Energy factories are casted in one piece using the pressure difference between and the pressure of one bar surrounding the mould in the production area. the pressure difference between the mould and the pressure R However, some blades produced at the Siemens Gamesa Renewable Energy factories are casted in one piece using the patented IntegralBlades technology (Stiesdal et al. of one bar surrounding the mould in the production area. Renewable Energy factories are casted in one piece using R the pressure difference between the mould and the pressure the patented IntegralBlades technology (Stiesdal et al. of one bar surrounding the mould in the production area. R of onetobar surrounding the mould inofthe production area. Renewable Energy factories are casted in one piece using R technology the patented IntegralBlades (Stiesdal et al. (2006)) instead of two separate half blades which are Due the inhomogeneous nature the flow inside the the patented IntegralBlades technology (Stiesdal et al. of onetobar surrounding the mould inofthe production area. R (2006)) instead of two separate half blades which are Due the inhomogeneous nature the flow inside the the patented IntegralBlades technology (Stiesdal et al. (2006)) instead of two separate half blades which are Due to the inhomogeneous nature of the flow inside the then united afterwards. One of the main disadvantages mould, risk of moulding defects such as dry spots (2006)) instead of two separate half blades which are Due to the inhomogeneous nature of the flow inside the then united afterwards. One of the main disadvantages mould, the inhomogeneous risk of moulding defects such as dry spots (2006)) instead of two separate half blades which are Due to the nature of the flow inside the then united afterwards. One of the main disadvantages mould, risk of moulding defects such as dry spots of this technology is that, it is not possible to visually and voids increases and this leads to deterioration of the then united afterwards. One of the main disadvantages mould, theincreases risk of and moulding defects such as dry ofspots of technology is it is possible to and voids this leads to deterioration the afterwards. One the main disadvantages mould, theincreases risk of and moulding defects such as dry of spots of this thisunited technology isof that, that, it of is not not possible to visually visually and voids increases and this leads to deterioration deterioration ofet the inspect the process transferring the epoxy resin to the mechanical properties of this the leads cast parts (Sreekumar al. then of this technology is that, it is not possible to visually and voids to the inspect the process of transferring the epoxy resin to the mechanical properties of the cast parts (Sreekumar et al. this technology is that, it is not possible to Gamesa visually and voids increases and leads toetdeterioration ofet inspect the process of transferring the epoxy resin to the mechanical properties of this the cast parts parts (Sreekumar et the al. of mould. Hence, the main motivation at Siemens (2007); Park et al. (2011); Matuzaki al. (2015)). Hence, inspect the process of transferring the epoxy resin to the mechanical properties of the cast (Sreekumar al. mould. Hence, the main motivation at Siemens Gamesa (2007); Park et al. (2011); Matuzaki et al. (2015)). Hence, inspect the process of transferring the epoxy resin to the mechanical properties of the cast parts (Sreekumar et al. mould. Hence, the main motivation at Siemens Gamesa (2007); Park et al. (2011); Matuzaki et al. (2015)). Hence, Renewable Energy is to develop a flow-front monitoring to prevent any potential fault from developing and for mould. Hence, the main motivation at Siemens Gamesa (2007); Park et al. (2011); Matuzaki et al. (2015)). Hence, Renewable Energy is to develop aa flow-front monitoring to prevent any potential fault from developing and for mould. Hence, the main motivation at Siemens Gamesa (2007); Park et al. (2011); Matuzaki et al. (2015)). Hence, Renewable Energy is to develop flow-front monitoring to prevent any potential fault from developing and for system which eventually aids in controlling the process better control of the process, it is necessary to monitor the Renewable Energy is to develop a flow-front monitoring to prevent any potential fault developing and the for system which eventually aids in controlling the process better control of the process, it isfrom necessary to monitor Energy is to develop flow-front monitoring to prevent fault developing and the for Renewable system which eventually aids trajectory in acontrolling controlling the process better control of potential the process, process, it necessary to monitor of mould filling. A controlled of the flow-front evolution ofany the flow-front inside the mould mould inmonitor a VARTM VARTM system which eventually aids in the process better control of the it is isfrom necessary toin the of mould filling. A controlled trajectory of the flow-front evolution of the flow-front inside the a system which eventually aids in controlling the process better control of the process, it is necessary to monitor the of mould filling. A controlled trajectory of the flow-front evolution of the the flow-front flow-front inside inside the the mould mould in in aa VARTM VARTM of inside the mould decreases the risk of areas with dry glass, process. mould filling. A controlled trajectory of the flow-front evolution of inside the mould decreases the trajectory risk of areas with dry glass, process. of mould filling. A controlled of the flow-front evolution of the flow-front inside the mould in a VARTM inside the mould risk of areas with dry glass, process. which eventually decreases the repair time and increase the inside the mould risk of areas with dry glass, process. which eventually decreases the repair time and increase the Sensors for monitoring the VARTM processes is an acinside the mould decreases the risk of areas with dry glass, process. which eventually eventually repair time and and increase increase the general quality of the produced blade. Sensors for monitoring the VARTM processes is an acwhich decreases the repair time the general quality of the produced blade. Sensors for monitoring the VARTM VARTM processes is an preactive field of research. Several types of sensors have Sensors for monitoring the processes is an acwhich eventually decreases the repair time and increase the general quality quality of of the the produced produced blade. blade. tive field of research. Several types of sensors have pregeneral Sensors for thein VARTM is an preac- To tive field of monitoring research. Several types of ofprocesses sensors have have preachieve this goal, this paper proposes a spatio-temporal viously been implemented VARTM processes including tive field of research. Several types sensors general quality of the produced blade. this this paper aa spatio-temporal viously been implemented in processes including tive field of research. Several types sensors have pre- To To achieve achieve this goal, goal,the thisdynamics paper proposes proposes spatio-temporal viously been implemented in VARTM VARTM processes including model to estimate of the inside permittivity sensors (Yenilmez (Yenilmez and of Sozer (2009)), presTo achieve this goal, this paper proposes aflow-front spatio-temporal viously been implemented in VARTM processes including model to estimate the dynamics of the inside permittivity sensors and Sozer (2009)), presTo achieve this goal, this paper proposes aflow-front spatio-temporal viously been implemented in VARTM processes including model to estimate estimate the the dynamics of the thewe flow-front inside permittivity sensors et (Yenilmez andand Sozer (2009)), presthe mould. To emulate real system, first model the sure sensors (Zhang al. (2011)), sensors based on model to the dynamics of flow-front inside permittivity sensors (Yenilmez and Sozer (2009)), presthe mould. To emulate the real system, we first model the sure sensors (Zhang et al. (2011)), and sensors based on to estimate the the dynamics of thewe flow-front inside permittivity sensors et (Yenilmez andand Sozer (2009)), presthe mould. To emulate real system, first model the sure sensors (Zhang et al. (2011)), (2011)), and sensors based on model spatio-temporal evolution of the flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. the mould. To emulate the real system, we first model the sure sensors (Zhang al. sensors based on spatio-temporal evolution of the flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. the mould. To emulate the real system, we first model the sure sensors (Zhang et al. (2011)), and sensors based on spatio-temporal evolution of the flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. partial differential equations (PDEs) based on Darcy’s evolution of the flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. spatio-temporal partial differential equations based on Darcy’s spatio-temporal evolution of the(PDEs) flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. partial differential equations (PDEs) based on Darcy’s Law (Whitaker (1986)). The PDE based model is then We thank the Manufacturing Academy of Denmark (MADE) for partial differential equations (PDEs) based on Darcy’s We thank the Manufacturing Academy of Denmark (MADE) for Law (1986)). The PDE based model is We thank partial differential equations based onevolution Darcy’s Law (Whitaker (Whitaker (1986)). The (PDEs) PDE based model is then then providing thethe financial support. We thank of theDenmark laboratory employees the Manufacturing Academy of (MADE) for simulated to generate different patterns of the Law (Whitaker (1986)). The PDE based model is then We thank Manufacturing Academy (MADE) for providing the financial support. We thank theDenmark laboratory employees simulated to generate different patterns of the evolution Law (Whitaker (1986)). The PDE basedof model is then We thank Manufacturing Academy (MADE) for at Siemens Gamesa Renewable Energy forDenmark their assistance and providing the financial support. thank the laboratory employees simulated to generate different patterns the evolution of the flow-front inside the mould. providing thethe financial support. We We thank of the employees simulated to generate different patterns of the evolution at Siemens Gamesa Renewable Energy for laboratory their assistance and of the flow-front inside the mould. providing theGamesa financialRenewable support. We thank the employees suggestions. at Siemens Energy for their and simulated to generate patterns of the evolution at Siemens Gamesa Renewable Energy for laboratory their assistance assistance and of inside the of the the flow-front flow-front inside different the mould. mould. suggestions. at Siemens Gamesa Renewable Energy for their assistance and suggestions. of the flow-front inside the mould. suggestions.
suggestions. 2405-8963 © © 2018 2018, IFAC IFAC (International Federation of Automatic Control) Copyright 383 Hosting by Elsevier Ltd. All rights reserved. Copyright © under 2018 IFAC 383 Control. Peer review responsibility of International Federation of Automatic Copyright © 2018 IFAC 383 Copyright © 2018 IFAC 383 10.1016/j.ifacol.2018.09.175 Copyright © 2018 IFAC 383
2018 IFAC SYSID 384 9-11, 2018. Stockholm, Sweden July
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
A PDE based model is generally not a preferred choice for online monitoring and control. Hence, the main contribution of this paper is to develop a spatio-temporal model of the flow-front dynamics using coupled stochastic differential equations (SDEs) (Øksendal (2010)). The validation of the coupled SDEs based model is performed for different patterns of flow-front dynamics generated using the PDE based model. Finally, an analysis of the level of measurement noise, sampling rate and the model complexity on the final estimated model is also presented. The paper is structured as follows: Section 2 briefly describes the VARTM process. In Section 3, the mathematical background and implementation of PDEs based flowfront simulation model is discussed. A model for estimation of the flow-front based on coupled SDEs is proposed in Section 4. Section 5 describes the results and finally the conclusions are given in Section 6. 2. VARTM PROCESS In this section, we briefly describe the VARTM process for blade manufacturing. Fig. 1 shows the work-flow of the VARTM process. Flow mesh
Vacuum foil
Resin inlet
Release film
Fibre reinforcement Vacuum hose
Aluminum plate
3. SIMULATING FLOW IN A POROUS MEDIUM To deal with the above mentioned challenges in practice, first it is important to completely understand the dynamics of the flow-front progression inside the mould. For that purpose, we first propose to simulate the flow-front progression in a VARTM process in a 80 cm × 90 cm rectangular casting for our simulation study. We have implemented Darcy’s law with relevant boundary conditions into a PDE solver software, FEniCS Project (Alnæs et al. (2015)), and use the output to create flow-front pattern data representative of the data which can be measured experimentally using commercially available sensors and instrumentation (Nauheimer et al. (2017)). In the section below, we briefly explain the mathematical background behind the flow-front modelling. 3.1 Flow-front simulation model Assumption 1: Epoxy is considered to be a Newtonian fluid. Remark 1: Uncured epoxy is not a Newtonian fluid. It is assumed that minimal or equal distributed physical force will be exerted on the epoxy within the mould and thereby it is assumed that it can be considered a Newtonian fluid for our application. Assumption 2: The increase in viscosity is small and constant enough to be omitted for shorter time scales. Remark 2: The curing of epoxy initiates as soon as it has been mixed therefore, as the epoxy cures, the viscosity increases. It is assumed that for the short time frame of our simulations, the increase in viscosity is small enough to be omitted.
Sealant
Fig. 1. The general material configuration and work-flow of a VARTM process. The resin enters the vacuum bag through the resin inlet and spreads into the pores of the fibre reinforcement from left to right. This process is only assisted by the vacuum generated through the vacuum hose. The VARTM process is a special variety of Resin Transfer Moulding (RTM) processes. As the name indicates the process is only assisted by the difference of the nearvacuum inside the mould and the surrounding pressure of approximately one bar. The flow of epoxy inside the mould is dependent on e.g. the viscosity of the epoxy, the permeability, the porosity of the reinforcement materials and the pressure gradient. Furthermore, the progression of the flow-front is highly dependent on both time and displacement of the front itself. This increases the challenge of estimating the flowfront. When epoxy is mixed it starts to cure immediately, increasing the viscosity over time. However, the continuous mixing of epoxy during the infusion of the mould decreases the cure rate. In addition, the viscosity of the epoxy is dependent on its temperature, with a higher temperature resulting in a lower viscosity but an increased cure rate. A decreased viscosity will result in an increased flow rate, but the epoxy curing process is an auto-catalytic exothermic process which will result in a fast increase of the temperature and thereby viscosity if heat cannot escape the blade mould. 384
Assumption 3: A flow model described by Darcy’s law in two spatial dimensions describes the general physics of a flow in three spatial dimensions. Remark 3: With decreasing thickness of the preform, a three-dimensional formulation can be approximated well by the two-dimensional formulation. The model of the flow-front propagation, can be formulated using PDEs determining the pressure as a function of space and time. We start with Darcy’s law in two spatial dimensions κH ∇p (1) q=− µ where q = q(x, y, t) is the flow rate of the fluid, κ = κ(x, y, t) is the permeability of the porous medium, H is a combined porosity and cross-sectional thickness measure, µ is the dynamic viscosity of the fluid, p = p(x, y, t) is the pressure, and ∇ = (∂x , ∂y ) is the spatial gradient. The flow rate is measured in m2 /s and is combined with the conservation of mass resulting in h˙ + ∇ · q = 0 (2) where h = h(x, y, t) ≤ H is the thickness of the fluid layer and h˙ indicates time derivative of h. The density can be eliminated by assuming a relationship p (3) h = min(H, ) ρg where ρ is the density of the fluid and g is the gravitational acceleration. Combining equations (2) and (3) results in,
While specifying the placement of glass fibre and auxiliary materials in a blade mould, it is ensured that the evolution of the flow-front during infusion process is homogeneous. However, sometimes due to manufacturing errors materials are handled and placed less carefully than specified, resulting in perturbations in the flow-front progression during the infusion process. Therefore, here both a homogeneous case, where the value of µκ is kept constant with respect to (w.r.t) time and displacement, and a heterogeneous case are simulated. In the heterogeneous case µκ is kept constant only w.r.t the time but along the spatial directions we choose to model it as κ c 0 (5) = 2πx µ 1 − A · cos Lx 1 − A · cos 2πy Ly
where A is a constant determining the relative decrease in permeability towards the middle, c0 being a correction constant and Lx = 0.8 m ; Ly = 0.9 m are the width and length of the casted area. This model describes a decrease of the permeability towards the centre of the cast area.
Fig. 2 shows a homogeneous flow-front progression with measurement noise added. The points are measurements at different, equally spaced, time steps along the y-axis where the preform is fully covered with epoxy. It is seen that the flow-front moves ahead evenly across the entire width of the preform. The time steps are evenly spaced showing that the velocity decreases with the progression 385
1.2 1.0 0.8 0.6 0.4 0.2 0.0
0.2
0.4
0.6
0.8
Flow front progression width (m)
Fig. 2. Simulated output data of a homogeneous flow-front from FEniCS with added measurement noise. The line colour changes gradually from red at the first time instance to green at the last time instance as the flowfront evolves with time. All time steps are equal. of the flow-front. Fig. 3 shows a heterogeneous flow-front
0.8
1.0
1.2
Measured flow
0.6
Assumption 4: Epoxy enters the mould evenly across a line perpendicular to the flow direction. Remark 4: The epoxy can only enter the mould through one or several inlet points. In an experimental setup a cavity can be created to ensure a homogeneous inlet across a line perpendicular to the flow direction.
0.0
0.4
The simulation model described above is used to emulate different flow-front patterns of the true system. The FEniCS software is configured such that the width is split into 64 cells and the length is split into 128 cells constructing a mesh of 64 × 128 cells and creating a total number of 65 × 129 vertices. The chosen amount of cells creates the opportunity for including a large number of lines to measure along, when trying to estimate a flow model and trying to determine the flow-front progression.
Measured flow
0.2
3.2 Simulated flow-front progression
385
0.0
κH dh p˙ = ∇ · ∇p (4) h˙ = dp µ This PDE is completed with boundary conditions: A pressure of zero at the outlet, a pressure of p0 = 1 bar at the inlet, and no-flux boundary conditions along the sides. The PDE is then solved with FEniCS software (Alnæs et al. (2015)) numerically by time marching, using a semi-implicit Euler step where the derivative dh dp is evaluated at the previous time step, the time derivative p˙ is approximated with a first order finite difference, and the right hand side is evaluated at the next time step. We note that dh dp vanishes in those parts of the domain that have already been filled with the fluid, so the system can be viewed as differential-algebraic. The numerical solution obtained from FEniCS is considered as our true data for this case study.
Flow front progression length (m)
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
Flow front progression length (m)
2018 IFAC SYSID July 9-11, 2018. Stockholm, Sweden
0.0
0.2
0.4
0.6
0.8
Flow front progression width (m)
Fig. 3. Simulated output data of a heterogeneous flow from FEniCS with added measurement noise. The line colour changes gradually from red at the first time instance to green at the last time instance as the flowfront evolves with time. All time steps are equal. progression with added measurement noise. In Fig. 3 it is seen how the decrease in the permeability also causes decreased flow rate towards the middle when comparing to the homogeneous flow shown in Fig. 2. From the discussion above it is seen that the dynamics of the flow-front can be formulated as spatio-temporal model describing the pressure as a function of time and space. This formulation is very useful to simulate the propagation of the flow-front inside the mould for different scenarios and operating conditions (e.g. permeability, temperature etc.) as demonstrated. Although these high-dimensional PDEs based models are good to understand and simulate a dynamical process they are not very useful when the final intended purpose of the
2018 IFAC SYSID 386 9-11, 2018. Stockholm, Sweden July
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
model is control and monitoring. Hence in the next section, we propose to use a coupled stochastic differential equation (SDE) (Øksendal (2010)) based modelling approach to model the flow-front dynamics. 4. LUMPED APPROXIMATION In this section, we introduce first briefly the concept of SDE and thereafter we show, how the coupled SDEs can be used to approximate the flow-front dynamics inside the mould along the spatial directions.
Y1,t
Y2,t
Y...
YN,t
4.1 Stochastic Differential Equations and Grey-box Models Grey-box models are typically SDEs based models (Kristensen et al. (2004)), where the structure of the model is built based on a combination of physical knowledge of the system, as in the white-box models, and on the statistical information based on the observations (measured data), as in the case of the black-box models. From a theoretical point of view, SDEs are the preferred choice to model stochastic, complex, and nonlinear systems where only a partial information about the system dynamics is available. An SDE based state-space model can be written as: dYt = f (Yt , Ut , t, θ)dt + σ(Ut , t, θ)dwt (6) Zk = h(Yk , Uk , tk , θ) + ek (7) The equations describing the dynamics of the states of the system, Yt , are formulated in continuous time and are separated in a drift term, f (Yt , Ut , t, θ), and a diffusion term, σ(Ut , t, θ). The observations, Zk , are linked to the states through the observation equation (7), which are typically formulated in discrete time and include the measurement error, ek . Ut represents the inputs and θ the parameters of the system. A clear separation of the residual error into diffusion and measurement noise results in a more correct description of the prediction error (Øksendal (2010)). Inclusion of the diffusion has another advantage, mainly related to the model building itself. By investigating the diffusion terms, one can retrieve information about how to improve an insufficient model (Kristensen et al. (2004)). Diffusion terms that are estimated to be relatively large indicate a model mismatch for the relevant part of the model. SDEs have been found to be useful within many areas of mathematical modelling of many dynamical systems (Øksendal (2010)). 4.2 Discretisation of the Spatial Domain The PDE model developed in the previous section is too complex for online monitoring and control purpose. Hence, we propose to use the coupled SDEs to model the flow-front dynamics for a limited amount of spatial discretisation points (as shown in Fig. 4) as compared to the PDE based model. First, note that one solution to the PDE in the homogeneous case, and for small H, is: y (8) p(x, y, t) = p0 · max(0, 1 − ), Yt where Yt is the position of the front, given by: κp0 1 dYt = , (9) dt µ Yt 386
Δx
Fig. 4. A pictorial representation, how the flow-front is approximated using the coupled SDEs. The red line shows the domain of spatial dimension x which is approximated using a SDE, whereas the dotted vertical lines show the points of discretisation for the spatial dimension x. and the material is flowing in the direction of the y-axis throughout the domain. Based on this, we discretise the x-axis using N grid points, and let Yi,t denote the position of the front at the corresponding value of x. No coupling between neighbouring grid points, results in N ordinary differential equations: κp0 1 dYi,t = , dt µ Yi,t
(10)
Since these equations assume perfect homogeneity, we introduce a noise term σi,t which parameterises all conceivable differences between the model and the true system. Finally, if the front is not perfectly aligned with the x-axis, the material flow will neither be perfectly aligned with the y-axis. This will stabilise the shape of the front, and this effect is included in the model through a diffusive term: κp0 1 Yi+1,t − 2Yi,t + Yi−1,t +D dt + σi,t dwt dYi,t = µ Yi,t (∆x)2 (11) For i = 2, . . . , N − 1, with similar equations for i = 1, N . For the actual estimation of a model, we choose a model where κpµ0 is combined to a single constant for each state D is treated as a single global equation, c0,i , and (∆x) 2 constant, D0 , for all state equations. This results in the following model being implemented c0,i + D0 · (Yi+1,t − 2Yi,t + Yi−1,t ) dt + σi,t dwt dYi,t = Yi,t (12) Except for the boundary cases where we model it as c0,1 + D0 · (Y2,t − Y1,t ) dt + σ1,t dwt (13) dY1,t = Y1,t We argue that these equations, while a stretch from the original physics, capture the essential dynamics of flowfront propagation and are in a form suitable for online flow-front estimation.
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
0.25
0.00 0.01 0.02 0.03 0.04 0.05 0.06
Sampling frequency (Hz)
0.08
Fig. 5. The RMSE values for the predicted states using models fitted for the homogeneous flow-front data. The data has been perturbed around the specified sampling frequencies for presentation.
0.07
Number of lines: 5 Number of lines: 8 Number of lines: 12 Noise σ: 0.10 Noise σ: 0.05 Noise σ: 0.01
5. RESULTS This section describes the results of using the Continuous Time Stochastic Modelling package in R to estimate models for the flow-front. For all model estimates the initial states are assumed to be known. The estimation model does not take the boundary conditions of the PDE into account. Therefore the data used for estimation only include the time steps where all the flow-front progression for all measurement lines are below 0.9 m. In the homogeneous case the value for the constants in equation (5) are A = 0.0 and c0 = 3.5 · 10−4 and for the heterogeneous we have chosen A = 0.9 and c0 = 1.5 · 10−4 . These have been chosen experimentally to create a flowfront that moves the distance from 0 meter to 0.9 meter to cover the entire area in approximately 600 seconds. For both the homogeneous and the heterogeneous case we perform the model estimation with 5, 8, and 12 line sensor measurements, each at sampling frequencies of 4, 10, and 20 seconds and noise levels with a standard deviation, σ of 0.01, 0.05, and 0.1. The root-mean-square errors (RMSE) are used to evaluate the estimated models. The resulting RMSE values are calculated as the mean of the RMSE values of the onestep ahead predictions of the states at each time step and these describe how well we estimate our true system. Fig. 5 and Fig. 6 show the relation between the RMSE values and the number of measurement lines, sampling frequency and measurement noise levels for the homogeneous and the 387
0.03
0.04
The corresponding value of the log-likelihood is the observed maximum likelihood value for that data set and model. All computations were done using the free statistical software, R (version 3.3.2), and the “CTSM-Rpackage” (Continuous Time Stochastic Modelling in R version 0.6.8-5, Juhl (2015)).
0.10
0.06
where θ is a set of parameters, ZN is the set of observations, L is the dimension of the observation space, and z0 is the initial condition. For a given set of parameters and initial states, k and Rk|k−1 are computed by a continuousdiscrete extended Kalman filter. The parameter estimates are found by maximising the log-likelihood: (15) θˆ = argmax log(L(θ; ZN )|z0 ) .
Number of lines: 5 Number of lines: 8 Number of lines: 12 Noise σ: 0.10 Noise σ: 0.05 Noise σ: 0.01
0.05
RMSE
The parameters, c0,i and D0 , of the coupled SDEs formulated above to model the evolution of the flow-front dynamics are estimated from the measured data simulated using the complex PDE based model. Because of the complex structure of SDEs, estimation of parameters is nontrivial except for some simple cases. Hence we propose to use a maximum likelihood method in combination with an extended Kalman filter method described in (Kristensen et al. (2004)) to estimate the parameters. The likelihood function is formulated using the one-step prediction errors, k , and the associated variances, Rk|k−1 : L(θ; ZN ) = p(ZN |θ) (14) N 1 T −1 exp − 2 k Rk|k−1 k = p(z0 |θ) √ det(Rk|k−1 )( 2π)L k=1
RMSE
4.3 Maximum Likelihood Estimation of SDEs
387
0.05
2018 IFAC SYSID July 9-11, 2018. Stockholm, Sweden
0.05
0.10
0.25
Sampling frequency (Hz)
Fig. 6. The RMSE values for the predicted states using models fitted for the heterogeneous flow-front data. The data has been perturbed around the specified sampling frequencies for presentation. heterogeneous cases, respectively. In the homogeneous case the RMSE value is constant w.r.t the number of measurement lines, however generally decreasing with decreasing measurement noise and increasing sampling frequency. For the heterogeneous flow the RMSE is decreasing with increasing number of measurement lines and is in general not affected by changes in measurement noise and sampling frequency. Fig. 7 shows an example of the estimated flowfront trajectories (solid lines) using an estimated flow-front model and the true flow-front trajectories, with measurement noise (dot-dashed), using the PDE model. 6. CONCLUSIONS AND FUTURE RESEARCH A good understanding of the flow-front dynamics is essential for designing a monitoring system for controlling the epoxy flow inside the mould in a VARTM process. In this paper, we have proposed a simulation based methodology
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
0.2
0.4
0.6
0.8
1.0
1.2
Estimated flow Measured flow
0.0
Flow front progression length (m)
2018 IFAC SYSID 388 9-11, 2018. Stockholm, Sweden July
0.0
0.2
0.4
0.6
0.8
Flow front progression width (m)
Fig. 7. The estimated flow-front with 12 lines (states), sampling frequency = 0.25 Hz, and the measurement noise σ = 0.1. The dot-dashed lines indicate the measurements from the true system and the solid lines indicate the estimated flow-front. The line colour changes gradually from red at the first time instance to green at the last time instance. to understand and model the flow-front dynamics. A PDE based model is proposed to simulate different patterns of the flow-front dynamics for various scenarios. Furthermore, a stochastic differential equations based modelling approach is proposed to estimate a grey-box model of the evolution of the flow-front dynamics inside the mould. The advantage of coupled SDEs based approach is that it gives the flexibility of combining the partial information about the physics behind the evolution of the flow-front inside mould and the measurements from various sensors (usually pressure or optical sensors in practice) into one parsimonious dynamical model. It can properly handle insufficient knowledge of the system and noise characteristics using the diffusion term in SDEs. The proposed methodology has been validated on simulation examples. The aim of our future research is to test the methodology for more complex flow patterns (e.g. a random change in permeability (w.r.t. displacement), increasing/stochastic viscosity (w.r.t time etc.) and to investigate its validity with real experimental data. REFERENCES Alnæs, M.S., Blechta, J., Hake, J., Johansson, A., Kehlet, B., Logg, A., Richardson, C., Ring, J., Rognes, M.E., and Wells, G.N. (2015). The FEniCS project version 1.5. Archive of Numerical Software, 3(100). doi:10.11588/ ans.2015.100.20553. Dominauskas, A., Heider, D., and Gillespie, J.W. (2003). Electric time-domain reflectometry sensor for online flow sensing in liquid composite molding processing. Composites Part A: Applied Science and Manufacturing, 34(1), 67 – 74. doi:https://doi.org/10.1016/ S1359-835X(02)00232-4. Juhl, R. (2015). CTSM-R: Continuous Time Stochastic Modelling in R: User’s Guide and Reference Manual. R package version 0.6.8-5. 388
Kristensen, N.R., Madsen, H., and Jørgensen, S.B. (2004). Parameter estimation in stochastic grey-box models. Automatica, 40(2), 225 – 237. doi:https://doi.org/10. 1016/j.automatica.2003.10.001. Kueh, S.R., Parnas, R.S., and Advani, S.G. (2002). A methodology for using long-period gratings and moldfilling simulations to minimize the intrusiveness of flow sensors in liquid composite molding. Composites Science and Technology, 62(2), 311 – 327. doi:https://doi.org/ 10.1016/S0266-3538(01)00217-2. Matsuzaki, R. and Shiota, M. (2016). Data assimilation through integration of stochastic resin flow simulation with visual observation during vacuum-assisted resin transfer molding: A numerical study. Composites Part A: Applied Science and Manufacturing, 84(Supplement C), 43 – 52. doi:https://doi.org/10.1016/j.compositesa. 2016.01.006. Matuzaki, R., Seto, D., Naito, M., Todoroki, A., and Mizutani, Y. (2015). Analytical prediction of void formation in geometrically anisotropic woven fabrics during resin transfer molding. Composites Science and Technology, 107(Supplement C), 154 – 161. doi:https:// doi.org/10.1016/j.compscitech.2014.12.013. Nauheimer, M., Relan, R., Madsen, H., Thygesen, U.H., Olesen, B., and Kirkeby, K. (2017). A combined experimental and simulation based approach to model the flow-front dynamics in the vacuum assisted resin transfer moulding process. In 2nd International Conference on Lightweight Design of Materials and Engineering Structures (LIMAS 2017), 28–35. Øksendal, B. (2010). Stochastic Differential Equations: An Introduction with Applications. Universitext. Springer Berlin Heidelberg. Park, C.H., Lebel, A., Saouab, A., Br´eard, J., and Lee, W.I. (2011). Modeling and simulation of voids and saturation in liquid composite molding processes. Composites Part A: Applied Science and Manufacturing, 42(6), 658 – 668. doi:https://doi.org/10.1016/j.compositesa. 2011.02.005. Sreekumar, P., Joseph, K., Unnikrishnan, G., and Thomas, S. (2007). A comparative study on mechanical properties of sisal-leaf fibre-reinforced polyester composites prepared by resin transfer and compression moulding techniques. Composites Science and Technology, 67(3), 453 – 461. doi:https://doi.org/10.1016/j.compscitech. 2006.08.025. Stiesdal, H., Enevoldsen, P., Johansen, K., Kristensen, J., N¨ortem, M., and Winther-Jensen, M. (2006). Method for manufacturing windmill blades. EP Patent 1,310,351. Whitaker, S. (1986). Flow in porous media i: A theoretical derivation of darcy’s law. Transport in Porous Media, 1(1), 3–25. doi:10.1007/BF01036523. Yenilmez, B. and Sozer, E.M. (2009). A grid of dielectric sensors to monitor mold filling and resin cure in resin transfer molding. Composites Part A: Applied Science and Manufacturing, 40(4), 476 – 489. doi:https://doi. org/10.1016/j.compositesa.2009.01.014. Zhang, F., Cosson, B., Comas-Cardona, S., and Binetruy, C. (2011). Efficient stochastic simulation approach for rtm process with random fibrous permeability. Composites Science and Technology, 71(12), 1478 – 1485. doi: https://doi.org/10.1016/j.compscitech.2011.06.006.
ScienceDirect
IFAC PapersOnLine 51-15 (2018) 383–388
A Stochastic Spatio-Temporal Model of the A Stochastic Spatio-Temporal Model of the A Stochastic Spatio-Temporal Model of the Flow-Front Dynamics in a Vacuum Assisted A Stochastic Spatio-Temporal Model of the Flow-Front Dynamics in a Vacuum Assisted Flow-Front Dynamics in a Vacuum Assisted Resin Transfer Moulding Process Flow-Front Dynamics in a Vacuum Assisted Resin Transfer Moulding Process Resin Transfer Moulding Process Resin Transfer Moulding Process ∗∗ Michael Nauheimer ∗∗ Rishi Relan ∗∗ ∗∗ Uffe Høgsbro Thygesen ∗∗
Michael Nauheimer Relan Thygesen ∗ ∗∗ ∗∗ ∗∗ ∗ Høgsbro ∗ ∗ Rishi ∗∗ Uffe ∗∗ Michael Nauheimer Rishi Relan Uffe Høgsbro Thygesen Henrik Madsen Olesen Kirkeby Michael Nauheimer Relan Høgsbro Thygesen ∗∗ Bendt ∗ Klaus ∗ ∗ Rishi ∗∗ Uffe ∗∗ Henrik Madsen Olesen Kirkeby ∗∗ Bendt ∗ Klaus ∗ Michael Nauheimer Rishi Relan Uffe Høgsbro Thygesen ∗∗ ∗ ∗ Henrik Madsen Madsen Bendt Bendt Olesen Olesen Klaus Klaus Kirkeby Kirkeby Henrik ∗∗ ∗ ∗ Henrik Madsen Bendt Olesen Klaus Kirkeby ∗ ∗ Siemens Gamesa Renewable Energy, Denmark ∗ Siemens Gamesa Renewable Energy, Denmark ∗ Siemens Gamesa Denmark (e-mail: [email protected], [email protected], Gamesa Renewable Renewable Energy, Energy, Denmark ∗ Siemens (e-mail: [email protected], [email protected], Siemens Gamesa Renewable Energy, Denmark (e-mail: [email protected], [email protected], [email protected]) (e-mail: [email protected], [email protected], [email protected]) ∗∗ (e-mail: [email protected], [email protected], [email protected]) Technical University of Denmark, Denmark [email protected]) ∗∗ DTU Compute, Technical University Denmark, ∗∗ DTU Compute, [email protected]) ∗∗ DTU Compute, Technical [email protected], University of of [email protected]) Denmark, Denmark Denmark (e-mail: [email protected], DTU Compute, Technical University of Denmark, Denmark ∗∗ (e-mail: [email protected], DTU Compute, Technical [email protected], University of [email protected]) Denmark, Denmark (e-mail: [email protected], [email protected], [email protected]) (e-mail: [email protected], [email protected], [email protected]) (e-mail: [email protected], [email protected], [email protected]) Abstract: With With more more emphasis emphasis on on the the use use of of green green energy, energy, the the size size of of the the turbines turbines and and blades blades Abstract: Abstract: With more emphasis on the use of green energy, the size of the turbines and blades in the wind turbines is continuously increasing. With increasing blade size, the casting process Abstract: With more emphasis on the use of green energy, the size of the turbines and blades in the wind With turbines is continuously increasing. With increasing blade the casting process Abstract: more on risk the use of green energy, the size of size, the turbines andmade blades in the turbines is continuously increasing. With increasing blade size, the casting process becomes more complicated and the the of faults faults increases. Production of such such blades, of in the wind wind turbines is emphasis continuously increasing. With increasing blade size, the blades, casting process becomes more complicated and risk of increases. Production of made of in the wind turbines is continuously increasing. With increasing blade size, the casting process becomes more complicated complicated and the the risk risk of faults faults increases. Production of such such blades, blades, made of fibre reinforced polymer composites, without the possibility of visual inspection of the infusion becomes more and of increases. Production of made of fibre reinforced polymer composites, without the possibility of visual inspection of the infusion becomes morefor complicated and the(possibly risk of faults Production ofprocess. such blades, of fibre reinforced polymer without the possibility of inspection of the infusion process calls a sensor sensorcomposites, system virtual) for monitoring monitoring the Therefore, in fibre reinforced polymer composites, without the increases. possibility of visual visualthe inspection ofTherefore, the made infusion process calls for a system (possibly virtual) for process. in fibre reinforced polymer composites, without the possibility of visual inspection of the infusion process calls for a sensor system (possibly virtual) for monitoring the process. Therefore, in this paper, we propose a stochastic spatio-temporal modelling approach to model the flow-front process calls for a sensor system (possibly virtual)modelling for monitoring the to process. Therefore, in this paper, we propose a spatio-temporal approach model flow-front process calls a sensor system (possibly for monitoring the to process. Therefore, this paper, wefor propose a stochastic stochastic spatio-temporal modelling approach to model the the flow-front dynamics inside the vacuum assisted resin virtual) transfermodelling moulding process, knowledge offlow-front which in is this paper, we propose a stochastic spatio-temporal approach model the dynamics inside the vacuum assisted resin transfer moulding process, knowledge of which this paper, we propose a stochastic spatio-temporal modelling approach model theof dynamics inside the vacuum vacuum assisted resinoftransfer transfer moulding process, to knowledge offlow-front which is is essential for determining the current state the infusion process. dynamics inside the assisted resin moulding process, knowledge which is essential for determining the current the infusion process. dynamics inside the vacuum assistedstate resinof moulding process, knowledge of which is essential determining the state of the process. essential for for determining the current current state oftransfer the infusion infusion process. © 2018, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. essential for determining the current state of the infusion process. Keywords: Wind Wind turbines, turbines, Partial Partial differential differential equations, equations, Stochastic Stochastic differential differential equations, equations, Keywords: Keywords: Likelihood, Wind turbines, turbines, PartialKalman differential equations, Stochastic differential equations, equations, Maximum Extended filter, Grey-box models Keywords: Wind Partial differential equations, Stochastic differential filter, Grey-box models Maximum Extended Keywords: Wind turbines, PartialKalman differential Maximum Likelihood, Likelihood, Extended Kalman filter, Grey-box models Maximum Likelihood, Extended Kalman filter,equations, Grey-boxStochastic models differential equations, Maximum Likelihood, Extended Kalman filter, Grey-box models 1. INTRODUCTION INTRODUCTION (2003)). In general these sensors are accurate but limited 1. (2003)). In these sensors are accurate but limited 1. INTRODUCTION INTRODUCTION (2003)). In general general these sensors are accurate butmoulded limited to measuring on or close to the surface of the 1. (2003)). In general these sensors are accurate but limited to measuring on or close to the surface of the moulded 1. INTRODUCTION (2003)). In general these sensors are accurate but limited to measuring measuring on or or closeinto tothethe the surface of the the based moulded parts. It has been shown past that sensors on Large scale composite shell structures like wind turbine to on close surface of moulded parts. It has been shown in the past that sensors based on Large scale composite shell structures like wind turbine to measuring on or close to the surface of the moulded parts. It has been shown in the past that sensors based on Large scale composite shell structures like wind turbine optical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer parts. It has been shown in the past that sensors based on Large scale composite shell structures like wind turbine optical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer parts. It has been shown in the past that sensors based on Large scale composite shell structures like wind turbine optical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer making it possible to measure through thickness of the structure. These structures are usually manufactured usoptical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer making it possible to measure through thickness of the structure. These structures are usually manufactured usoptical fibres (Kueh et al. (2002)) can be cast into parts blades consist of a glass fibre reinforced thermoset polymer making it possible to measure through thickness of the structure. These structures are usually manufactured usshell. Recently, (Matsuzaki and Shiota (2016)) proposed ing the vacuum assisted resin transfer moulding (VARTM) making it possible to measure through thickness of the structure. Theseassisted structures are usually manufactured us- shell. Recently, (Matsuzaki and Shiota (2016)) proposed ing the vacuum resin transfer moulding (VARTM) it using possible to measure thickness ofdata the structure. These structures usually manufactured us- making Recently, (Matsuzaki and through Shiota (2016)) (2016)) proposed ing the vacuum vacuum assisted resin transfer moulding (VARTM) ashell. method two-sided visual observations for 4D process. The VARTM process involves infusing a therthershell. Recently, (Matsuzaki and Shiota proposed ing the assisted resinare transfer moulding (VARTM) ashell. method using two-sided visual observations for 4D data process. The VARTM process involves infusing a Recently, (Matsuzaki and Shiota (2016)) proposed ing the vacuum assisted resin transfer moulding (VARTM) a method using two-sided visual observations for 4D data process. The VARTM process involves infusing a therto accurately reconstruct the 3D flow moset polymer polymer into a fibre fibre reinforcement preform.aEpoxy Epoxy method using two-sided visual observations forresin 4D data process. The VARTM process involves infusing ther- aassimilation to accurately reconstruct the 3D flow moset into reinforcement preform. aand method using two-sided observations forresin 4D data process. The VARTM process involves infusing aEpoxy ther- assimilation assimilation to accurately reconstruct the 3D resin flow moset polymer into aaatofibre fibre reinforcement preform. Epoxy permeability filed of aavisual fibre preform. resin is transferred the mould through several inlets. assimilation to accurately reconstruct the 3D resin flow moset polymer into reinforcement preform. and permeability filed of fibre preform. resin is transferred to the mould through several inlets. assimilation to accurately reconstruct the 3D resin flow moset polymer into a fibre reinforcement preform. Epoxy and permeability filed of a fibre preform. resin is transferred to the mould through several inlets. Progression of the epoxy inside the mould is only driven by and permeability filed of a fibre preform. resin is transferred to the mould through several inlets. However, some blades produced at the Siemens Gamesa Progression of the epoxy inside the mould is only driven by and permeability filed of a fibre preform. resin is transferred to the mould through several inlets. blades produced at Gamesa Progression of the epoxy epoxy inside the the mould mould is is only driven by However, the pressure pressure of difference between and thedriven pressure Progression the inside only by However, some some blades produced at the theinSiemens Siemens Gamesa Renewable Energy factories are casted one piece using the difference between the mould mould is and the pressure some blades produced at the Siemens Gamesa Progression of the epoxy inside the only driven by However, Renewable Energy factories are casted in one piece using the pressure difference between and the pressure of one bar surrounding the mould in the production area. the pressure difference between the mould and the pressure R However, some blades produced at the Siemens Gamesa Renewable Energy factories are casted in one piece using the patented IntegralBlades technology (Stiesdal et al. of one bar surrounding the mould in the production area. Renewable Energy factories are casted in one piece using R the pressure difference between the mould and the pressure the patented IntegralBlades technology (Stiesdal et al. of one bar surrounding the mould in the production area. R of onetobar surrounding the mould inofthe production area. Renewable Energy factories are casted in one piece using R technology the patented IntegralBlades (Stiesdal et al. (2006)) instead of two separate half blades which are Due the inhomogeneous nature the flow inside the the patented IntegralBlades technology (Stiesdal et al. of onetobar surrounding the mould inofthe production area. R (2006)) instead of two separate half blades which are Due the inhomogeneous nature the flow inside the the patented IntegralBlades technology (Stiesdal et al. (2006)) instead of two separate half blades which are Due to the inhomogeneous nature of the flow inside the then united afterwards. One of the main disadvantages mould, risk of moulding defects such as dry spots (2006)) instead of two separate half blades which are Due to the inhomogeneous nature of the flow inside the then united afterwards. One of the main disadvantages mould, the inhomogeneous risk of moulding defects such as dry spots (2006)) instead of two separate half blades which are Due to the nature of the flow inside the then united afterwards. One of the main disadvantages mould, risk of moulding defects such as dry spots of this technology is that, it is not possible to visually and voids increases and this leads to deterioration of the then united afterwards. One of the main disadvantages mould, theincreases risk of and moulding defects such as dry ofspots of technology is it is possible to and voids this leads to deterioration the afterwards. One the main disadvantages mould, theincreases risk of and moulding defects such as dry of spots of this thisunited technology isof that, that, it of is not not possible to visually visually and voids increases and this leads to deterioration deterioration ofet the inspect the process transferring the epoxy resin to the mechanical properties of this the leads cast parts (Sreekumar al. then of this technology is that, it is not possible to visually and voids to the inspect the process of transferring the epoxy resin to the mechanical properties of the cast parts (Sreekumar et al. this technology is that, it is not possible to Gamesa visually and voids increases and leads toetdeterioration ofet inspect the process of transferring the epoxy resin to the mechanical properties of this the cast parts parts (Sreekumar et the al. of mould. Hence, the main motivation at Siemens (2007); Park et al. (2011); Matuzaki al. (2015)). Hence, inspect the process of transferring the epoxy resin to the mechanical properties of the cast (Sreekumar al. mould. Hence, the main motivation at Siemens Gamesa (2007); Park et al. (2011); Matuzaki et al. (2015)). Hence, inspect the process of transferring the epoxy resin to the mechanical properties of the cast parts (Sreekumar et al. mould. Hence, the main motivation at Siemens Gamesa (2007); Park et al. (2011); Matuzaki et al. (2015)). Hence, Renewable Energy is to develop a flow-front monitoring to prevent any potential fault from developing and for mould. Hence, the main motivation at Siemens Gamesa (2007); Park et al. (2011); Matuzaki et al. (2015)). Hence, Renewable Energy is to develop aa flow-front monitoring to prevent any potential fault from developing and for mould. Hence, the main motivation at Siemens Gamesa (2007); Park et al. (2011); Matuzaki et al. (2015)). Hence, Renewable Energy is to develop flow-front monitoring to prevent any potential fault from developing and for system which eventually aids in controlling the process better control of the process, it is necessary to monitor the Renewable Energy is to develop a flow-front monitoring to prevent any potential fault developing and the for system which eventually aids in controlling the process better control of the process, it isfrom necessary to monitor Energy is to develop flow-front monitoring to prevent fault developing and the for Renewable system which eventually aids trajectory in acontrolling controlling the process better control of potential the process, process, it necessary to monitor of mould filling. A controlled of the flow-front evolution ofany the flow-front inside the mould mould inmonitor a VARTM VARTM system which eventually aids in the process better control of the it is isfrom necessary toin the of mould filling. A controlled trajectory of the flow-front evolution of the flow-front inside the a system which eventually aids in controlling the process better control of the process, it is necessary to monitor the of mould filling. A controlled trajectory of the flow-front evolution of the the flow-front flow-front inside inside the the mould mould in in aa VARTM VARTM of inside the mould decreases the risk of areas with dry glass, process. mould filling. A controlled trajectory of the flow-front evolution of inside the mould decreases the trajectory risk of areas with dry glass, process. of mould filling. A controlled of the flow-front evolution of the flow-front inside the mould in a VARTM inside the mould risk of areas with dry glass, process. which eventually decreases the repair time and increase the inside the mould risk of areas with dry glass, process. which eventually decreases the repair time and increase the Sensors for monitoring the VARTM processes is an acinside the mould decreases the risk of areas with dry glass, process. which eventually eventually repair time and and increase increase the general quality of the produced blade. Sensors for monitoring the VARTM processes is an acwhich decreases the repair time the general quality of the produced blade. Sensors for monitoring the VARTM VARTM processes is an preactive field of research. Several types of sensors have Sensors for monitoring the processes is an acwhich eventually decreases the repair time and increase the general quality quality of of the the produced produced blade. blade. tive field of research. Several types of sensors have pregeneral Sensors for thein VARTM is an preac- To tive field of monitoring research. Several types of ofprocesses sensors have have preachieve this goal, this paper proposes a spatio-temporal viously been implemented VARTM processes including tive field of research. Several types sensors general quality of the produced blade. this this paper aa spatio-temporal viously been implemented in processes including tive field of research. Several types sensors have pre- To To achieve achieve this goal, goal,the thisdynamics paper proposes proposes spatio-temporal viously been implemented in VARTM VARTM processes including model to estimate of the inside permittivity sensors (Yenilmez (Yenilmez and of Sozer (2009)), presTo achieve this goal, this paper proposes aflow-front spatio-temporal viously been implemented in VARTM processes including model to estimate the dynamics of the inside permittivity sensors and Sozer (2009)), presTo achieve this goal, this paper proposes aflow-front spatio-temporal viously been implemented in VARTM processes including model to estimate estimate the the dynamics of the thewe flow-front inside permittivity sensors et (Yenilmez andand Sozer (2009)), presthe mould. To emulate real system, first model the sure sensors (Zhang al. (2011)), sensors based on model to the dynamics of flow-front inside permittivity sensors (Yenilmez and Sozer (2009)), presthe mould. To emulate the real system, we first model the sure sensors (Zhang et al. (2011)), and sensors based on to estimate the the dynamics of thewe flow-front inside permittivity sensors et (Yenilmez andand Sozer (2009)), presthe mould. To emulate real system, first model the sure sensors (Zhang et al. (2011)), (2011)), and sensors based on model spatio-temporal evolution of the flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. the mould. To emulate the real system, we first model the sure sensors (Zhang al. sensors based on spatio-temporal evolution of the flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. the mould. To emulate the real system, we first model the sure sensors (Zhang et al. (2011)), and sensors based on spatio-temporal evolution of the flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. partial differential equations (PDEs) based on Darcy’s evolution of the flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. spatio-temporal partial differential equations based on Darcy’s spatio-temporal evolution of the(PDEs) flow-front dynamics using electrical time-domain reflectometry (Dominauskas et al. partial differential equations (PDEs) based on Darcy’s Law (Whitaker (1986)). The PDE based model is then We thank the Manufacturing Academy of Denmark (MADE) for partial differential equations (PDEs) based on Darcy’s We thank the Manufacturing Academy of Denmark (MADE) for Law (1986)). The PDE based model is We thank partial differential equations based onevolution Darcy’s Law (Whitaker (Whitaker (1986)). The (PDEs) PDE based model is then then providing thethe financial support. We thank of theDenmark laboratory employees the Manufacturing Academy of (MADE) for simulated to generate different patterns of the Law (Whitaker (1986)). The PDE based model is then We thank Manufacturing Academy (MADE) for providing the financial support. We thank theDenmark laboratory employees simulated to generate different patterns of the evolution Law (Whitaker (1986)). The PDE basedof model is then We thank Manufacturing Academy (MADE) for at Siemens Gamesa Renewable Energy forDenmark their assistance and providing the financial support. thank the laboratory employees simulated to generate different patterns the evolution of the flow-front inside the mould. providing thethe financial support. We We thank of the employees simulated to generate different patterns of the evolution at Siemens Gamesa Renewable Energy for laboratory their assistance and of the flow-front inside the mould. providing theGamesa financialRenewable support. We thank the employees suggestions. at Siemens Energy for their and simulated to generate patterns of the evolution at Siemens Gamesa Renewable Energy for laboratory their assistance assistance and of inside the of the the flow-front flow-front inside different the mould. mould. suggestions. at Siemens Gamesa Renewable Energy for their assistance and suggestions. of the flow-front inside the mould. suggestions.
suggestions. 2405-8963 © © 2018 2018, IFAC IFAC (International Federation of Automatic Control) Copyright 383 Hosting by Elsevier Ltd. All rights reserved. Copyright © under 2018 IFAC 383 Control. Peer review responsibility of International Federation of Automatic Copyright © 2018 IFAC 383 Copyright © 2018 IFAC 383 10.1016/j.ifacol.2018.09.175 Copyright © 2018 IFAC 383
2018 IFAC SYSID 384 9-11, 2018. Stockholm, Sweden July
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
A PDE based model is generally not a preferred choice for online monitoring and control. Hence, the main contribution of this paper is to develop a spatio-temporal model of the flow-front dynamics using coupled stochastic differential equations (SDEs) (Øksendal (2010)). The validation of the coupled SDEs based model is performed for different patterns of flow-front dynamics generated using the PDE based model. Finally, an analysis of the level of measurement noise, sampling rate and the model complexity on the final estimated model is also presented. The paper is structured as follows: Section 2 briefly describes the VARTM process. In Section 3, the mathematical background and implementation of PDEs based flowfront simulation model is discussed. A model for estimation of the flow-front based on coupled SDEs is proposed in Section 4. Section 5 describes the results and finally the conclusions are given in Section 6. 2. VARTM PROCESS In this section, we briefly describe the VARTM process for blade manufacturing. Fig. 1 shows the work-flow of the VARTM process. Flow mesh
Vacuum foil
Resin inlet
Release film
Fibre reinforcement Vacuum hose
Aluminum plate
3. SIMULATING FLOW IN A POROUS MEDIUM To deal with the above mentioned challenges in practice, first it is important to completely understand the dynamics of the flow-front progression inside the mould. For that purpose, we first propose to simulate the flow-front progression in a VARTM process in a 80 cm × 90 cm rectangular casting for our simulation study. We have implemented Darcy’s law with relevant boundary conditions into a PDE solver software, FEniCS Project (Alnæs et al. (2015)), and use the output to create flow-front pattern data representative of the data which can be measured experimentally using commercially available sensors and instrumentation (Nauheimer et al. (2017)). In the section below, we briefly explain the mathematical background behind the flow-front modelling. 3.1 Flow-front simulation model Assumption 1: Epoxy is considered to be a Newtonian fluid. Remark 1: Uncured epoxy is not a Newtonian fluid. It is assumed that minimal or equal distributed physical force will be exerted on the epoxy within the mould and thereby it is assumed that it can be considered a Newtonian fluid for our application. Assumption 2: The increase in viscosity is small and constant enough to be omitted for shorter time scales. Remark 2: The curing of epoxy initiates as soon as it has been mixed therefore, as the epoxy cures, the viscosity increases. It is assumed that for the short time frame of our simulations, the increase in viscosity is small enough to be omitted.
Sealant
Fig. 1. The general material configuration and work-flow of a VARTM process. The resin enters the vacuum bag through the resin inlet and spreads into the pores of the fibre reinforcement from left to right. This process is only assisted by the vacuum generated through the vacuum hose. The VARTM process is a special variety of Resin Transfer Moulding (RTM) processes. As the name indicates the process is only assisted by the difference of the nearvacuum inside the mould and the surrounding pressure of approximately one bar. The flow of epoxy inside the mould is dependent on e.g. the viscosity of the epoxy, the permeability, the porosity of the reinforcement materials and the pressure gradient. Furthermore, the progression of the flow-front is highly dependent on both time and displacement of the front itself. This increases the challenge of estimating the flowfront. When epoxy is mixed it starts to cure immediately, increasing the viscosity over time. However, the continuous mixing of epoxy during the infusion of the mould decreases the cure rate. In addition, the viscosity of the epoxy is dependent on its temperature, with a higher temperature resulting in a lower viscosity but an increased cure rate. A decreased viscosity will result in an increased flow rate, but the epoxy curing process is an auto-catalytic exothermic process which will result in a fast increase of the temperature and thereby viscosity if heat cannot escape the blade mould. 384
Assumption 3: A flow model described by Darcy’s law in two spatial dimensions describes the general physics of a flow in three spatial dimensions. Remark 3: With decreasing thickness of the preform, a three-dimensional formulation can be approximated well by the two-dimensional formulation. The model of the flow-front propagation, can be formulated using PDEs determining the pressure as a function of space and time. We start with Darcy’s law in two spatial dimensions κH ∇p (1) q=− µ where q = q(x, y, t) is the flow rate of the fluid, κ = κ(x, y, t) is the permeability of the porous medium, H is a combined porosity and cross-sectional thickness measure, µ is the dynamic viscosity of the fluid, p = p(x, y, t) is the pressure, and ∇ = (∂x , ∂y ) is the spatial gradient. The flow rate is measured in m2 /s and is combined with the conservation of mass resulting in h˙ + ∇ · q = 0 (2) where h = h(x, y, t) ≤ H is the thickness of the fluid layer and h˙ indicates time derivative of h. The density can be eliminated by assuming a relationship p (3) h = min(H, ) ρg where ρ is the density of the fluid and g is the gravitational acceleration. Combining equations (2) and (3) results in,
While specifying the placement of glass fibre and auxiliary materials in a blade mould, it is ensured that the evolution of the flow-front during infusion process is homogeneous. However, sometimes due to manufacturing errors materials are handled and placed less carefully than specified, resulting in perturbations in the flow-front progression during the infusion process. Therefore, here both a homogeneous case, where the value of µκ is kept constant with respect to (w.r.t) time and displacement, and a heterogeneous case are simulated. In the heterogeneous case µκ is kept constant only w.r.t the time but along the spatial directions we choose to model it as κ c 0 (5) = 2πx µ 1 − A · cos Lx 1 − A · cos 2πy Ly
where A is a constant determining the relative decrease in permeability towards the middle, c0 being a correction constant and Lx = 0.8 m ; Ly = 0.9 m are the width and length of the casted area. This model describes a decrease of the permeability towards the centre of the cast area.
Fig. 2 shows a homogeneous flow-front progression with measurement noise added. The points are measurements at different, equally spaced, time steps along the y-axis where the preform is fully covered with epoxy. It is seen that the flow-front moves ahead evenly across the entire width of the preform. The time steps are evenly spaced showing that the velocity decreases with the progression 385
1.2 1.0 0.8 0.6 0.4 0.2 0.0
0.2
0.4
0.6
0.8
Flow front progression width (m)
Fig. 2. Simulated output data of a homogeneous flow-front from FEniCS with added measurement noise. The line colour changes gradually from red at the first time instance to green at the last time instance as the flowfront evolves with time. All time steps are equal. of the flow-front. Fig. 3 shows a heterogeneous flow-front
0.8
1.0
1.2
Measured flow
0.6
Assumption 4: Epoxy enters the mould evenly across a line perpendicular to the flow direction. Remark 4: The epoxy can only enter the mould through one or several inlet points. In an experimental setup a cavity can be created to ensure a homogeneous inlet across a line perpendicular to the flow direction.
0.0
0.4
The simulation model described above is used to emulate different flow-front patterns of the true system. The FEniCS software is configured such that the width is split into 64 cells and the length is split into 128 cells constructing a mesh of 64 × 128 cells and creating a total number of 65 × 129 vertices. The chosen amount of cells creates the opportunity for including a large number of lines to measure along, when trying to estimate a flow model and trying to determine the flow-front progression.
Measured flow
0.2
3.2 Simulated flow-front progression
385
0.0
κH dh p˙ = ∇ · ∇p (4) h˙ = dp µ This PDE is completed with boundary conditions: A pressure of zero at the outlet, a pressure of p0 = 1 bar at the inlet, and no-flux boundary conditions along the sides. The PDE is then solved with FEniCS software (Alnæs et al. (2015)) numerically by time marching, using a semi-implicit Euler step where the derivative dh dp is evaluated at the previous time step, the time derivative p˙ is approximated with a first order finite difference, and the right hand side is evaluated at the next time step. We note that dh dp vanishes in those parts of the domain that have already been filled with the fluid, so the system can be viewed as differential-algebraic. The numerical solution obtained from FEniCS is considered as our true data for this case study.
Flow front progression length (m)
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
Flow front progression length (m)
2018 IFAC SYSID July 9-11, 2018. Stockholm, Sweden
0.0
0.2
0.4
0.6
0.8
Flow front progression width (m)
Fig. 3. Simulated output data of a heterogeneous flow from FEniCS with added measurement noise. The line colour changes gradually from red at the first time instance to green at the last time instance as the flowfront evolves with time. All time steps are equal. progression with added measurement noise. In Fig. 3 it is seen how the decrease in the permeability also causes decreased flow rate towards the middle when comparing to the homogeneous flow shown in Fig. 2. From the discussion above it is seen that the dynamics of the flow-front can be formulated as spatio-temporal model describing the pressure as a function of time and space. This formulation is very useful to simulate the propagation of the flow-front inside the mould for different scenarios and operating conditions (e.g. permeability, temperature etc.) as demonstrated. Although these high-dimensional PDEs based models are good to understand and simulate a dynamical process they are not very useful when the final intended purpose of the
2018 IFAC SYSID 386 9-11, 2018. Stockholm, Sweden July
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
model is control and monitoring. Hence in the next section, we propose to use a coupled stochastic differential equation (SDE) (Øksendal (2010)) based modelling approach to model the flow-front dynamics. 4. LUMPED APPROXIMATION In this section, we introduce first briefly the concept of SDE and thereafter we show, how the coupled SDEs can be used to approximate the flow-front dynamics inside the mould along the spatial directions.
Y1,t
Y2,t
Y...
YN,t
4.1 Stochastic Differential Equations and Grey-box Models Grey-box models are typically SDEs based models (Kristensen et al. (2004)), where the structure of the model is built based on a combination of physical knowledge of the system, as in the white-box models, and on the statistical information based on the observations (measured data), as in the case of the black-box models. From a theoretical point of view, SDEs are the preferred choice to model stochastic, complex, and nonlinear systems where only a partial information about the system dynamics is available. An SDE based state-space model can be written as: dYt = f (Yt , Ut , t, θ)dt + σ(Ut , t, θ)dwt (6) Zk = h(Yk , Uk , tk , θ) + ek (7) The equations describing the dynamics of the states of the system, Yt , are formulated in continuous time and are separated in a drift term, f (Yt , Ut , t, θ), and a diffusion term, σ(Ut , t, θ). The observations, Zk , are linked to the states through the observation equation (7), which are typically formulated in discrete time and include the measurement error, ek . Ut represents the inputs and θ the parameters of the system. A clear separation of the residual error into diffusion and measurement noise results in a more correct description of the prediction error (Øksendal (2010)). Inclusion of the diffusion has another advantage, mainly related to the model building itself. By investigating the diffusion terms, one can retrieve information about how to improve an insufficient model (Kristensen et al. (2004)). Diffusion terms that are estimated to be relatively large indicate a model mismatch for the relevant part of the model. SDEs have been found to be useful within many areas of mathematical modelling of many dynamical systems (Øksendal (2010)). 4.2 Discretisation of the Spatial Domain The PDE model developed in the previous section is too complex for online monitoring and control purpose. Hence, we propose to use the coupled SDEs to model the flow-front dynamics for a limited amount of spatial discretisation points (as shown in Fig. 4) as compared to the PDE based model. First, note that one solution to the PDE in the homogeneous case, and for small H, is: y (8) p(x, y, t) = p0 · max(0, 1 − ), Yt where Yt is the position of the front, given by: κp0 1 dYt = , (9) dt µ Yt 386
Δx
Fig. 4. A pictorial representation, how the flow-front is approximated using the coupled SDEs. The red line shows the domain of spatial dimension x which is approximated using a SDE, whereas the dotted vertical lines show the points of discretisation for the spatial dimension x. and the material is flowing in the direction of the y-axis throughout the domain. Based on this, we discretise the x-axis using N grid points, and let Yi,t denote the position of the front at the corresponding value of x. No coupling between neighbouring grid points, results in N ordinary differential equations: κp0 1 dYi,t = , dt µ Yi,t
(10)
Since these equations assume perfect homogeneity, we introduce a noise term σi,t which parameterises all conceivable differences between the model and the true system. Finally, if the front is not perfectly aligned with the x-axis, the material flow will neither be perfectly aligned with the y-axis. This will stabilise the shape of the front, and this effect is included in the model through a diffusive term: κp0 1 Yi+1,t − 2Yi,t + Yi−1,t +D dt + σi,t dwt dYi,t = µ Yi,t (∆x)2 (11) For i = 2, . . . , N − 1, with similar equations for i = 1, N . For the actual estimation of a model, we choose a model where κpµ0 is combined to a single constant for each state D is treated as a single global equation, c0,i , and (∆x) 2 constant, D0 , for all state equations. This results in the following model being implemented c0,i + D0 · (Yi+1,t − 2Yi,t + Yi−1,t ) dt + σi,t dwt dYi,t = Yi,t (12) Except for the boundary cases where we model it as c0,1 + D0 · (Y2,t − Y1,t ) dt + σ1,t dwt (13) dY1,t = Y1,t We argue that these equations, while a stretch from the original physics, capture the essential dynamics of flowfront propagation and are in a form suitable for online flow-front estimation.
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
0.25
0.00 0.01 0.02 0.03 0.04 0.05 0.06
Sampling frequency (Hz)
0.08
Fig. 5. The RMSE values for the predicted states using models fitted for the homogeneous flow-front data. The data has been perturbed around the specified sampling frequencies for presentation.
0.07
Number of lines: 5 Number of lines: 8 Number of lines: 12 Noise σ: 0.10 Noise σ: 0.05 Noise σ: 0.01
5. RESULTS This section describes the results of using the Continuous Time Stochastic Modelling package in R to estimate models for the flow-front. For all model estimates the initial states are assumed to be known. The estimation model does not take the boundary conditions of the PDE into account. Therefore the data used for estimation only include the time steps where all the flow-front progression for all measurement lines are below 0.9 m. In the homogeneous case the value for the constants in equation (5) are A = 0.0 and c0 = 3.5 · 10−4 and for the heterogeneous we have chosen A = 0.9 and c0 = 1.5 · 10−4 . These have been chosen experimentally to create a flowfront that moves the distance from 0 meter to 0.9 meter to cover the entire area in approximately 600 seconds. For both the homogeneous and the heterogeneous case we perform the model estimation with 5, 8, and 12 line sensor measurements, each at sampling frequencies of 4, 10, and 20 seconds and noise levels with a standard deviation, σ of 0.01, 0.05, and 0.1. The root-mean-square errors (RMSE) are used to evaluate the estimated models. The resulting RMSE values are calculated as the mean of the RMSE values of the onestep ahead predictions of the states at each time step and these describe how well we estimate our true system. Fig. 5 and Fig. 6 show the relation between the RMSE values and the number of measurement lines, sampling frequency and measurement noise levels for the homogeneous and the 387
0.03
0.04
The corresponding value of the log-likelihood is the observed maximum likelihood value for that data set and model. All computations were done using the free statistical software, R (version 3.3.2), and the “CTSM-Rpackage” (Continuous Time Stochastic Modelling in R version 0.6.8-5, Juhl (2015)).
0.10
0.06
where θ is a set of parameters, ZN is the set of observations, L is the dimension of the observation space, and z0 is the initial condition. For a given set of parameters and initial states, k and Rk|k−1 are computed by a continuousdiscrete extended Kalman filter. The parameter estimates are found by maximising the log-likelihood: (15) θˆ = argmax log(L(θ; ZN )|z0 ) .
Number of lines: 5 Number of lines: 8 Number of lines: 12 Noise σ: 0.10 Noise σ: 0.05 Noise σ: 0.01
0.05
RMSE
The parameters, c0,i and D0 , of the coupled SDEs formulated above to model the evolution of the flow-front dynamics are estimated from the measured data simulated using the complex PDE based model. Because of the complex structure of SDEs, estimation of parameters is nontrivial except for some simple cases. Hence we propose to use a maximum likelihood method in combination with an extended Kalman filter method described in (Kristensen et al. (2004)) to estimate the parameters. The likelihood function is formulated using the one-step prediction errors, k , and the associated variances, Rk|k−1 : L(θ; ZN ) = p(ZN |θ) (14) N 1 T −1 exp − 2 k Rk|k−1 k = p(z0 |θ) √ det(Rk|k−1 )( 2π)L k=1
RMSE
4.3 Maximum Likelihood Estimation of SDEs
387
0.05
2018 IFAC SYSID July 9-11, 2018. Stockholm, Sweden
0.05
0.10
0.25
Sampling frequency (Hz)
Fig. 6. The RMSE values for the predicted states using models fitted for the heterogeneous flow-front data. The data has been perturbed around the specified sampling frequencies for presentation. heterogeneous cases, respectively. In the homogeneous case the RMSE value is constant w.r.t the number of measurement lines, however generally decreasing with decreasing measurement noise and increasing sampling frequency. For the heterogeneous flow the RMSE is decreasing with increasing number of measurement lines and is in general not affected by changes in measurement noise and sampling frequency. Fig. 7 shows an example of the estimated flowfront trajectories (solid lines) using an estimated flow-front model and the true flow-front trajectories, with measurement noise (dot-dashed), using the PDE model. 6. CONCLUSIONS AND FUTURE RESEARCH A good understanding of the flow-front dynamics is essential for designing a monitoring system for controlling the epoxy flow inside the mould in a VARTM process. In this paper, we have proposed a simulation based methodology
Michael Nauheimer et al. / IFAC PapersOnLine 51-15 (2018) 383–388
0.2
0.4
0.6
0.8
1.0
1.2
Estimated flow Measured flow
0.0
Flow front progression length (m)
2018 IFAC SYSID 388 9-11, 2018. Stockholm, Sweden July
0.0
0.2
0.4
0.6
0.8
Flow front progression width (m)
Fig. 7. The estimated flow-front with 12 lines (states), sampling frequency = 0.25 Hz, and the measurement noise σ = 0.1. The dot-dashed lines indicate the measurements from the true system and the solid lines indicate the estimated flow-front. The line colour changes gradually from red at the first time instance to green at the last time instance. to understand and model the flow-front dynamics. A PDE based model is proposed to simulate different patterns of the flow-front dynamics for various scenarios. Furthermore, a stochastic differential equations based modelling approach is proposed to estimate a grey-box model of the evolution of the flow-front dynamics inside the mould. The advantage of coupled SDEs based approach is that it gives the flexibility of combining the partial information about the physics behind the evolution of the flow-front inside mould and the measurements from various sensors (usually pressure or optical sensors in practice) into one parsimonious dynamical model. It can properly handle insufficient knowledge of the system and noise characteristics using the diffusion term in SDEs. The proposed methodology has been validated on simulation examples. The aim of our future research is to test the methodology for more complex flow patterns (e.g. a random change in permeability (w.r.t. displacement), increasing/stochastic viscosity (w.r.t time etc.) and to investigate its validity with real experimental data. REFERENCES Alnæs, M.S., Blechta, J., Hake, J., Johansson, A., Kehlet, B., Logg, A., Richardson, C., Ring, J., Rognes, M.E., and Wells, G.N. (2015). The FEniCS project version 1.5. Archive of Numerical Software, 3(100). doi:10.11588/ ans.2015.100.20553. Dominauskas, A., Heider, D., and Gillespie, J.W. (2003). Electric time-domain reflectometry sensor for online flow sensing in liquid composite molding processing. Composites Part A: Applied Science and Manufacturing, 34(1), 67 – 74. doi:https://doi.org/10.1016/ S1359-835X(02)00232-4. Juhl, R. (2015). CTSM-R: Continuous Time Stochastic Modelling in R: User’s Guide and Reference Manual. R package version 0.6.8-5. 388
Kristensen, N.R., Madsen, H., and Jørgensen, S.B. (2004). Parameter estimation in stochastic grey-box models. Automatica, 40(2), 225 – 237. doi:https://doi.org/10. 1016/j.automatica.2003.10.001. Kueh, S.R., Parnas, R.S., and Advani, S.G. (2002). A methodology for using long-period gratings and moldfilling simulations to minimize the intrusiveness of flow sensors in liquid composite molding. Composites Science and Technology, 62(2), 311 – 327. doi:https://doi.org/ 10.1016/S0266-3538(01)00217-2. Matsuzaki, R. and Shiota, M. (2016). Data assimilation through integration of stochastic resin flow simulation with visual observation during vacuum-assisted resin transfer molding: A numerical study. Composites Part A: Applied Science and Manufacturing, 84(Supplement C), 43 – 52. doi:https://doi.org/10.1016/j.compositesa. 2016.01.006. Matuzaki, R., Seto, D., Naito, M., Todoroki, A., and Mizutani, Y. (2015). Analytical prediction of void formation in geometrically anisotropic woven fabrics during resin transfer molding. Composites Science and Technology, 107(Supplement C), 154 – 161. doi:https:// doi.org/10.1016/j.compscitech.2014.12.013. Nauheimer, M., Relan, R., Madsen, H., Thygesen, U.H., Olesen, B., and Kirkeby, K. (2017). A combined experimental and simulation based approach to model the flow-front dynamics in the vacuum assisted resin transfer moulding process. In 2nd International Conference on Lightweight Design of Materials and Engineering Structures (LIMAS 2017), 28–35. Øksendal, B. (2010). Stochastic Differential Equations: An Introduction with Applications. Universitext. Springer Berlin Heidelberg. Park, C.H., Lebel, A., Saouab, A., Br´eard, J., and Lee, W.I. (2011). Modeling and simulation of voids and saturation in liquid composite molding processes. Composites Part A: Applied Science and Manufacturing, 42(6), 658 – 668. doi:https://doi.org/10.1016/j.compositesa. 2011.02.005. Sreekumar, P., Joseph, K., Unnikrishnan, G., and Thomas, S. (2007). A comparative study on mechanical properties of sisal-leaf fibre-reinforced polyester composites prepared by resin transfer and compression moulding techniques. Composites Science and Technology, 67(3), 453 – 461. doi:https://doi.org/10.1016/j.compscitech. 2006.08.025. Stiesdal, H., Enevoldsen, P., Johansen, K., Kristensen, J., N¨ortem, M., and Winther-Jensen, M. (2006). Method for manufacturing windmill blades. EP Patent 1,310,351. Whitaker, S. (1986). Flow in porous media i: A theoretical derivation of darcy’s law. Transport in Porous Media, 1(1), 3–25. doi:10.1007/BF01036523. Yenilmez, B. and Sozer, E.M. (2009). A grid of dielectric sensors to monitor mold filling and resin cure in resin transfer molding. Composites Part A: Applied Science and Manufacturing, 40(4), 476 – 489. doi:https://doi. org/10.1016/j.compositesa.2009.01.014. Zhang, F., Cosson, B., Comas-Cardona, S., and Binetruy, C. (2011). Efficient stochastic simulation approach for rtm process with random fibrous permeability. Composites Science and Technology, 71(12), 1478 – 1485. doi: https://doi.org/10.1016/j.compscitech.2011.06.006.