Optimization of the Laser Assisted Tape Winding Process using an Inverse Kinematic-Optical-Thermal model

Laser-assisted tape winding/placement (LATW/P) is a process in which fiber-reinforced thermoplastic prepregs are heated by a laser source and in-situ consolidated by a compaction roller. Maintaining a constant temperature along the prepreg width prior to the nip point is utmost necessary to manufacture fiber-reinforced thermoplastic composites with proper bonding quality. The heat flux distribution on the incoming prepreg tape and already placed substrate is affected by the local geometry and variation in process parameters during the process. In order to maintain the processing temperature distribution at the desired conditions, a new process optimization approach for the laser power distribution is presented. A laser source with variable power distribution is considered in the optimization scheme which can improve the bonding quality of the final product by achieving the desired constant nip point temperature distribution. The variable laser power distribution can be realized practically with the Vertical-Cavity Surface-Emitting Laser (VCSEL) technology. An inverse optical and thermal process models are developed in which an ideal laser power distribution is calculated based on the surface temperature distributions on the tape and substrate described as an input to the process model. Two different optimization studies are performed based on two different input temperature profiles namely stepwise and linearly ramp profile. First, the inverse thermal model is used to calculate the required heat flux distribution which is forwarded as an input to the inverse optical model. The obtained heat flux distribution is transformed into the laser power distribution using the optical model in which the ray-tracing method is used. The relation between the surface temperature distribution and laser power distribution is investigated. The obtained laser power distribution using the inverse optical model is compared with the calculated laser power distribution as an input for the already developed optical-thermal model. 2020 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license https://creativecommons.org/licenses/by-nc-nd/4.0/) Peer-review under responsibility of the scientific committee of the 23rd International Conference on Material Forming.


Preface
I would like to thank Ismet Baran for supervising and for giving me the opportunity to work on the presented topic, Mohammad Hosseini for his guidance during the assignment and for sharing his time and effort in co-writing the papers, Amin Zaami for helping me with understanding all the theoretical aspects of the described manufacturing process. I would like to thank the members of the committee, Remko Akkerman and Semih Perdahcioǧlu, for being available during the holiday period to attend at the defence. Furthermore, I would like to thank the University of Twente, especially for allowing me to participate in the ESAFORM conference, the whole Production Technology group that gave feedback on the assignment during the PT meetings, and of course friends and family for their support throughout the Master's assignment.
iii iv Summary For years, the use of composite materials has grown in many industries. Since producing composite materials can be laborious and expensive, faster and more reliable manufacturing methods are being used. One of these manufacturing methods is Laser Assisted Tape Winding (LATW), an additive manufacturing method that uses a laser source to heat up fibre reinforced thermoplastic tapes that are compressed onto a mandrel. The fibre reinforced material consolidated during the process, requiring no curing stage. The use of robotic assistance in this process results in a fast and accurate method of producing cylindrical composite products.
Although promising, this manufacturing process does not achieve the desired bond quality between the different layers, resulting in poor performance of the final product due to varying heating conditions during the process. The limitations were shown by investigations using numerical and experimental approaches. The use of homogeneous laser sources of constant power during the LATW process does not allow for the required flexibility that is needed to stay within optimal process conditions such that a constant nip point temperature is achieved.
The introduction of vertical cavity surface emitting laser (VSCEL) allows for different optimization approaches that use custom heat distribution in the area heated by the laser. Using this laser source leads to more control over heating of fibre reinforced materials, so that ideal bonding conditions are achieved and maintained. Using an inverse approach, the ideal settings for this laser source can be determined for a given nip point temperature. However, a complete inverse model calculating these settings does not yet exist.
A three dimensional transient inverse Kinematic-Optical-Thermal model (IKOT) was constructed in order to determine the ideal laser settings. A grid-based laser source was introduced, as this laser source was assumed to be more versatile in optimization of three dimensional geometries compared to a VSCEL module. The IKOT model calculates the required laser power of each cell in this laser grid individually, based on a nip point temperature and temperature profile in the heated area.
Three case studies were performed. A steady state analysis was performed to study the accuracy and convergence of the developed IKOT model using different temperature input profiles. Power distributions were found and validated using the Kinematic-Optical-Thermal (KOT) model as the output temperatures matched the input temperature profiles.
Transient analyses were performed on hoop and helical winding processes. Desired nip point temperatures were achieved and maintained. For the hoop winding process, a decrease in total laser power was found for consecutive added layers negating the heat built up during standard hoop winding processes. For the helical winding process, in which a single layer was wound on the dome of a pressure vessel, the IKOT model was able to maintain a near constant nip point temperature whilst process conditions were changing as a result of curvatures in the substrate geometry.
In both transient analyses, different laser sources were simulated. Best results were obtained by the grid structured laser source, as it was able to achieve the desired nip point temperatures and temperature profiles, whilst requiring the lowest total laser power, compared to array structured laser sources or sources with a homogeneous laser power.
Although validated numerically, experimental results are required for exact validation. This has yet to be performed as a grid-based laser source is currently not available. However, it is expected that the technology will be available at some point in the future.

Introduction
The use of fibre reinforced Polymer (FRP) materials is constantly increasing and expanding into different industries. The advantage of the weight to strength ratio and rigidity of these materials is useful for many applications, for example in the aircraft industry in which weight saving is a critical aspect. A less common reason for using FRP materials is the behavior failure mechanics of the material. The ability of FRP products to absorb energy upon fracture results in less dangerous situations compared to traditional materials such as metals. However, manufacturing products using FRP materials can be laborious and expensive.
Although there are manufacturing methods for producing FRP products, there is always a need for faster, cheaper and more reliable methods of manufacturing. Many laborious tasks, for example the tape laying process, are taken over by robots that can produce parts with great accuracy. However, in most cases a curing stage or post process is required, slowing down the process. A process that does not require the curing stage is the Laser Assisted Tape Winding (LATW) process.

Laser assisted tape winding
The LATW process is an additive manufacturing process in which layers of fibre reinforced thermoplastic tape are placed onto a rotating substrate. An complete overview of the process is visualized in Fig. 1.1a. A close up of the nip point geometry is seen in Fig. 1  During the process, the tape and substrate are heated locally near the nip point. This is the point where the tape meets the substrate. Heat is applied by a Near Infrared laser source that melts a thin layer of thermoplastic material in front of the nip point. A rubber roller is used to firmly press the tape 3 onto the substrate, applying a constant pressure whilst the tape and substrate are cooling. The roller flattens under compression, giving a larger area of contact with the substrate. This allows the molten thermoplastic material to mix at the interface during compression so that a bond is formed during cooling.
The tape spool, compaction roller and laser source are forming the Tape Laying Head (TLH), which is fixed to a robotic arm. This arm allows for the TLH to move following the desired trajectories over the product geometry so that tape is deposited with precision, whilst a pressure perpendicular to the substrate surface is applied by the roller.
LATW is an ideal process for forming cylindrical or tubular shapes due to the rotation of the mandrel resulting in a fast coverage. Some example products that are manufactured by the LATW process are shown in Fig. 1.2. The LATW process can operate with many different thermoplastic and fibrous materials. The main advantage of this process is that the FRP tapes solidify during compression. Therefore, curing or other finishing steps are not needed, which speeds up the process [1].

Motivation
The quality of a finished product is one of the most important aspects in a manufacturing process. For this, the complete process needs to be optimal from start to finish in order to achieve the desired quality. This has yet to be achieved by the LATW process, as there are many difficulties that arise during the process [2,3]. These difficulties include: • Positioning gaps between applied tapes.
• Formation of voids during compression.
• Thermal degradation of the thermoplastic matrix.
• Not reaching the melting temperature of the thermoplastic matrix.
• Changes in laser irradiation due to changes in geometry.
• Changes in heating time due to changes in velocity.
• Heat accumulation during lay up of multiple layers.
A lot of these difficulties are related to heating the FRP materials during the process. Proper heating is one of the most important factors for achieving the correct conditions under which a proper bond can form [4,5]. A high quality LATW process thus requires full control over the heating process, yet the current process is lacking the flexibility that allows the required control.

Problem statement and goals
Most LATW processes use a laser source with homogeneous heat distribution on the geometry surfaces that are irradiated by the laser. Parameter estimation or iterative optimization methods are commonly used to optimize the the LATW process for individual process variables with use of numerical models and other tools [6,7,8,9,10,11].
For the LATW process, there is minimal optimization of the heating stage. Inverse thermal models are used in order to find ideal heating conditions for a set nip point temperature [12,13]. These patterns are approached by iterative approaches using a Vertical Cavity Surface-Emitting Laser (VCSEL) source, seen in Fig. 1.3 [14,15], since a three dimensional transient numerical model that is able to accurately predict the required laser power and settings to achieve a desired nip point temperature is not yet present. Therefore the goal of this thesis is: Development of an inverse numerical model, that is able to constantly achieve a desired temperature along the width of the nip point.
In order to make a simulation model that is competitive with models that perform similar tasks, requirements to the numerical model were set. These requirements of the inverse model are: • The model must be able to achieve a constant temperature along the width of the nip point for three dimensional geometries.
• The model must work based on any given temperature input profile.
• The model must solve the required power to deliver the constant nip point temperature.
• The model must be as fast or faster than numerical models that perform similar tasks.
• The results must be achieved for different composite materials, process properties, and laser settings.
For this, an Inverse Kinematic-Optical-Thermal (IKOT) model will be constructed that can calculate the required power distribution based on a temperature input. The IKOT model will be constructed using the three dimensional KOT model developed by Grouve [4]. The thermal model was improved by Hoeksema [16] and the optical model used was created by Reichardt [17].
As there is a huge potential in the VSCEL source Fig. 1.3, a laser source similar to this type is introduced. A grid-based laser source is used in the simulations instead of an array based source such as the VCSEL. The grid structured laser source allows for the two dimensional optimization of the heating zone required for a three dimensional transient simulation of the LATW process. The construction of the IKOT model is done according to the following steps: • Construction of an inverse thermal model.
• Validation of the inverse thermal model. a Source: Trumpf. Retrieved from: https://www.Trumpf.com/ • Integration of the grid-based laser source into the optical model.
• Determination of the approach for calculating of the required laser power.
• Validation of the model using a 2.5D simulation.
• Construction of a transient three dimensional model and validation of the complete IKOT model.

Outline
This thesis is split into three parts. In Part I, is an extended summary of the two written papers. These written papers are found in Part II. Part III contains the appendixes.

Part I
In Chapter 2 the fundamental theory required for numerical simulation of a LATW process is given. Three research cases are introduced in Chapter 3. The results of these cases are summarized in Chapter 4 followed with a discussion on the performance of the developed model. The final chapter, Chapter 5, includes the conclusion and recommendations.

Part II
Two papers were written regarding the developed numerical model. 2. Optimization of laser-assisted tape winding process using inverse kinematic-optical-thermal model The topic of the first paper is the development of a static inverse optical model as a proof of concept on the first case study. An improved numerical model is introduced in the second paper, that is able to perform the desired optimization approach for a transient three dimensional model. These papers include the the full numerical analysis performed on two case studies.

Part III
The appendix contains information regarding to changes to the optical model that were not discussed in the paper. A schematic of the process is presented as well.

Theoretical Background
Over the years, several numerical and analytical models have been developed in order to understand the physical behavior of tape winding and placement processes [4,6,17,18,19]. These models allow for optimization of many process parameters before manufacturing of the products. An schematic of a KOT model is seen in Fig. 2.1. As is visible, the numerical model is made out of several sub models, which are described in the sections below. Depending on the optimization approach and required input or output, some sub models may be left out.

Kinematic model
The kinematic model is used to simulate the movements of the TLH. A path is plotted over the surface of the simulated product. This path is split up in smaller sections representing the distance traveled during a step in time, known as time step. For each time step, the position of the laser head, compression roller and nip point are determined as well as the area on which tape is deposited. This information forms the input for the optical model.
Choosing the correct reference frame for the kinematic model is important as it determines the kinematics and fundamental thermal equations required to model the system. If the tape-roller combination is taken as a reference frame, the mandrel must rotate (Fig. 2.2a), whereas if the mandrel is used as the reference frame, the tape-roller combination is rotating around the mandrel (Fig. 2.2b). The choice of the correct reference frame depends on the goal of the optimization model.

Optical model
The optical model was required as it became clear that heat from used laser sources did not have a uniform heat distribution on the modeled surfaces [4]. Several optical models have been developed with increasing complexity in order to simulate the exact behavior of the laser light [17]. In this section, the sub models making up the optical model are described.

Ray tracing
The ray tracing procedure simulates a significant amount of light rays coming from the surface of the laser source. These rays are evenly distributed over the laser source, and have a predetermined angle with respect to this surface. The trajectory of these rays is determined and intersection points with the geometry surfaces are calculated. With this method, parallel or diverging light is simulated.
At the intersection points of the rays with the surface, the reflection angles and energy absorbed by the surface are calculated. Depending on the desired number of reflections, this process is continued for reflected rays. Rays that do not intersect with any surface are neglected. Reflections are often not used, but it is needed for calculating the reflected absorbed energy.

Laser Roller
Fibre material Ray Absorbed Energy

Reflection models
For simple optical models, a specular reflection model is used. However the surfaces of fibre reinforced materials are non-specular. This means that the incoming light does not reflect uniformly on the surface, but that it is scattered. As a result, very distinct reflection patterns are formed depending on the orientation of the fibres and the angle of incidence of the laser light at the substrate surface.
Complex reflection models have been developed that simulate the microstructure of the fibre reinforced material so that a good prediction of the reflective behavior is generated. By comparing the difference in energy of the incoming and outgoing laser rays, the energy absorbed by the fibre reinforced material is calculated.

Energy distribution
In order to calculate the spatial temperature distribution in the fibre material, the substrate is modeled by a mesh consisting of nodes over the surface and through the thickness of the tape and substrate. In order to provide the necessary heat flux input from the optical model, the absorbed energy provided by the laser light must be distributed over the surface nodes. There are several methods that can be used to distribute the absorbed energy at the ray intersection points over the thermal nodes. A brute force method such as a binning method can be applied, in which the energy of each ray is assigned to the closest surface node. A more sophisticated approach is interpolation, in which the energy is distributed over several nodes close to the intersection point. The distribution of energy is determined by the distance between the nodes at the intersection point.
A trade off is present between the two models. The interpolation method requires more calculation, but is more accurate for a lower number of used rays. If a significant amount of rays are used for simulation, the binning method is more suitable.

Thermal model
The thermal model forms the main component in simulation of the LATW process. For the most basic simulations, only this sub model is used. The thermal model calculates the temperature distribution of used composite materials for a given heat input, either supplied by a gas torch or a Near Infra Red (NIR) laser source. Every thermal model is determined by the specific case for which it is used. Different thermal models can be used to simulate the thermal behavior of the tape and substrate. In this section, several important aspects are described that form the basis needed to set up the correct numerical thermal model. The first step in producing a thermal model is determining the frame of reference and modeling approach.

Eulerian and Lagrangian approach
Selecting the correct approach is important in constructing the thermal model. This often goes hand in hand with the chosen reference frame. Two approaches can be used, namely the Eularain and Lagrangian approach.
In context of the LATW process, the Lagrangian approach follows a piece of fibre reinforced tape whilst it is heated up in front of the nip point. When reaching the nip point, the current temperature as well as the history of that piece of material is known. The Eulerian method only focuses on the nip point line, noting the temperature of a piece of material as it passes by. Only the temperature at that point in time is known.
For the substrate, the Lagrangian approach is often used in combination with the reference frame seen in Fig. 2.2b. The velocity of the TLH is already incorporated by the movement of the kinematic model. As a result the heat equation seen in Eq. 2.1 can be used, assuming that no heat is generated or dissipated within the material. In this equation, ρ is the density, C p the specific heat capacity and κ the thermal conductivity of the simulated material.
The Eulerian approach can use different reference frame for example the reference frame found in Fig. 2 If the temperature evolution of the material is important, the Lagrangian approach is the only option to use, but this requires the storage of the temperature of each node for all generated time steps, resulting in a lot of data for longer simulations.
The formulated equations are solved using numerical approximations. For both the time derivative and Laplacian operator, different approaches can be used.

Time derivative
Analytical models evaluate a steady-state solution, making the time derivative obsolete, therefore it is left out. However, for transient numerical models, this derivative must be solved. Several methods can be applied which are often based on the (Forward) Euler equation in Eq. 2.3.
The Euler method is not suitable for the heat equation as the time step is heavily constrained in order to maintain stability. An improvement to this equation is the backward (implicit) Euler equation in Eq. 2.4.
This equation is stable at the cost of extra calculation needed to estimate the temperatures at the next time step. More sophisticated and accurate methods are present, for example the Runge-Kutta method, which are often derived from the Euler method. In most cases the backward Euler method has a sufficient accuracy in approximating the solution.

Second-order derivative
There are multiple methods for solving the second order differentials. The most used method is the second-order central finite difference method seen in Eq. 2.5. This method approximates the derivative by calculating the temperature in a node based on the values of the neighbouring nodes.
These equations allow for the approximation of the internal temperature distribution of the simulated material. Although, an input is required in the form of boundary equations so that the equation can be solved.

Boundary equations
Boundary equations are required in order to solve the differential equation as they form the initial condition and input of the thermal model. These equations accommodate heat put into the system or heat loss either through convection, conduction or radiation. Conduction is the heat transfer between two different solid materials. Convection is the heat loss from solid to gas or liquid materials. Radiation heat loss is loss to the environment through infrared radiation.
The general boundary condition in the heating zone is found in Eq. 2.6, in which h is the convection coefficient, the radiation coefficient and σ the Boltzmann constant. This boundary equation is also used as the input for the cooling stage after compression. In that situation, the incoming heat of the model is equal to zero.
Heat transfer between different layers is through conduction, for which the boundary equation is seen in Eq. 2.7. This boundary equation is used between substrate and mandrel, tape and roller, and between the applied tape and substrate during compression.

Crystallinity model
The bonding quality of the finished simulation is evaluated with the crystallinity model. This model simulates healing of the used thermoplastic material in order to estimate the strength of the formed bond. Healing is the reformation of thermoplastic crystals at the bond interface. Heating and cooling cycles have a very large effect on the rate of the healing process [4]. Therefore, an optimization process that can result in the best healing conditions is highly sought after.

Case studies
Three different case studies were analyzed. In this chapter, these case studies are introduced and their significance for the analysis is discussed.

Case study 1
In the first case study, a laser assisted tape placement process was analyzed using a steady state approach. This case study was used, as it will functioned as a proof of concept of the developed IKOT model. A similar case study was used by Stokes and Compston [7]. Using this case, it was possible to validate the new model and resulting power distributions. A schematic of this process geometry is found in Fig. 3 in Paper 1.

Case study 2
A hoop winding process was the subject of the second case study. Experimental and numerical analysis were showing a significant increase in surface temperature during LATW manufacturing of a FRP pressure vessel. This behavior was also found for simple hoop winding processes [3,20], in which layers were consecutively stacked on top of each other. A schematic of this process is seen in Fig. 3  thermal properties and dimensions of the substrate. Using this case, it was possible to optimize the total laser power required in order to maintain a constant nip point temperature, negating the heat built up during hoop winding processes. Optimization was performed on three different laser types, namely a homogeneous laser source, A VCSEL source with a 28 arrays and an introduced grid based source with 308 individual laser cells.

Case study 3
In the third case study, analysis was performed on the dome of a pressure vessel. This geometry is seen in Fig. 3.1. It was shown, with experimental and numerical analysis [8], that significant temperature variations were present during the helical winding process on the dome of a gas tank. These variations were caused by changes in dome curvature that impacted the overall shape of the heated zone. This case study was selected to test the impact of these variations in geometry on the IKOT model. Two different laser types, namely a VCSEL array and grid-structured laser source were used.

Results and Discussion
In this chapter, the results of the analysed case studies are summarized. The complete case studies are found in Paper 1 and 2. After the results, a discussion found regarding the development and performance of the IKOT model.

Steady state performance
In the first case study, a steady-state process was evaluated. Focus was on the accuracy of the IKOT model using different temperature input profiles and process variables. A stepped and linear ramp temperature profile were used as the input for the IKOT model. These profiles are found in Fig.4.2a.
The resulting laser power distribution acquired with the IKOT model, found in Fig.4.1, were then used as the input in the KOT model. The validated temperatures from the KOT, seen in Fig. 4.2b, were compared to the input profiles.  As is seen in Fig. 4.2, there was a significant difference between the input and validated temperatures for the stepped profile. This profile was heavily influenced by heat flux limitations, described in Section 2.2 in Paper 1. Still, the desired nip point temperature or 500 • C was reached along the width of the nip point for both input profiles. A better solution was found when smaller mesh elements were used. The linear ramp profile was showing minor differences and was therefore used as the preferred input profile in the succeeding case studies. The convergence of the IKOT model was studied and conditions for convergence were described.

Hoop winding optimization
In the second case study, the IKOT model was used to optimize a hoop winding process. Four cycles of five layers each were simulated, using different laser grids, namely a 1 × 1, a 28 × 1, and a 28 × 11 grid. The total laser power, seen in Fig.4.3, was decreasing for all used grids. The resulting power distributions for different laser grids are found in Fig. 13 and 14 in Paper 2.  The average process temperatures per layer are found in Fig.4.4. Optimization of the 1 × 1 laser grid was resulting in large variations in surface temperatures for tape and substrate. It is seen in Fig. 4.4, that the tape surface temperatures were decreasing and substrate temperatures were increasing in every cycle. Best results were found with the 28 × 11 laser grid. A temperature close to the desired nip point was achieved whilst the lowest total laser power was required. Similar temperatures were found with the 28 × 1. However, a significant larger total laser power was required.

Helical winding optimization
[h] A helical winding process was the subject of the third case study. The numerical analysis was performed with the KOT model. An optimization analysis of the helical winding process was performed with the IKOT model using a 22 × 1 and a 22 × 11 laser grid. Temperature profiles were achieved on the substrate surface during transient analysis. The nip point temperature was achieved for different heating lengths, seen in Fig. 4.5. The resulting surface temperature histories of the tape and substrate acquired with the KOT and IKOT models are visualized in Fig.4.6. The surface temperatures found with the KOT model, depicted grey in Fig. 4.6, showed large deviations from the desired nip point temperature of 140 • C. These large deviations were not present after the optimization process. A maximum deviation of 20 • C was found, although the deviation was much lower for most of the process. Although a similar performance, less power was required on average for the 22/times11 grid source compared to the 22 × 1, as a smaller area of the laser surface was used for heating.

Computing time and accuracy
There was a difference in simulation time between the IKOT and KOT models. For a steady state analysis, the simulation time difference was not noticeable as was mentioned in the results in Paper 1, however simulation times were significantly larger for the transient analysis in the second and third case studies performed in Paper 2. For the used sub models, similar calculation times were found. Both models use the same kinematic model and almost identical optical and thermal models. The difference in calculation time was the result of the power distribution model that was added to the IKOT model. The calculation time of this sub model was dependent on two major factors, namely the number of laser cells used and the size of the mesh elements used for the inverse thermal model. The power was calculated per laser cell, using the values from the irradiated nodes by this cell. The calculation time was increasing significantly for an increase in the number of cells or mesh elements used.
The accuracy of the models used in the IKOT model, except for the power distribution model, were already proven. The inverse thermal model was validated using the thermal model it was developed from. No large deviations between these models were found. The accuracy of the power distribution model is mostly dependent on the used mesh on both substrates and on the laser source. This was clearly shown in Paper 1. Paper 2 showed that the resulting power distributions are heavily influenced by movements in between time steps during transient analysis. This is clearly visible in the results between the tape and substrate found in Section 5.2 in Paper 2. For the tape, a static mesh was used. This made it easier to match the results.

Laser reflections
Although laser reflections can be calculated by the optical model, they were not used in the case studies described in the papers. In the first model, reflections were not for simplifications, as the fist case study was a proof of concept for the IKOT model. Analysis with reflections were performed during the second case study. However, this resulted in similar or higher power distributions compared to simulations without reflections. As the analysis without reflections was giving desirable solutions it was decided not to use reflections.

Non-iterative laser power solver
The fundamental idea of the power distribution model in the IKOT model was finding the power so that the heat flux generated by the optical model (A) was matching with the heat flux generated by the inverse thermal model (B) as is seen in Fig. 2 in Paper 1. Standard iterative solvers were applied in order to find a solution to this problem. However, these solvers were not able to find a consistent and logical power distribution.
A model was set up to find the lowest power for which A was equal to B in the nodes of the thermal model. The main issue is that there were a lot more constraints compared to the number of laser cells. Furthermore, the constraint values in the node were influenced by multiple laser cells. However, there is no predetermined relation between laser cells as they were operating independently, thus the ratio of power in the constraints was also unknown. As a result, a non-uniform distribution was generated, with large variations in cell power along the length and width of the laser source and the solution was not consistent at all. For that reason, work on this method was abandoned in favour of a more simple and reliable method.

Inverse optical model
In paper 1, the developed model was named the Inverse Optical Model (IOM). Although an inverse approach was applied in terms of having a nip point temperature input and a laser power distribution output, there was no true IOM model constructed. A more suitable name would have been the Inverse Thermal Optical Model.
The optical model was required in order to know the area irradiated by the laser source. Reversing this process was not possible, as the ray tracing procedure was required, which could not be reversed. The irradiated area was required as the temperature input profiles were constructed within the heated zone. Furthermore, the normalized heat flux was also used to constraint the thermal model by limiting the maximum heat flux as no laser source is able to generate infinite power. An updated flowchart of this process is visualized in Fig. 6 found in Paper 2. A flowchart of the updated optical model with the added grid-structured laser source is found in Appendix A.

Conclusions and recommendations
In this chapter, conclusions based on the overall development and performance of the developed numerical model are made after which the recommendations are found.

Case 1:
• The inverse thermal model created is as accurate as the standard numerical model as a result of similarities between the models. The same holds for the optical model as these models were validated prior to implementation into the IKOT model.
• A laser power distribution was found using the IKOT model with different temperature patterns as an input. These distributions were validated with the KOT model and compared with the input profiles.
• Several parameters that influenced the accuracy of the solved power distributions were determined and convergence was achieved.

Case 2:
• The IKOT model was able to negate the effect of heat build up by decreasing the laser power with each added layer during the hoop winding process, whilst maintaining the desired temperature distribution on the surface. • The numerical model was able to find a solution for different given input materials and temperatures.
• The IKOT model is not suitable for optimization of a homogeneous laser source.

Case 3:
• A numerical model was created that was able to calculate the required laser power for a given nip point temperature for steady-state and transient cases.
• A desired nip point temperature and temperature profile were achieved, independent of the size of the zone irradiated by the simulated laser source.

Accuracy and model implementation
• Most potential was found for the grid structured laser source as the least amount of power was required.
• The model has difficulties calculating the required power to heat the nip point, due to influences from required power from tape and substrate in combination with low heat flux from the laser that resulted in a maximum required laser power for this area.
• The resulting power distributions were less accurate for the transient analysis. The results were difficult to achieve due to the triangular mesh used for the thermal models in combination with poorly chosen time steps.
• Most of the calculation time was required by the optical model. However, this was matched by the power distribution model during transient analysis. The kinematic, optical and thermal sub models used by both models require similar calculation times.
• For more complex simulations with reflections an iterative solver is required. However, this increases the complexity of the model and resulting computing time.
• Nip point temperatures were achieved close to the desired temperature along the width of the nip point.

Recommendations
As the optical model is a bottleneck in terms of simulation time, most improvements can be made to the total simulation time by optimizing the optical model.
Furthermore, the whole simulation process can be optimized. Simulation is performed for all time steps in the order of the sub models. For example, all the ray tracing is performed before the thermal model is used. This is useful for research purposes of individual sub models, although a new approach could simulate use the models consecutively per time step. In between, comparison with the results of the previous time step could be made. Based on smart programming, some time steps may be skipped using this approach.
The power distribution calculates the power from scratch each time step. The previous time step could be used as an initial guess to see if it fits the solution, possible skipping steps if the process is in steady state, for example a hoop winding process in which there is only a change at the start of a new layer.
A substitute mesh could be used in front of the nip point. This will give a more accurate prediction for the laser power distribution. The results can then be interpolated onto the substrate.
For inclusion of reflections, an iterative solver is the best method, for which the current model can be used for the initial guess during simulation.
If an iterative approach is used, it is wise to first solve the power of laser cells irradiating the nip point, then to solve the power required the heat further away from the nip point.
The calculated power distribution could be more accurate if the irradiated area of a node is known. This will prevent very high required powers at the edges of the area irradiated by the laser source.

Introduction
In the field of thermoplastic composites manufacturing, laser-assisted tape winding/placement (LATW/P) processes are becoming more popular due to their potential for automation. Fiber-reinforced prepreg tapes are wound around a mandrel or liner during the LATW process. A laser source is used to heat the incoming tape and already wound substrate prior to the consolidation by the compaction roller. The incoming tape and substrate are bonded at the nip point. A schematic overview of this process is depicted in Fig. 1. The in-situ consolidation of the tape and substrate at the nip point makes the LATW process relatively fast. Improper temperature distribution on the tape and substrate surfaces together with process uncertainties and variation in geometrical parameters make the process difficult to control [1,2].
Several process models have been introduced to predict the material behavior during LATW/P processes [3][4][5][6]. Recently, the optimization of the LATW/P receives more attention to develop control strategies. Optimum heat flux distribution on the tape and substrate surfaces were obtained by using inverse

Introduction
In the field of thermoplastic composites manufacturing, laser-assisted tape winding/placement (LATW/P) processes are becoming more popular due to their potential for automation. Fiber-reinforced prepreg tapes are wound around a mandrel or liner during the LATW process. A laser source is used to heat the incoming tape and already wound substrate prior to the consolidation by the compaction roller. The incoming tape and substrate are bonded at the nip point. A schematic overview of this process is depicted in Fig. 1. The in-situ consolidation of the tape and substrate at the nip point makes the LATW process relatively fast. Improper temperature distribution on the tape and substrate surfaces together with process uncertainties and variation in geometrical parameters make the process difficult to control [1,2].
Several process models have been introduced to predict the material behavior during LATW/P processes [3][4][5][6]. Recently, the optimization of the LATW/P receives more attention to develop control strategies. Optimum heat flux distribution on the tape and substrate surfaces were obtained by using inverse 2 Author name / Procedia Manufacturing 00 (2019) 000-000 thermal modeling of the manufacturing process [7][8][9][10]. The optimum heat flux distributions were obtained based on the given temperature input. The optimum distribution of the laser power was not studied in literature to the best knowledge of the authors which would require an optical model-based process optimization strategy. The recent development of Vertical-Cavity Surface-Emitting diode Laser (VCSEL) technology allowed new optimization approaches for varying laser power distribution [10][11][12]. Modules with multiple VCSEL chips allow for individual control over the chips such that custom heating patterns can be created, giving more control over the heating process in LATW/P processes. The difference between the optical power profile generated by ordinary laser sources and VSCEL sources is that the latter allows for changes in the power profile during the process. This is important as it allows for small changes in the power profile to compensate for small differences in material and geometrical anomalies during the process such that nearconstant nip point temperatures can be held. The effectiveness of VSCEL lasers is investigated by Weiler et al. [11]. Laser chips are arranged in arrays acting as a single laser which limits the adaptability of the module. Ideally, for optimal control over the laser, a grid type module is desired, which can be adapted to several heating profiles.
This paper formulates an inverse optical-thermal modeling approach to find the optimal power distribution needed to achieve a desired heat flux distribution. An inverse thermal model is first developed to generate the heat flux distribution based on the given temperature distributions. The obtained heat flux distributions are subsequently used to determine the optimum power distribution with the inverse optical model in which a ray-tracing approach is implemented. A LATP setup is used as a basis for the model implementation. Two different temperature profiles are used as the desired temperature distribution near the nip point.

Optimization approach
The optimization approach consisted of multiple models as depicted in Fig. 2 as a flow chart. The direct optical and thermal models developed in [5,13] were implemented in this work. In the direct optical model, the laser power distribution is known as well as the heat flux distribution of the output. Similarly, for the direct thermal model the surface heat flux is the input and the temperature distribution is the output. Initially, the percentage of power absorbed by the heating surfaces of tape and substrate was estimated by using the direct optical model (A in Fig. 2) with the uniform power distribution of the laser source. Subsequently, the output of the direct optical model is employed as an input for the power distribution model. The absorbed power was also used in the inverse thermal model for checking that the temperature profile was within the laser illuminated area or not. This will be necessary for winding along curved trajectories in which the heating of the substrate is not uniform. The given temperature profile was used in the inverse thermal model to calculate the required heat flux distribution (B in Fig. 2) over the same thermal nodes. The optimum laser power distribution which was defined as the output was obtained by using the power distribution model, for which the heat flux distribution obtained from the inverse model (B in Fig. 2) was the input. The input is a combination of all system settings. The output is an optimal laser power distribution to these input settings.

Optical model
The optical model used in this work was based on the already developed model in [8][9][10][11]. The system setup for numerical implementation is illustrated in Fig. 3. A laser grid was made, visible on the right side of Fig. 3, in which each rectangular cell can be customized to preferred dimensions. A ray-tracing procedure was used to calculate the intersecting points on the surfaces in three-dimensional (3D) space. This data was stored for each ray as well as the laser cell from which the ray originated. The total absorbed heat flux was the sum of the energy absorbed from all the rays. The absorbed energy of each ray was calculated as a percentage of the original cell power. The energy of each ray was distributed over the closest surface thermal nodes using bilinear interpolation as seen in Fig. 4 and a bin method, in which all the energy is collected in the closest node only. The bilinear interpolation gives a more accurate representation of the energy distribution when a minimal number of rays was used compared to the bin method. However, the bin method is computationally less expensive as the ray energy is only represented by one node on one of the surfaces. Fig. 4. Schematic of the distribution of the energy from a ray intersection point (blue point) over four corner nodes (p0, p1, p2, p3) in the thermal model. The energy (E) is distributed by dividing the total area (A=A0+A1+A2+A3) over the area opposing the nodes. To illustrate, the energy at node p2 is estimated as E×A/A2 as the node is closest to the intersection point, it has the largest amount of energy.
The difference in distribution obtained by implementing bilinear interpolation and bin method can be seen in Appendix A for a different number of rays. When a relatively large number of rays was used, both models generated almost the same heat flux distribution of the total absorbed energy on the substrate surfaces (see Fig. A1), therefore, the bin method was chosen in the present work as it reduced the total calculation time significantly.

Inverse thermal model
The required laser heat flux based on the desired temperature distribution was calculated by using the inverse thermal model. The already developed 1D thermal model in the through-thickness direction solved by an implicit finite difference (FD) scheme was used for temperature prediction [5,14,15]. The 1D model was almost as accurate as a 2D model, as the effect of heat conduction in placement and width direction was negligible [15]. A temperature profile can be generated automatically for the inverse thermal model using the heat flux distribution area calculated in the optical model. This is useful if the laser position and angle are fixed. Also, a custom temperature profile can be constructed on the surfaces for prescribed areas. The latter one is useful for determining a laser position for a set angle and allows the use of the inverse thermal model without the use of the optical model.
The flowchart of the implemented inverse thermal model is presented in Fig. 5. The inverse thermal model calculated the heat flux distribution for the tape and substrate surfaces according to the defined temperature profile. By calculating the temperatures in the material thickness direction, the boundary equation (eq. 1) was used to calculate the required surface heat flux: where , ℎ, , and ∞ were the thermal conductivity in the thickness direction, ambient heat convection coefficient, surface emissivity, Stefan-Boltzmann constant with the value of 5.67×10 -8 W/m 2 -K 4 and ambient temperature, respectively.
The inverse thermal solver was limited to a maximum allowable heat flux (qmax) of 1.67×10 6 W/m 2 [16] (Table B1 in Appendix B). When exceeded in a timestep during simulation, the direct thermal model was used to calculate the temperature using the maximum heat flux. The new temperature was used for the next timestep in the simulation.

Power distribution model
Using the required heat flux obtained from the inverse thermal model and the absorbed heat flux from the direct optical model as discussed in section 2.1 and 2.2, the optimum laser power distribution can be calculated. For each ray, the ideal corresponding laser cell power can be calculated directly by using Eq. 2. Of each ray originating from laser cell , the percentage of the total heat flux on the intersecting node ( % ( ) calculated from the direct optical model) was used to find the contribution of ray on the energy of the thermal node ( ( )). This was multiplied with the required heat flux ( ( ) calculated from the inverse thermal model), the laser cell area ( ( )) and the number of launched rays from the laser cell ( ( )) to give the power that cell needs to deliver for ray . 4 Author name / Procedia Manufacturing 00 (2019) 000-000 For laser cell , all the calculated ( ) values of all corresponding rays were averaged. The outlier values were neglected if ( ) is not within the range of standard deviation of all . Furthermore, cells irradiating on the border of the heating zone were set to zero when more of half of the total rays did not hit within the required heating zone, such that boundary conditions were clearly defined. The power of cells which irradiated outside the defined heating zone was automatically solved with this approach. These cases were treated as: However, the following criteria should be satisfied: The calculated heat flux distribution was used as an input for the direct thermal model for validation. If the calculated nip point temperature and temperature distribution were in line with the desired input values, a good optimal solution was assumed to be found.
The solution method was able to approach an inverse optical model, as inversing the ray-tracing procedure was very difficult. Calculation times were kept short by using the optical, inverse and validation model performed once to maintain a solution. This approach was only viable if no reflections were used. With reflections, multiple cells gave a heat flux on a single node on the surface. As the power ratio between laser cells was unknown, the percentage of a ray with respect to the total flux could not be determined and thus the power distribution could not be determined.

Results and discussion
As a case study, a numerical LATP setup was considered. The material used was a carbon-reinforced PEEK. The process, geometrical and material parameters used can be found in Table B1 in Appendix B. A uniform temperature of 500°C was set for the nip point along the width. This value was chosen as an arbitrary value above the melting temperature. Two temperature profiles for the inverse thermal model have been set along the placement direction and are represented in Fig. 6. The heated length of 30 and 20 mm on the substrate and tape surfaces were considered, respectively. A constant temperature before the nip point as well as across the width was desired for both profiles. In order to study the performance and accuracy of the inverse optical model three different enmeshments were considered:  Three meshes with element sizes of 1, 0.25 and 0.0625 mm 2 were used for the substrate and tape surfaces in the thermal model.  These cases were simulated with generating 1000, 10000 and 50000 rays in the optical model.  The calculated required heat flux by the inverse thermal model is depicted in Fig. 7. The trends of the heat fluxes were found to agree with the results presented in [3]. Fig. 7 clearly shows a problem for the power optimizer. As can be seen for both temperature profiles, at the nip point itself, two different ranges of heat flux were required, which mostly must be generated by a single row of laser cells. A laser cell can only be optimized to either one of two heat fluxes at the nip point for which the step profile differed a factor of ten, i.e. approximately 0.7×10 6 W/m 2 for the substrate and 0.07×10 6 W/m 2 for the tape. The same holds for the linear ramp profile although the values lie closer to each other at the nip point, i.e. approximately 0.8×10 6 W/m 2 for the substrate and 0.3×10 6 W/m 2 for the tape.
The surface temperatures calculated by the inverse thermal model are presented in Fig. 8 (refer to the flowchart of Fig.5). The stepwise temperature profile was not accurately following the intended temperature on the surface as the model was limited by the maximum laser power given in Table B1 (Appendix B). However, the required temperature was reached within a reasonable amount of distance from the start of the heating zone, acquiring and maintaining the goal temperature. The calculated temperatures exactly followed the input profile for the linear temperature case providing a more ideal and predictive heating solution as compared with the step profile.  To tackle this challenge either a temperature profile must be determined which results in an equal heat flux on the nip point for both surfaces, or no heat must be applied to this region at all. The latter found to be more favorable for obtaining the desired nip point temperature.
The resulting power distribution generated by the stepwise temperature profile is represented in Fig. 9. Although each cell was solved independently, the distribution of the power per cell was found to be similar along the width of the laser (3-9 mm region). As expected only the central laser cells were activated to cover the tape and substrate width (6.35 mm) which is smaller than the laser spot width (12 mm). The top of the laser cell, between 12 to 16 mm, mostly illuminated the tape whilst the bottom, 2 to 9 mm, illuminated the substrate. The section with low power in the middle was the intersection point between the tape and substrate i.e. nip point. Note that according to Fig. 7, both tape and substrate surfaces required an equal maximum heat flux, i.e. 1.67×10 6 W/m 2 . However, this is not directly seen in the power distribution in Fig. 9 due to the fact that the tape and substrate had different laser incident angles which resulted in different absorption behavior of the laser rays according to the Fresnel equation. More specifically, the tape had a larger incident angle than the substrate therefore maximum power for the tape was found to be lower than the substrate after the inverse optical procedure. The linear ramp profile is represented in Fig 10. The illuminated areas were the same as the step profile. The heating was spread more effectively and gradually over the area with the highest heating power in the middle of the laser. The sum of the power distribution for each mesh setting is given in Tables 1 and 2. The results show that the total laser power converged to approximately 325 W for the linear ramp profile and 395 W for the step profile. These results were best achieved by the smallest laser cells as they fit the small changes in the required heat flux in a better way.
Using 1000 rays was found to be not adequate for finding accurate results, especially when a small surface mesh was used. Sufficient amount of rays were needed such that all the nodes requiring heat flux were associated with a ray. As more rays were used, it was not needed to use the interpolation heat flux distribution on the surfaces for the laser optimization process. The use of finer surface mesh was only useful when more than 50.000 rays were used. However, the same results can be obtained with a coarse surface mesh and fewer rays resulting in less calculation time. The influence of the laser cell size was only noticeable for the 4 mm 2 cells. These tend to overestimate the power distribution as the cells were too large for representing a finer power distribution. Therefore, the optimal settings as a benchmark were the 1 mm 2 surface mesh and laser cell size, combined with the usage of 10.000 rays, as they were a good representation of the result with the lowest computing time. The obtained power distributions were used as the input for the direct thermal model to compare the input temperature profile with the optimized one. In this model, the temperature across the width at the nip point was calculated as shown in Fig. 11. It is seen that a uniform temperature on the nip point width was achieved. However, some deviation from the desired temperature (500°C) was observed. The tape achieved a temperature very close to the set temperature of 500°C however the substrate temperature was slightly higher for all the simulations averaging between 505 and 560 degrees. In some cases, this trend was the other way around, with the substrate achieving temperatures close to the set temperature and the tape temperature being too high. This corresponds to the behavior described earlier in this section, in which a single laser cell can only be optimized for either the tape or substrate heat flux near the nip point. Fig. 12 illustrates the output temperature distribution on the tape and substrate surfaces calculated by the direct thermal model by applying the optimum laser power distribution. The output temperatures closely followed the desired temperature profiles except for the substrate in the step profile. This temperature evolution was overestimated for all step temperature profile simulations with the only difference being the height of the peak temperature. With the largest laser cells, the temperature peak was found to be 560 °C and for the lowest cell size, the temperature reached a maximum of 510°C. The reason for that was the overlapping of the laser cell with the rapid drop in required heat flux around 20 mm before the nip point (Fig. 7). Using an even finer laser cell size, this behavior can be eliminated at the cost of computing time. For all simulations, the results were deviating when getting close to the nip point, giving further evidence for the behavior near the nip point due to the laser cells illuminating the nip point being optimized to either the substrate or the tape. The step profile needed more total power as the temperature must be maintained at a high level for a longer period compared to the linear ramp profile. Therefore, the linear approach was more favorable as a heating pattern since it was more accurate and required less power.

Conclusion
An inverse optical model was developed for the optimization of the nip point temperature in LATP/W processes that calculated an optimal laser power distribution. The model approach was described using a single step and linear ramp temperature profile as the input. Results showed that a desired constant nip point temperature along the width was reached. The inverse thermal model was performing fast and accurate. The inverse optical model yielded in a uniform distribution of the laser rays over the surface. The results were consistent for all simulations and the convergence of the results was established.
A simple, fast and effective method was described and used to find the optimal power distribution for multiple temperature input profiles using the inverse optical model. Optimal power distributions were found with converging results for a stepwise and linear ramp input profile. Optimal simulation settings and potential problems were determined.
Future work will include a more sophisticated solving method for the inverse optical solver and optimization of heating near the nip point. After this, more complex cases during LATW/P will be analyzed.

Appendix A
The difference in energy distribution by bilinear interpolation and binning using a minimal number of rays is depicted in Figure A1. . Energy distribution on the tape. (a) 1000 rays using bilinear distribution, (b) 10000 rays using bilinear distribution, (c) 1000 rays using bin distribution, (d) 10000 rays using bin distribution. The x-axis represents the distance from the nip point on the right side. The y-axis represents the distance from the middle of the tape.

Appendix B
Process setting parameters and material properties for carbonreinforced PEEK are depicted in Table B1.  Laser-assisted tape winding

Introduction
The laser-assisted tape winding (LATW) process is gaining more traction as a highly-automated manufacturing method for tubular fiber-reinforced thermoplastic composite products. The LATW process, in which a rotating mandrel or liner is utilized, allows the production of large parts without the need for a curing step, enhancing the production rate. Similar to the laser-assisted tape placement 5 (LATP) process, in which the tooling is stationary, the incoming pre-impregnated (prepreg) tape and substrate are heated to the melting temperature with the use of a diode laser before touching each Preprint submitted to Elsevier August 13, 2020  other at the 'nip point'. The tape is compressed onto a mandrel or liner substrate by means of a deformable compaction roller to establish the in-situ consolidation which is the main mechanism to form the final product as seen in Fig. 1. The process temperature, speed, compaction pressure 10 [1, 2, 3, 4], and fiber path trajectory [5,6,7,8] are the important factors for developing a proper bonding quality at the interface between the tape and substrate to manufacture defect-free products with high mechanical performance. In addition, optimal ranges of the tape winding and placement parameters such as tape tension and roller force were provided to reduce the void content [9, 10] and residual stress [11,12] and to increase the tensile strength [9,11], interlaminar shear strength [13,14], 15 peel strength [15,16], and production rate [14,15,12] The nip point temperature as an indicator for the process temperature might vary during the LATW process resulting in a substantial variation in the consolidation quality. As an instance, the nip point temperature increased significantly with the increase of substrate thickness during the continuous multi-layer hoop winding process due to the heat accumulation in the the previously wound layers 20 [17,18]. Besides, an increase in the winding angle resulted in higher temperatures for the tape and lower ones for the substrate during LATW process in [19]. Maintaining the nip point temperature within the processing temperature is a key challenge for the current LATW technology. Therefore, efforts were performed to optimize the process settings. An inverse thermal model was developed to calculate the required heat fluxes for a target temperature profile on the substrate surface for an LATP 25 process with a flat tooling [20]. In another attempt for tape lay-up process [21], the unknown heat transfer coefficients were determined by achieving the consistency between the measured and predicted temperatures which is known as inverse heat conduction problems (IHCPs). In [22], the laser power was adjusted in order to ensure a constant target temperature for an AFP process on a flat tooling using the proper generalized decomposition (PGD) method. 30 The current LATW technology is mostly limited by the use of laser sources that mainly generate uniform power distribution. The introduction of the vertical-cavity surface-emitting laser (VCSEL) paved the road for optimization procedures by allowing real-time control in the laser power distribution 2 [23,24,25]. To this end, there has been a limited research focusing on finding the optimum laser power distribution during LATW and LATP processes. An inverse analytical thermal model was developed to 35 tailor the power of the individual emitters of VCSEL which was validated by a simplified test setup in [24]. The transient total power of VCSEL for a uniform process temperature was introduced during the AFP process of a single-curved tooling [23]. In a recent study [26],  [27]. Two case studies are considered based on the previous works which described the temperature variation during i) multi-layer hoop winding [18] and ii) single-layer helical winding of a doubly-curved pressure vessel [27]. For the former case, the presented in [27,18].

Problem definition
In this work, two different winding cases are defined to be optimized: i) multi-layer hoop winding process on a cylindrical mandrel [18] and ii) helical winding process on the dome section of a pressure vessel [27]. The variation of the nip point temperature in these winding cases were studied numerically 60 and experimentally.
The evolution of the process temperature was predicted by using the KOT model in which the kinematics of the process was coupled with the optical-thermal model. The kinematics of the winding process presented in [18]    time-dependent process settings [27].
The hoop winding process was performed in 4 cycles of 5 continuous layers starting with a 5-layers-thick substrate. The used unidirectional carbon fibre reinforced prepreg tape was TC1200 AS-4/PEEK provided by TenCate. The material properties and process settings can be found in [18].
In order to reduce the computational cost, only a sixth of the cylinder circumference was simulated (see  can be found in [18]. Regarding the helical winding of the pressure vessel's dome section, the CELSTRAN ® glass fiber GF70-01/ high-density polyethylene (HDPE) [29] prepreg material was used. The material properties and process settings can be found in [27]. Fig. 4 represents the history of the substrate curvature, 105 winding angle, and the velocity of the tape laying head during the helical winding. Temperatures of the tape and substrate were acquired at 5 mm prior to the nip point are shown in Fig. 5. The time-dependent process velocity presented in [27] was neglected to simplify the helical process even further with a constant velocity of 50 mm s −1 . Therefore, the results in Fig. 5 are not the same as the results shown in [27] due to the applied simplifications, yet following the same trend.   A constant temperature of 125 • C on the tape surface was achieved whereas large fluctuations in the substrate temperature were found mainly due to the change in the local tooling curvature. The maximum temperature on the substrate was slightly over 225 • C mainly due to the lower incident of the laser beams leading to a higher energy absorption according to the Fresnel law. In order to minimize the variations in the process temperature and keep the process temperature within a target 115 temperature range, a process optimization framework is developed based on the inverse KOT model which is explained in the following section.

Inverse Kinematic-Optical-Thermal (IKOT) model
A new inverse kinematic-optical-thermal (IKOT) model was developed by coupling the KOT model with the inverse optical-thermal (IOT) model introduced in [26]. The IOT model developed in [26] was The kinematic model provided the history of substrate thickness distribution (A1) and the tape laying head (TLH) position with respect to a surface mesh (A2). The TLH consisted of the laser source, compaction roller, and incoming tape which followed a given input fiber path discretized according to the distance that TLH can travel in 3D space depending on the time-dependent tape feeding rate 130 during the constant time step (∆t). For each time step, the TLH location (A1) and the thickness distribution (A2) were passed into the optical and the IT model, respectively. A normalized heat flux distribution (q f w (X, t)) on the tape and substrate surfaces (B) was calculated using the optical model as function of time (t) and the 3D space (X). The thickness distribution (A2) and theq f w (X, t) (B) were used in the IT model to determine the required heat flux distribution (q iv (X, t)) on the tape and 135 substrate surfaces (C) needed to maintain the nip point temperature as close as possible to the target temperature. Finally, the magnitude of the required heat flux (C) and the normalized heat flux distribution (B) were used to calculate the optimized laser power distribution by using the laser power distribution (LPD) model. The details of the IT mode and LPD are explained in the following sections.

140
The tape and substrate temperature distributions in the KOT model were obtained by using two independent thermal models which were coupled at the nip point. Multiple slices of 2D thermal domains were considered along the winding direction and through the thickness direction for the tape while multiple triangulated through-thickness 1D thermal domains were used at the substrate mesh locations as described in [27]. In the thermal model the objective was to calculate the temperature 145 evolution as the output, while in the IT model the required heat flux distribution was calculated on the heating surface to achieve the target temperatures at the tape and substrate nip point (T nip target ).
As discussed in section 2, the substrate temperature variations in the helical process was due to the changing heating flux distribution e.g. heating length, due to the local variation in the surface curvature. Therefore, it was not trivial to define a unique target temperature profile on the substrate 150 surface. To this end, linear temperature profiles for each of the triangulated 1D substrate thermal domains on the heating surface were defined. No matter how the temperature profile evolved, the aim was to reach T nip target at the nip point (t nip ). To illustrate the variation in the required linear target temperature profile, an example is given in Fig. 7 in which the fibre path orientation was changing in one case and was constant in another case. A schematic view of heating zones at three time stamps (t1 155 to t3) is represented by the highlighted areas. The inner and outer edge of the fibre path trajectory depict the region on the substrate on which the tape was deposited. For both cases, the top surface of the highlighted thermal domains i.e. triangular elements P1 and P2, were located on the nip line at t = t3 at the inner and outer edges, respectively. Initially, all the thermal domains had the same temperature of T i target at t = t i . At t = t1, the element P1 on the inner edge was out of the heating 160 zone while element P2 was heated as seen in Fig. 7. while both elements were heated during t2 to t3.
Besides, P1 was partially heated and P2 was fully heated. Therefore, element P1 had a different heat flux history with a shorter heating duration (t2 to t3) to achieve the T nip target compared to element P2 (t1 to t3) which was on the outer edge. As a result, a higher temperature gradient and therefore higher laser power was required for the inner edge of the curved trajectory to maintain an even temperature 165 along the width of the nip point. On the other hand, both P1 and P2 elements were illuminated for the same amount of time from t1 to t3 in the case of constant fiber path orientation. Therefore, the linear target temperature profile for both points were the same as depicted in Fig. 7. The linear target temperature profiles for the heated substrate thermal domains can then be interpreted as a time-dependent target temperature distribution on the mesh surface such that each triangular element 170 on the top of 1D thermal domain follows the defined target temperature profile during the process.
The linear target temperature profile was determined by using Eq. 1 as also visualized in Fig. 7. Eq. 1 explains the slopes/gradient definition of the linear temperature profiles.
where T nip target was the target temperature at the nip point, T i target was the initial temperature at the beginning of the heating zone at t = t i i.e. t i = t1 and t i = t2 for P2 and P1, respectively. t nip was the time at which the nip point was located on the desired triangular element i.e. end of the heating zone. ∆m = t nip − t i was the heating duration to achieve T nip target which depended on the transient heat flux 175 distribution. In each time step, the illuminated and un-illuminated triangular elements of the substrate were known by the KOT model which used to calculated the ∆m. It was assumed that the T nip target was achieved one time step (∆t) prior to the nip point as was done also in [26] and shown in Fig. 7.
Once the target temperature profile for the substrate was defined on the top surface (T target ), the required heat flux profile to reach T target was calculated based on the heat transfer governing equation  for 1D thermal domains shown in Fig. 9.
where ρ, c p and k z are the density, specific heat capacity, and thermal conductivity in the thickness direction, respectively. Eq. 2 was solved using implicit control volume based finite difference (CV/FD) 180 scheme which was implemented in MATLAB 2019b. There were three general types of nodes in the 1D domain i.e., the interior nodes and two boundary nodes. The finite difference expressions of these nodes are given in Eq. 3 for the interior nodes and in Eqs. 4 and 5 for the boundary nodes.
where F o = α · ∆t/∆z was the dimensionless Fourier number, α = k z /(ρc p ) was the diffusivity of the material, and Bi = h · ∆z/k z was the dimensionless Biot number. ∆z was the mesh size in thickness 185 9 direction. h a and h t were the heat transfer coefficients for the surfaces in contact with air and tooling at the corresponding temperature of T a (for cv = 1) and T t (for cv = K), respectively. T cv=1 was the temperature of the control volume at the heating surface and T cv=K was the temperature at the substrate-tooling interface. q n+1 iv was the unknown required absorbed heat flux on the surface induced by the laser source.

190
In the IT model the condition T n+1 1 = T n+1 target was defined as the temperature at the top surface should follow the target temperature profile. In order to find unknown q n+1 iv , therefore Eq. 4 is rewritten for the cv = 1 as: The linear system of equations out of the discretized governing equations is shown in Eq. 7 using Eqs. 3, 5, and 6.
with matrix d equal to the right hand side of Eqs. 3 to 5.
The temperature distribution through the thickness at t = n + 1 and the required heat flux were calculated using Eq. 7.
It should be noted that the maximum power a laser cell can generate was limited. As a result, the maximum allowable heat flux (q max ) was also restricted to be delivered to a thermal node. If the 195 calculated q n+1 iv exceeded the q max , another calculation step is required to find the maximum allowable temperature distribution using Eq. 7 and q max . In this case the T n+1 1 was less than T n+1 target , and therefore it is needed to be compensated in the next time step.
The flowchart of the IT model for the substrate is shown in the Fig. 8.
For the tape, the same IT model as the [26] was used since the temperature profile reached to a 200 steady-state phase due to the constant TLH geometry and the consequent constant heating zone during the process. The target temperature profile was then constant during the process and considered as a linear profile from the room temperature at the beginning of the heating zone to T nip target at the nip point as shown also in Fig. 7.

205
The LPD model was based on the arrays of laser cells which are independently controllable to optimize the time-dependent laser power distributions during the process as shown in Fig. 9 for an exemplary case. The reflections of the laser rays were neglected to simplify the LDP model. As seen in Fig. 6, the inputs of LDP model was the time-dependent absorbed heat flux which was normalized with respect to laser power intensity (q f w (X, t)) and the required absorbed heat flux (q iv (X, t)) which 210 were generated by the KOT and IT models, respectively.
The normalized absorbed surface heat flux (q j f w (i, n)) by the thermal domain 'i' which was illuminated by laser cell 'j' for the time stamp 'n' is schematically shown in Fig. 9. It is seen that also the laser cell 'j+1' had a contribution to the heat flux of thermal top surface of the domain 'i'. Therefore, the total heat flux on the substrate thermal node 'i' (q f w (i, n)) was determined by summing up all the contributing laser cells as described in Eq. 8 which was the procedure for the modified KOT model.
with P j equaled to the laser cell power, A j was the cell area, and N c was the total number of laser cells heating the thermal node 'i'.
The power of each laser cell was determined based on the required heat flux (q iv (i, n)) of the illuminated triangulated elements on the substrate. In the Fig. 9, up to 8 triangular thermal elements were illuminated on the substrate surface by the laser cell 'j' as an exemplary case which depends on the substrate surface mesh and laser cell sizes. As a result, a weight-averaging approach was consider to estimated the required laser power in Eq. 9.
Nc j=1q j f w (i, n) (9) in which N t was the total number of triangular thermal elements illuminated by the laser cell 'j' and R j (i, n) was the ratio of normalized heat flux from laser cell 'j' divided by the total normalized heat 215 flux for the top surface of thermal domain 'i'.
By increasing the area of the laser cell, more triangular elements were illuminated by a single laser cell resulting in larger number of data to be weight-averaged in Eq. 9. Therefore, the exact required heat flux by each thermal elements was subjected to deviation since the required heat flux of other illuminated thermal elements should also be satisfied by the common laser cell. Therefore, illumination 220 of one thermal element by one laser cell provides the optimum condition foe IKOT model ideally. The major advantage of the IKOT optimization technique is the small added computational cost to KOT model. This is because the optical and kinematical models were the same for IKOT and KOT models.
The IT solver had the same size of system of equations to be solved as the thermal model in KOT model. Finally, the extra LPD model was not computationally expensive since only the analytical Eq. 225 9 needed to be solved. However, the gradient-based optimization approaches in which the initial guess converges to the target temperature needs to iteratively solve the optical and thermal models requiring much higher computational cost.

Absorption ratio
The laser cell absorption ratio is introduced to calculate the efficiency of the illuminated surfaces to absorb the irradiated heat by the laser source with uniform power distribution. Besides, it helps to identify the laser cells responsible for heating the tape or substrate. The absorption ratio (R j abs ) is defined by Eq. 10.
The I k,j abs term was the absorbed energy of a the ray 'k' launched from laser cell 'j' as also indicated in 230 Fig. 9. This term was calculated by the KOT model considering the Fresnel equation and the geometry of the setup which determined the incident angle of the laser rays. I 0 was the initial energy of the ray 'k' which was the same for all the rays due to the uniform power distribution pf the laser source, and N r was the number of rays coming from the laser cell 'j'.

Process optimization studies 235
The developed IKOT model was used to optimized the laser power distribution for the cases studied explained in Section 2.
• Case-1: Three laser grids settings were considered to evaluate the performance of the IKOT model for multi-layer hoop winding process, namely a 1 × 1, 28 × 1 and 28 × 11 grids. As the laser spot size was 28 mm × 11 mm, the individual laser cells sizes were 1 mm × 11 mm and 240 1 mm × 1 mm for 28 × 1 and 28 × 11 grids, respectively. The target nip point temperature for the substrate and tape was determined as 400 • C for the hoop winding of C/PEEK tapes [30]. The  initial temperature of the layer 5 was higher than the layer 1 due to the accumulated heat during winding of layers 1 to 4 which was obtained from the KOT model.

Hoop winding
The distribution of the laser cell absorption ratio is plotted on the laser surface for 28 × 1 and 28 × 11 270 laser grids in Fig. 12. The laser cell rows 1 to 11 and 13 to 28 were considered to irradiate the tape and substrate, respectively. For the 28 × 11 grid, the first and last three columns of the tape region were pointing to the roller. The nip point geometry which was heated by the laser row 12 had the least ability to absorb the laser heat (R abs = 0.59) due to the lower incident angles of laser rays. The tape had lower R abs = 0.82 on average compare to the substrate R abs = 0.87. The roller, however, had the highest R abs = 0.96 due to the lower refractive index compare to the composite surface. The absorption ratio for the 1 × 1 grid was 0.82 which included all the tape, substrate, and roller components.
The optimized laser power distributions for 28 × 1 and 28 × 11 laser grids required to achieve the target temperature profile during hoop winding are shown in Figs. 13 and 14, respectively. The laser power distribution of the area which heated the tape i.e. rows 1 to 11, remained the same for each 280 added layer as expected. The power of laser cells in columns 1, 2, 10, and 11 of the tape region was zero for the 28 × 11 laser grid which means that the corresponding laser cells were not active. By excluding the deactivated laser cells, the average power (intensity) of all laser cells were 12.8 W (1.16 W mm −2 ) and 1.17 W (1.17 W mm −2 ) for 28 × 1 and 28 × 11 laser grids, respectively. It is seen that the power of the cells that were used to heat the substrate was gradually decreasing with each added layer to 285 compensate for the accumulated heat in the substrate. The non-uniform power distribution across the width of the laser source for 28 × 11 laser grid was due to the unstructured substrate mesh elements with respect to the laser source as depicted earlier in Fig. 9. The cells pointing at the nip points of the tape and substrate i.e. row 12 had the largest power to compensate for the relatively low absorption ratio up to 45 and 4.4 W for 28 × 1 and 28 × 11 laser grids, respectively. The heat flux and temperature history of a single thermal domain on the substrate in the steady-state phase is depicted in Fig. 15 for layers 1 and 5. The optimized laser power distribution was utilized as the input for the KOT model to generate the optimized heat flux and temperature histories.
Although only 5 time steps were considered in the heating region a good match between the target temperature profile and optimized results generated by the IKOT model was achieved. It should be 295 noted that the IT results were actually the same as the target temperature profile. On the other hand, the obtained temperatures from the IKOT model at the nip point, i.e. at 0.32 s, were found to be higher than the target temperature by 19 • C (4.7% deviation) and 25 • C (6.2% deviation) for layers 1 and 5. This relatively small deviation was mainly due to the mismatch in the used grid distribution and the triangular thermal domains.

300
The optimized heat flux and temperature distributions along the centerline of the substrate at the steady-state phase during winding of layer 1 and 5 are plotted in Fig. 16 for 28 × 11 laser grid. In addition, the heat flux and temperature distributions obtained from the IT model are also shown in Fig.   16 in order to evaluate the effectiveness of the IKOT model. Note that the results illustrated in Fig. 16 were extracted from the triangular elements at the top surface of the substrate at a specific process 305 time, i.e. the time at which the nip point temperature reached steady-state phase [27]. A step-wise temperature distribution was predicted because of the relatively large time step of 0.08 s which resulted in the heat flux movement of 8 mm per time step (v × ∆t = 100 mm s −1 × 0.08 s). This is different than  the target temperature profile in Fig. 10 which was defined in time domain for each triangular thermal element whereas the temperature distribution was obtained at a specific time stamp of the process.

310
In general, a very good match between the obtained heat flux from the IKOT model by using the optimized laser power distribution and the required heat flux from the IT model (q iv (X)) was achieved.
However, the IKOT heat fluxes were slightly larger than IT results around 23, 15, and 8 mm prior to the nip point which led to the overestimation of the temperature distribution as a consequence as seen in Fig 16b. The required heat flux was increasing in the range of 33 to 8 mm away from the nip point 315 up to 1.1 W mm −2 and 0.78 W mm −2 for layers 1 and 5, respectively. However, a decreasing trend was predicted in the region of 8-0 mm since the target temperature profile was in the plateau region. This trend was also observed in Figs. 13 and 14 in which the laser power distribution reached to a peak at row 17 and then the required laser power reduced for rows 15 to 13 just before the nip point.
The target temperature profile was also followed by the IT and IKOT models with relatively good 320 match. However, around the nip point and 29, 21, 13 mm away from the nip point the IKOT temperature was higher than IT model due to the larger heat fluxes. the optimized substrate temperature during consecutive winding of 5 layers using the 28 × 1 and 28 × 11 laser grids varied in the range of 410 • C to 420 • C (up to 5% mismatch) and 405 • C and 395 • C (up to 1.25% mismatch), respectively. Therefore, the 28 × 11 laser grid provided a process temperature evolution closer to the target temperature than the 28 × 1 laser grid.

335
The evolution of the optimized total power for each layer which was obtained by summing up the power of each optimized laser cells is depicted in Fig. 18. The required total power to maintain the desired nip point temperature showed a similar pattern for all the laser grids.
The total power decreased during consecutive winding of 5 layers. The lowest average total power was obtained for the 28 × 11 laser grid configuration as approximately 300 W which was 25% lower 340 than the total power of the non-optimized case. The average total power for the 28 × 1 and 1 × 1 laser grid configurations was found to be relatively higher as compared with the 28 × 11 laser grid. Therefore, 28 × 11 laser grid was found to be effective based on the obtained process temperature evolution and the required average total laser power.

345
The 22 × 1 and 22 × 11 laser grids were used to optimize the helical winding process of the dome section of the pressure vessel. The 1 × 1 laser grid was not used here since it did not converge to the desired process temperature in the hoop winding optimization. The distributions of the absorption ratio for the laser sources are shown in Fig. 19 as an exemplary case when the nip pint was located at facet C at 6.5 s of the winding time.   125 • C with a 2% deviation as compared to optimized temperature of 238 • C.
The resulting total laser power for the tape and substrate are shown in Fig. 25. The corresponding laser cells heating the tape and substrate were determined based on the Fig. 19. Both laser grids 400 resulted in a more or less similar pattern for the total power. The tape power was almost constant at 204 W and 93 W for the 22 × 1 and 22 × 11 laser grids, respectively. However, the required power for heating the substrate was fluctuating due to the changing substrate curvature and winding angle. It is noteworthy that for both laser grids the tape needed higher power as compared with the substrate as the non-optimized tape temperature increased to reach to target temperature while the substrate 405 temperature reduced to be optimized. Fig. 25 shows that 22 × 11 laser grid required an averaged power of 160 W which was 240 W lower than the non-optimized case. The required power for 22 × 11 laser grid was also half of 22 × 1 laser grid, however, both configurations resulted in almost the same optimized temperature history according to Fig. 24. predicted lower power consumption while resulting almost the same optimized temperature output.
The current IKOT dismissed the reflections of the laser rays which is needed to be considered for a more accurate laser power distribution. The overall resulting laser power is expected to reduce due to 430 the absorption from the reflections by the opposite surface. Furthermore, the IKOT model needs to be 24 validated with actual experimental tests by using the VCSEL heat source for the multi-layer hoop winding process as the future work. Finally, to improve the accuracy of the IKOT model a finer substrate mesh and laser cell size at the vicinity of the nip point are suggested to limit the deviations from the target temperature.

435
Appendix A

Changes made to the optical model
In the two written papers, the modification to the thermal models and the new power distribution model were discussed. These were not the only changed models. Changes to the optical model were made and discussed, yet on a superficial level. Changes were made to the model developed by Reichardt [17]. A schematic of the new optical model is presented in Fig. A  Most functions of the optical model were unchanged. A grid structured laser source was introduced. The new source was functioning the same as a normal laser source. For the input variables, the total dimensions of the laser source were kept. Extra parameters were added that divides the laser source in rows and columns, which was done in the construct laser grid sub model. Each cell of the laser grid was treated as a single laser source, with a set number of rays and an equal ray distribution. The total number of rays used as the input was divided over the number of laser cells and rounded to integers. The ray distribution was only calculated once. These values were used for all laser cells.
The standard input variables, for example divergence could be used for each cell individually. This was also true for other advanced options such as the complex surface reflection models. The model was further altered so that for each ray the cell of origin was known and stored in memory. This data was used to construct the output of the optical model, which was the normalized heat flux per laser cell, for each node on the irradiated surfaces. Furthermore, the model was slightly altered such that the results could be visualized in the developed model. Calculations on ray tracing and reflection model were calculated per generated ray as usual.