Previousresearchhashintedonfurtherimprovementsofthebucklingbehaviourofvariable-stiﬀness laminates by incorporating overlaps, resulting in a variable thickness proﬁle that is non-linearly cou- pled to the steering angles. The present study compares two modelling strategies to consider the variable thickness distribution: 1) as-manufactured with discrete thicknesses; and 2) smoothed with a continuous thickness distribution. The as-manufactured discrete thickness created by overlapping tows is obtained by means of virtually manufactured laminates. The smeared approximation is much simplertoimplement,wherebythelocalthicknessisanon-linearfunctionofthelocalsteeringangle. Linear buckling analyses are performed by means of fast semi-analytical models based on the Ritz methodusinghierarchicalpolynomialsandclassicalplateformulation. Byassumingasmoothmanu- facturingmouldononeside,aone-sidedthicknessvariationisproduced,resultinginnon-symmetric laminates for which the mid-plane surface is varied accordingly. Modelling guidelines are provided regarding the use of the smeared model in a study covering a wide range of geometries, loading and boundary conditions. With these guidelines, one can apply the smeared thickness technique in semi- analytical models to reach a correlation within ,5~compared to a costly discrete-thickness ﬁnite element model. 1. Introduction In the aerospace industry, weight reduction has always been an important goal in order to increase eﬃciency. This hascausedashiftinthematerialsfromwoodtometal,until latelywithcomposites,suchasintheBoeing787orAirbus A350. Owing to their speciﬁc strength and stiﬀness, com- bined with their mechanical tailoring properties, composite materials have been employed progressively more. The lat- est composite developments involve the emergence of Tow- Steered Laminates (TSL), also called variable-stiﬀness or Variable-Angle Tow (VAT) laminates, increasing structural capabilities by incorporating a higher degree of tailoring. Earlyuseoflocalﬁbresteeringwasperformedinorderto increaseopenholestrength[47]. Lateron,Biggersetal.[3] strategicallypositioneduni-directionalstripsinordertoim- prove the buckling capabilities of composite plates, locally varyingthestiﬀnessandbeneﬁtingfromloadredistribution. The addition of these strips resulted in Variable Stiﬀness Laminates (VSL). Building upon this idea of varying the stiﬀness distribution locally, VAT, also called steered-tow orﬁbre-steeredlaminateswerecreated. VATlaminatescon- sist in the spatial variation of the ﬁbre orientations within a ply, as opposed to straight ﬁbres used in conventional Con- stant Stiﬀness (CS) laminates. This diﬀerence is illustrated in Figure 1. Since, the beneﬁts of VAT have been used and re- searchednumericallyintoincreasingmechanicalproperties, such as stiﬀness [19], buckling [18, 20, 29, 46] and post- buckling [45, 44, 36] capabilities for laminates alone, or <Corresponding authors, lander.vertonghen@onera.fr, S.G.P.Castro@tudelft.nl ORCID(s): 0000-0002-3915-6323 (L. Vertonghen); 0000-0001-9711-0991 (S.G.P. Castro) Figure 1: Conventionalconstant-stiﬀness(CS)laminateonthe left, and variable-angle tow (VAT) laminate on the right [7]. along the interaction of stiﬀeners in a bay design [15, 33] and linear buckling of VAT ﬁlament-wound cylinders [40]. These beneﬁts have been proven in experiments to increase propertiesindiﬀerentapplications,suchasbending[38]and buckling for cylinders [24, 25] and plates [27, 30]. Theproductionofﬁbresteeringrequiresautomatedman- ufacturingmethods. AutomatedFibrePlacement(AFP)isa consolidated process, largely used in the manufacturing of compositeparts. Theroboticarmscanbeconﬁguredtobend thetowswhilefollowingareferencecurvedpath[6]. Novel manufacturingtechniquessuchasContinuousTowShearing (CTS) produces curved ﬁbres by means of shearing [23]. DuringtheAFPprocess,overlapsappearduetothenon- parallel path course of variable tow angle, as illustrated in Figure 2. If these overlaps are to be avoided to keep a more uniform thickness, the distance between two subse- quenttowscanbeincreasedtoobtaingapsinstead[4]. Such gaps become resin rich pockets and may become predomi- nant locations for damage and failure [26]. Subsequent cut- and-restart of adjacent tows is another possibility to mini- mize the presence of overlaps and gaps, at the price of in- creasing the manufacturing costs [31]. Instead of gaps or overlaps, the CTS method produces a continuously chang- ing thickness proﬁle [23]. L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 1 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 2: VAT manufacturing results: AFP with gaps (a) and AFP with overlaps (b) [5]. When tows are bent in the AFP manufacturing process, ﬁbres are prone to wrinkling and twisting [2, 7], as a di- rect result of compressive strains applied on the inner side of the curved ﬁbres. To prevent this negative eﬀect, it is a commonpractisetoconstraintthepathcurvature,withearly assumptions using a minimum radius of 635 mm[28], re- gardlessofthematerialpropertiesormanufacturingparam- eters;whereasrecentstudiesclearlystatethatthemaximum curvature depends on deposition rate, temperature and ma- terials [14]. Clancy et al. [14] used laser-assisted AFP of thermoplastictowsandpresentedprocessingparametersthat allowed a minimum steering radius as low as 400 mm. On theotherhand,theCTSprocessallowsaminimumsteering radius of 50 mm[23]. The presence of gaps and overlaps inﬂuences the ideal response of a VAT laminate. The eﬀects of gaps has been incorporated by the defect layer method, where the stiﬀ- ness matrix is reduced based on the relative gap percentage in a Finite Element Method (FEM) element [16]. Instead, Mishraetal.[31]adaptedthisdefectlayermethodbymeans of a homogenization of the stiﬀness matrix as a function of the ﬁbre angular distortion. However, despite showing higher strength and stiﬀness than conventional laminates, gap designs are less stiﬀ and have a lower ﬁrst-ply failure than the overlap design counterparts [27], even when mass- normalized,inbending[38]andbucklingexperiments[30]. However,littleliteraturehasinvestigatednumericallytheef- fects of overlaps on the buckling behaviour of VAT lami- nates. Wu[41]analyzedtheeﬀectofavaryingthicknessproﬁle due to overlaps created by AFP on a cylinder. After repli- cating the tow courses and overlap locations, as shown in Figure 3, Wu used this information in a FEM analysis. He deﬁnedeachelement’sstiﬀnesswiththethicknessandlayup informationofthevirtualmodelclosesttotheelement’scen- troid and perform a non-linear buckling analysis. This ap- proach showed to be in good comparison with experiments of the manufactured cylinder with overlaps [42]. Upon the variable thickness phenomena, Castro et al.[10] proposed a smeared thickness formulation to ap- Figure 3: Discretised laminate thickness [41]. proximate the local thickness due to overlaps in AFP- manufactured laminates. The formulation is based on the towwidthapproximationfromBlometal.[5],whilstapply- ing a constant-volume constraint. Groh and Weaver [18] investigated the inﬂuence of a varying thickness proﬁle in CTS manufactured laminates, modellingthevariationbothasa2Dplateandasa3Dcurved shell, to assess which one represents best the 3D behaviour determinedbymeansofamodelusinglinearbrickelements. They concluded that a curved shell formulation represents better the thickness inﬂuence, and also saw the possibility of using double curved thickness variation to increase the buckling performance. Alternatively,Irisarrietal.[22]andPeetersandAbdalla [35] used variable thickness by means of ply drops along curvedﬁbrestoretrievethestackingsequencesafteralinear buckling optimization using lamination parameters. Thepresentstudyfocusesontheuseofcontinuoustows in AFP manufacturing allowing the formation of overlaps. Virtually-manufactured laminates are used to obtain the lo- cationsoftheoverlapsandlayer-wisethicknessdistribution, according to a given ﬁbre orientation and tow width. The resultinginformationisadiscretethicknessdistributionthat canbeusedeitherinaFEMorasemi-analyticalmodel. The potential of the smeared thickness approximation by Castro et al. [10] is assessed for representing the inﬂuence of over- laps on the buckling behaviour of the AFP laminates herein studied, and compared with the modelling approach using discrete overlaps. Finally, guidelines on the use of the sim- plersmearedthicknessmodellingaregivenforawiderange of design possibilities and manufacturing parameters. 2. Virtual Manufacturing This section addresses the idealization of variable stiﬀ- nessbydeﬁningtheﬁbrevariation,followedbyspecialcon- siderations regarding the overlaps. Finally, the procedure to obtain the discrete thickness distribution of a laminate and the way to incorporate manufacturing constraints are dis- cussed. L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 2 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 4: Global.x;y/coordinate system and manufacturing .x¨;y¨/coordinate system rotated by the angle . Fibre angles are deﬁned using the angle , measured about the manufac- turing coordinates. 2.1. Fibre Steering Two coordinate systems are introduced: a global and a manufacturing one; both illustrated in Figure 4. The man- ufacturing coordinate is rotated about the global coordinate system through the ply rotation angle . The global coordi- natesystem( x;y)isusedtorepresentthelaminateproperties and associated loading conditions, whereas the manufactur- ing coordinate system ( x¨;y¨) is used to deﬁne the ﬁbre ori- entation and the overlap locations on a ply basis. Aone-dimensionalﬁbrevariationisdeﬁnedintheman- ufacturing coordinate system along the x¨axis by means of a Lagrange interpolation scheme, as suggested by Wu et al. [46] and shown in Equation 1, where the manufacturing angle is shown in Figure 4. This interpolation scheme al- lows for complex angle deﬁnitions using only a few control points(CP)perply. InEquation1,thereare Mcontrolpoints located at ﬁxed positions x¨ k. .x¨/=MÉ m=0 mMÇ kmx¨*x¨ k x¨ m*x¨ k(1) The manufactured ﬁbre orientation .x¨/must be trans- lated to the global coordinate angle , also deﬁned as posi- tive in the counter-clockwise direction and measured about the x-axis, according to Equation 2 and illustrated in Fig- ure 4. The exact method and considerations to rotate a ply are treated in Section 2.2. .x;y/= .x¨/+ (2) A grid in the global coordinate system is created by Ns xNssampling points to keep track of the thickness and stiﬀnessdistribution,necessaryforthebucklingsimulations. For each point in this grid, the following parameters are stored for every ply: •Manufacturing angle, •Global ﬁbre orientation, Table 1 Diﬀerent thickness proﬁle for the same global ﬁbre orientation. CP 1 CP 2 [°] 1[°]Thickness ratio [-] 2[°]Thickness ratio [-] Case 1 0 0 1 45ù 2 Case 2 45 -45ù 2 0 1 •Fibre steering curvature,  •Local thickness (either discrete or smeared), tlocal 2.2. Ply Rotation Onecanconstructthesameglobalﬁbreorientationwith diﬀerentplyrotationangles,asexempliﬁedbythetwocases ofTable1. Twocontrolpoints(CP)areused,creatingalin- ear variation of the orientation angle . Note that in both Case 1 and Case 2 the ﬁnal global ﬁber orientation is the same,asgivenby = + ofEquation2. However,thenon- linear thickness variation due to overlaps is only a function of the manufacturing angle , as discussed in Section 2.3 and 2.4, resulting in diﬀerent thickness proﬁles for Case 1 and Case 2. Until now, the literature only assumed uniform thickness at any location, such that the ply rotation angle could simply be added afterwards [19, 20, 26, 31, 43, 46]. By means of this example, it is clear that when overlaps are considered,therotationangle mustbejudgedconcurrently with the manufacturing angles . The methodology herein suggested, taking the rotation angleintoaccount,isbasedontheactualproductionprocess: a manufacturing area is assumed, according to the manu- facturingcoordinatesystem .x¨;y¨/,whichencompassesthe complete laminate, represented by the black-dotted area in Figure 5. This manufacturing area can be viewed as the re- gion where tows, such as the one in Figure 8, are being de- positednexttooneanotherbytherobotheaduntilfullycov- ering the laminate. The manufacturing angles are interpolated in the manufacturing coordinates system according to Equation 1. However, as the dimension along the x¨axis varies depend- ing on the rotation angle , the CP locations are expressed as percentages rather than ﬁxed positions. In this way, the CP can always be placed on the edges of the manufacturing area,whichareatleastcornersofthelaminate,suchthatan interpolation is possible over the complete domain. Thereafter, the global laminate mesh is remapped to the manufacturing coordinate system by means of Equation 3, which can then be used in the interpolation of Equation 1. Thenewx¨stillservesaslocationfortheinterpolation,how- ever,they¨partisalsoimportantinordertoknowtheexact overlap location, as there could only be partial overlapping dependingonthemanufacturingangle,asexplainedinSec- tion 2.3. L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 3 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 5: Manufacturing area encompassing a plate, with as- sociated rotated CP. < x¨ y¨= =4 cos. / *sin. / sin. / cos. /5< x y= (3) Finally,apositiveplyrotationyieldsavaryingmanufac- turing angle pattern from the top-left corner to the bottom- right corner, and the opposite for a negative rotation. The rotation of a ply can serve two purposes: either increasing the thickness for a given global ﬁbre orientation; to achieve aglobalorientationnotpossibleduetotheconstraintonthe manufacturing angles, as discussed later; or a combination thereof. The rotation strategy, and subsequent steps in this methodology, are demonstrated on a virtual manufacturing example for one ply, whose characteristics are given in Ta- ble 2. Note that the CP locations are given as percentages along the manufacturing coordinate system .x¨;y¨/and that for each CP there is a corresponding angle. Four CPs are used, leading to a third-order interpolation of the angle , according to Equation 1. In Figure 6, the ﬁbre orientation is plotted in the manufacturing coordinate system with the platerotated,wheretheCPvaluesof *30°and*50°areob- tained at the top-left and bottom-right corners. When the plate is viewed in its original global coordinate system, the obtained interpolation pattern is the one shown in Figure 7, where the ply rotation has also been applied to obtain the angle distribution . 2.3. Overlap Location To evaluate the overlap locations, every tow path is de- terminedandplotted[21]usingagivenopacitylevel. When theopaquetowssuperimpose,thedistinctopacitylevelscan berelatedtodiscretethicknessvalues. Thisprocedureisap- plied for each individual ply, whose separate thickness can thenbeaddedtogethertoobtainthetotallaminatethickness distribution.Table 2 Virtual manufacturing example characteristics. Variable Value Length 50 [mm] Width 50 [mm] -30 [ °] 1at 0% -30 [ °] 2at 30% 15 [ °] 3at 60% 65 [ °] 4at 100% -50 [ °] Figure 6: Virtual manufacturing example. Interpolated ﬁbre orientation in the manufacturing coordinate system. Figure 7: Virtual manufacturing example. Interpolated and rotated ﬁbre orientation in the global coordinate system. BymeansoftheLagrangeinterpolationschemeofEqua- tion1,thetowpathisrelatedtotheﬁbreorientationthrough Equation 4 [5] in the manufacturing coordinate system for a given ply. The tow path consists of a sequence of points x¨;y¨, in whichx¨is obtained using ﬁxed increments x¨, andy¨is obtained using forward Euler integration, where y¨is calculated using Equation 5. The step size x¨is de- ﬁned by the length of the encompassed manufacturing do- mainalongthe x¨axisdividedby Ns*1. Oncethetowpath L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 4 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 8: Virtual manufacturing example: ﬁbre orientation interpolation with corresponding tow path. for a single tow is obtained and plotted, the process is re- peated by changing the starting point according to the shift direction, until the complete plate is covered. The starting pointisshiftedbyadistanceequaltothetowwidth,assum- ing that for any ply a manufacturing angle of 0°will always bepresent,suchthatnogapsarepresent. Asingletowofthe virtual example is shown in Figure 8 in the manufacturing coordinatesystem,withthecorrespondingﬁbreorientation. dy¨ dx¨=tan (4) y¨=tan. .x¨//x¨(5) According to the aforementioned shifting strategy, the manufacturing angles = 0°result in no overlaps or gaps. However, as soon as 0the greater eﬀective width we, given by Equation 6 [5], will lead to overlapping ad- jacent tows. The limit case of a single full overlap is ob- tained for an eﬀective tow width of we=2w, covering half of the previous and subsequent tows. This corresponds to ð ð= 60°, according to Equation 6. Even more overlaps are appearing as ð ðincreases. The manufacturing angles and subsequent overlap formations are shown in Figure 9, with darker regions representing more thickness and over- laps. The exact amount of overlaps is determined by means ofpre-establishedthresholdsoftheopacity value. Theto- tal local thickness can be obtained by means of Equation 7, meaningthelowestnominaltowthicknessisobtainedwhere nooverlapsarepresent,andincreasesasanintegermultipli- cation for each given overlap. weùw cos (6) tlocal=.Noverlaps+1/ttow (7)As the overlap patterns are periodic along the shifting axis since no variation is present in that direction, only a restrained plotting surface using the full x¨length could be used in determining the exact overlap locations. Thereafter, the overlaps could be mapped to the complete rotated plate inthemanufacturingarea,basedontheirrelativecoordinate to this repeat surface. For each manufacturing grid point withtransformed x¨andy¨coordinates,theplythicknessac- cordingtoEquation7alongtheﬁbreorientationobtainedby the interpolation scheme of Equation 2 is stored. The dis- cretethicknesspatternofthevirtualmanufacturingexample is shown in Figure 10 and 11. 2.4. Smeared Thickness Castroetal.[10]proposedasmearedthicknessformula- tionforAFP-manufacturedlaminates,wherethelocalthick- nesstlocalisapproximatedbyEquation8. Thethicknesshas anon-lineardependencywiththemanufacturingsteeringan- gle deﬁnedinFigure12, tlocalandttowbeingrespectively the local smeared thickness and nominal tow thickness. tlocal.x¨;y¨/=ttow cos .x¨;y¨/(8) Thisapproximationhasseverallimitations,asﬁrstlyitis onlyvalidfor ð ð<90°,sincemanufacturinganglesoutside thisrangewillyieldanegativethicknesswhichisphysically not possible. Besides that, angles outside this range would meantowsarepointinginthe negativex¨directionandhence will not cover the remaining part of the ply. Furthermore, as ð ðapproaches90°, the smeared thick- ness of a nominal tow behaves asymptotically to inﬁnite as shown in Figure 12. In reality, the case =90°would pro- duceamajorthicknessbuild-upequaltothenumberofsub- sequent tows, which is neither desired nor will it result in inﬁnitethickness. Nevertheless,thesmearedapproximation forthecurrentstudyisrestrictedto ð ð<=60°alongsidethe considerations discussed in Section 3.1, which corresponds to an increase of 3% in local thickness for a 1°increment of the manufacturing angle , as illustrated in Figure 12. The smeared thickness distribution of the virtual manufacturing example shown in Figure 11 is given in Figure 13, where the smeared distribution results in the same volume as the discrete distribution. 2.5. Manufacturing Constraints The manufacturing process of variable stiﬀness lami- nates is signiﬁcantly constrained by the minimum steering radius, especially for AFP. In order to achieve feasible tow paths, their curvature must be lower than 1.575 1_m, as discussed in Section 1. To evaluate this constraint, the tow path radius can be precisely calculated by means of Equa- tion 9 [17], which is related to the curvature as=1_rs. rs=4 1+ dy¨ dx¨253 2 ðd2y¨ dx¨2ð(9) L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 5 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 9: Overlap formation depending on the manufacturing angle . Figure 10: Number of overlapping tows in the manufactur- ing coordinate system for the virtual manufacturing example showing the overlap locations. Applying this expression in the manufacturing coordi- nate system, the ﬁrst derivative of the ﬁbre path ( dy¨_dx¨) is the tangent of the ﬁbre orientation as shown in Equa- tion4[5],where .x¨;y¨/isdirectlycalculatedusingEqua- tion 1. The second derivatived2y¨ dx¨2becomes: d2y¨ dx¨2=1 cos2. /d dx¨(10) in which the derivative d _dx¨can be constructed in a closed form using Equation 1, leading to the expression given in Equation 11. d dx¨=MÉ m=0 mMÉ i=0;ijL 1 x¨ j*x¨ iMÇ m=0;m.i;j/x¨*x¨ m x¨ j*x¨ mM (11) The curvature constraint is evaluated for each ply at ev- ery grid point in the laminate mesh, remapped to the man- ufacturing coordinate system. The ﬁbre variation and steer- inghappensinthiscoordinatesystem,andissimplyrotated Figure 11: Discretethicknessdistribution, givenbythenumber of overlapping tows, in the global coordinate system for the virtual manufacturing example showing the overlap locations. Figure 12: Local tow thickness approximation and thickness sensitivitytoa 1°variationofthe angle, basedonEquation8. when viewed in the global coordinate system: the rotation doesnotalterthemanufacturingconstraints. Assoonasany mesh location does not satisfy the curvature limit, then that particularplycannotbemanufacturedandadiﬀerentvalue fortheCPlocationorCPanglehastobedeterminedinorder to render a feasible design. L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 6 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 13: Smeared thickness distribution, given in [ mm], in the global coordinate system for the virtual manufacturing ex- ample. 3. Linear Buckling Analysis Inthissection,thesimulationassumptionsarediscussed, followed by the derivation of the buckling problem for the semi-analytical problem, and ﬁnally the FEM implementa- tion. 3.1. Laminate Stiﬀness The VAT laminate mechanical properties will be sim- ulated on a mesoscale by means of the ABD stiﬀness ma- trix[37],whosegenericformulationisgiveninEquation12. Here,Qijis the ply stiﬀness [37], rotated by the local angle in the global coordinate system; zdeﬁnes the position of theplythrough-the-thickness. Toconsiderthevaryingthick- ness proﬁle, each ply thickness zk*zk*1is given by Equa- tion7andincaseofthesmearedapproach,theplythickness is given by Equation 8. Aij;Bij;Dij=³n k=1QijL zk*zk*1; z2 k*z2 k*1 2; z3 k*z3 k*1 3M (12) Duetotheoverlapconsiderations,thecrosssectiondoes not remain symmetric, as during manufacturing the plies are added on only one side, assuming that at the other side there is a smooth mould. Therefore, the geometrical mid planevariesthroughoutthelaminate,creatingavaryingoﬀ- set. This phenomenon was simulated by Groh et al. [18], whereavariablemidplanewasusedinthecurvedshelland 2Dplateformulations. Forthecurvedshellavaryingradius ofcurvaturewasusedthroughoutthelaminatetogetherwith thevaryingoﬀset,whereasforthe2Dplateonlythevarying oﬀset for the mid plane was used. Both approaches allowed a geometrically accurate representation of the ply thickness distributionovertheVATlaminate,butthe2Dplateformu- lation does not consider the coupling between out-of-plane displacements and in-plane strains due to the presence of a shell curvature created by the thickness buildup. Both2Dand3Doptionsprovedtobeingoodagreement for angles up to 60°, beyond which only the shell represen-tation renders a good accuracy [18]. The proposed method- ology of the present study uses the full ABD matrix in the buckling problem derivation, and limits the manufacturing angle to60°, alongside the thickness increase considera- tions of Section 2.4. Hence, a 2D plate model with correct mid-planeoﬀsetvaluesisselectedtorepresentthebehaviour of the VAT laminates with varying thickness. 3.2. Semi-analytical Approach The general buckling problem is derived based on the neutral equilibrium criterion of the total potential energy given in Equation 13. Following the derivation in [12, 9], the general form of Equation 14 can be obtained. 2=0 (13) .KKK+KgKgKg/=0 (14) whereKKKistheconstitutivestiﬀnessmatrixbasedonthe geometry and stiﬀness properties of the laminate; KgKgKgis the geometricstiﬀnessmatrixcalculatedusingthepre-buckling in-plane stresses. Their detailed derivation can be found in Castroetal. [9,11,12,13],andtheimplementationispub- licly available in a Python package [8]. Both KKKandKgKgKgare integrated numerically, given that the constitutive stiﬀness represented by the the ABD matrix spatially changes due to ﬁbre steering. The pre-buckling in-plane stresses are calculated by meansofastaticanalysis. Thisanalysisstemsfromtheﬁrst variation of the total potential energy functional: =0U*V=0 (15) Using the global coordinate system x;y, the in-plane displacements u.x;y/;v.x;y/andout-of-planedisplacement w.x;y/can be approximated as: u.x;y/=SuSuSu.x;y/ccc v.x;y/=SvSvSv.x;y/ccc w.x;y/=SwSwSw.x;y/ccc(16) where matrices SuSuSu;SvSvSv;SwSwSwcontain the shape functions oftheapproximateddisplacements,detailednext;and ccccon- tain the corresponding Ritz coeﬃcients. Thus, the variation of the strain energy Ucan be written as: U=cccTKKKccc (17) The work of the external forces V, considering only a force vector fffacting on the boundary ) of the two- dimensional plate domain , can be written as: V=cccTfff (18) L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 7 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Replacing these expressions in Equation 15 leads to: cccT.KKKccc*fff/=0 (19) which holds true for any ccc. Hence, the pre-buckling static equation is obtained: KKKccc=fff (20) Ifthein-planeforcesattheboundaryaregivenbyknown functionsfx.x;y/;fy.x;y/,fffcan be calculated as: fff=Ê) 4 SuSuSu SvSvSv5T<fx fy= d) (21) For displacement-controlled problems, the shape func- tions based on the Legendre polynomials allow a precise partition of the Ritz coeﬃcient’s vector cccsuch that the dis- placementsandrotationsofeachedgecanbeindependently deﬁned. For instance, setting flagt1andflagt2foruandv inEquation23. Inthiscase,thepartitionof fffcorresponding to the unknown Ritz coeﬃcients can be calculated from the pre-buckling static equation, and the remaining system can be solved for the unknown Ritz coeﬃcients. Displacement- controlledproblemsarenotincludedinthepresentworkand are left as a topic for future research. Theunknowndisplacementﬁeldissolvedbymeansof theRitzmethod,usingInthepresentworkthedisplacement ﬁeld is approximated using shape functions based on Leg- endrehierarchicpolynomials,asderivedbyRodrigues[34]. These polynomials have the advantage to allow a full con- troloftheboundaryconditionsfortheout-of-planedisplace- ment using the 4 ﬁrst terms given in Equation 22, whereas thehigher-ordertermsgiveninEquation24areusedforen- riching the inner domain. For the in-plane displacements, the boundary conditions are controlled using the 2 terms of Equation 23 instead of Equation 22, assuming that Su;v i=3= Su;v i=4= 0. For all cases, the inclusion of ﬂags allows (ﬂag = 1) or prohibits (ﬂag = 0) translations ( t1andt2) and rota- tions (r1andr2) at the plate edges, enabling the method to work with any combination of simply-supported, clamped and free boundary conditions [9]. Sw i=1.or/=.1 2*3 4+1 43/flagt1 Sw i=2.or/= 1 8*1 8*1 82+1 83 flagr1 Sw i=3.or/= 1 2+3 4*1 43 flagt2 Sw i=4.or/= *1 8*1 8+1 82+1 83 flagr2(22) Su;v i=1.or/= **1 2+1 flagt1 Su;v i=2.or/= .+1/ 2 flagt2(23)Table 3 Shape function terms for each ﬁeld variable. 1st term 2nd term 3rd term 4th term 5th term 6th term u Su;v i=1Su;v i=2Si=5Si=6Si=7Si=8 v Su;v i=1Su;v i=2Si=5Si=6Si=7Si=8 w Sw i=1Sw i=2Sw i=3Sw i=4Si=5Si=6 Si>4.or/=i_2É p=0.*1/p.2i*2p*7/@@ 2pp@.i*2p*1/@.i*2p*1/(24) Table 3 shows how the shape functions of Equation 22 - 24 are ordered for each ﬁeld variable. The order of the approximation polynomial along each direction ;is con- trolled separately and the ﬁnal approximation function is shown in Equation 25. Details regarding the implementa- tion are given in Castro and Donadon [9]. S.;/=³m i=1³n j=1Si./Sj./ (25) 3.3. FEM The four-node shell element S4 available in ABAQUS [39] was adopted for the ﬁnite element buckling analysis, using full numerical integration as recommended in [5, 46, 43]. The element properties are assigned based on the layup information of the manufacturing grid point closest to the element’s centroid, as performed in [41]. 3.4. Test cases Three diﬀerent cases are used to investigated the use of the smeared thickness approximation in buckling simula- tions. Allthecaseshavefourcontrolpoints(CP)foreachply, located respectively at 0%, 33%, 66% and 100%. A generic layup is applied to each case, as shown in Table 4, accord- ing to the notation [ ê 1ð 2ð 3ð 4ë] suggested by Guerdal et al. [19]. The orthotropic material properties are Exx= 181e3 GPa, Eyy= 10.3e3 GPa, xy= 0.28,Gxy= 7.2e3 GPa andttow= 0.127 mm [46], with a tow width wtowof 6.35 mm, or 1/4 inch. Case 1 Case 1, represented in Figure 14, is an equal bi-directional loading case for a square plate, simply supported on three sidesandpinnedonthelast,withfourpliesanddimensions of 1500x1500 mm. Case 2 Case2,showninFigure15,representsaﬂange,whereaplate of higher aspect ratio is simulated, with three edges simply supported and one long edge free. It is only loaded in the xdirection, and has a length of 1000 mm and width of 200 mm. The laminate consists of four plies. L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 8 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 14: Deﬁnition of test case 1. Figure 15: Deﬁnition of test case 2. Figure 16: Deﬁnition of test case 3. Table 4 Test cases layup. Case 1[[*4ê45ð7ð13ð6ë];[*59ê59ð5ð37ð*25ð*19ë]; [16ê27ð*25ð*8ð40ë];[1ê*39ð53ð15ð*14ë]] Case 2[[42ê*25ð48ð54ð*8ë];[*86ê*42ð9ð*40ð24ë]; [14ê*42ð*25ð*7ð30ë];[4ê16ð*52ð*16ð*35ë]] Case 3[[*21ê*8ð*48ð12ð40ë];[*60ê*39ð*24ð*4ð42ë]; [*42ê*43ð54ð*29ð*3ë];[*68ê2ð*48ð*45ð*33ë]; [45ê23ð*5ð43ð*13ë];[47ê53ð*25ð13ð10ë]; [52ê3ð*7ð49ð*43ë];[98ê*33ð*29ð*10ð*39ë]] Case 3 Case3,representedinFigure16,isamultiloadingcasebe- ing applied to 2000 x 1000 mm geometry. It is loaded by a runningloadinthe ydirection30%ofthedistributedloadin thexdirection and a shear load applied at the opposite side oftheclampwithatotalloadhalfthatofthexdirection. The laminate consists of eight plies. Figure 17: Tow path radius comparison between the virtual manufacturing and the analytical expression of Equation 26. 4. Virtual Manufacturing Veriﬁcation The veriﬁcation of the proposed virtual manufacturing approach for VAT laminates is done in two stages: 1)start- ing with the veriﬁcation of the radius of curvature of the towpaths; 2)comparisonbetweenthethicknesspatternsob- tained with the discrete and smeared approaches. The radius calculated with Equation 9 is veriﬁed for the virtualmanufacturingexampleofFigure8andTable2. The method herein proposed is compared to the analytical ex- pressiongiveninEquation26[20],withtheresultshownin Figure 17. The graph has been cut oﬀ at rs= 100mmfor better visualization and both methods give the same result, giving conﬁdence to use the approach herein proposed for calculating manufacturability constraints of each VAT ply. r=a sin 1*sin 0(26) The thickness distribution of the obtained patterns for both discrete and smeared simulation are compared for all three cases of Section 3.4 and given in Figures 18 - 20. There is a good visual agreement in the thickness patterns forallthreecaseshereininvestigated. Comparingthevalues of the maximum thickness tmaxin Table 5 it can be noted that, for the discrete case, higher thicknesses are achieved duetothelocallyoverlappingtows;whereasinthesmeared casetheseoverlapsareaveragedout,resultinginasmoother thickness distribution and therefore a lower value for tmax. The quantitative accuracy of the obtained thickness pattern is performed by means the same volume requirement [10], heretranslatedtoanaveragethicknessrequirement. Table5 presents the relative average thickness error calculated us- ing the discrete and smeared conﬁgurations, showing that a maximum error of 0:14~occurred for Case 1. 5. Linear Buckling Results Before comparing the results of the diﬀerent simula- tions, the convergence of the models is checked: this is L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 9 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 18: Discrete (top) and smeared (bottom) thickness proﬁle of case 1. Figure 19: Discrete (top) and smeared (bottom) thickness proﬁle of case 2. Table 5 Particular value comparison between the discrete and smeared thickness. tmax[mm] taverage[mm] Relative error of the average Discrete Smeared Discrete Smeared thickness [%] Case 1 1.02 0.66 0.57 0.57 0.14 Case 2 1.02 0.71 0.62 0.62 0.07 Case 3 1.91 1.38 1.22 1.22 0.10 assessed in the FEM model by changing the number of elements, whereas in the semi-analytical model, the con- vergence is controlled by the amount of integration points and the order of the Legendre polynomials constituting the shapefunctions. Bydefault,200x200virtualmanufacturing grid points Nsare used to generate the thickness and Figure 20: Discrete (top) and smeared (bottom) thickness proﬁle of case 3. layup information for the FEM model, augmented with the integration point locations in case of the semi-analytical model,alongthematerialpropertiespresentedinSection3.4 5.1. Smeared Thickness Approximation TheconvergenceoftheFEMresultswithsmearedthick- ness is given in Appendix A, and these results are used predominately for verifying the semi-analytical implemen- tation. Semi-Analytical Modelling Figure 21 shows the convergence of the ﬁrst eigenvalue for allthreecasesusingthesameamountofintegrationpointsin eachdirection(N)equaltotheorderoftheLegendrepolyno- mials (cf. Table 3); with one extra integration point (N+1); and with three extra integration points (N+3). The idea of evaluating over-integrated schemes came from the possibil- itythatthevariablestiﬀnessrequiresmorespatialdiscretiza- tion than the displacement ﬁelds being approximated. The results,however,showthatusing12shapefunctionsandan amountofintegrationpointsequaltoNsuﬃcesinobtaining a converged result. 5.2. Discrete Thickness Semi-Analytical Modelling For the semi-analytical simulation of the discrete thickness case, Figure 22 shows again the convergence results for N, N+1 and N+3 integration points; with N being the order of the Legendre polynomials used as shape functions. Note in Figure 22 that there is more oscillations in the convergence of the discrete thickness models, when com- pared to the convergence of the smoothed thickness model given in Figure 21. Increasing the amount of integration points for low-order shape functions, makes the result os- cillate up to 20%. This is expected since the discrete thick- L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 10 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 21: Semi-analytical model convergence of the smeared thickness approach for N, N+1 and N+3 integration points. Figure 22: Semi-analytical model convergence of the discrete thickness approach for N, N+1 and N+3 integration points; whereNistheorderoftheLegendrepolynomialsusedasshape functions. nessgeneratesadiscontinuousstiﬀness,whichispoorlyap- proximated by quadrature rules original developed for con- tinuousfunctions[32]. Forhigher-ordershapefunctionsthe solution converges, whilst showing smaller oscillations for the same reason. Future studies could attempt diﬀerent nu- merical integration schemes, such as presented in Abedian and Düster [1], to evaluate and mitigate the observed con- vergence issues. Nevertheless, for the present work, the re- quired number of shape functions and integration points for aconvergedsemi-analyticalsimulationcannotbepredicted withcertainty,andthereforethediscretethicknessapproach is not advised in the context of semi-analytical modelling. Finite Element Method The FE buckling eigenvalue convergence is shown in Fig- ure 31 for all three cases. It does not suﬀer from oscilla- tions,incontrasttothesemi-analyticalformulation. Usinga 1%convergencecriterion,convergenceisattainedacrossall three cases with 50 elements in each direction. Figure 23: FE convergence using the discrete thickness ap- proach. Table 6 FE reﬁned and converged eigenvalue comparison for the dis- crete thickness. Reﬁned model [-] Converged model [-] Diﬀerence [%] Case 1 0.74 0.74 0.11 Case 2 8.49 8.52 0.33 Case 3 24.92 24.66 -1.04 The converged simulations are ﬁnally compared to a much more reﬁned virtual manufacturing model and simulation, to assess the inﬂuence of the thickness grid deﬁning the element’s properties. These reﬁned models are run with a 1000x1000 virtual manufacturing grid and 100elementsineachdirection,withtheresultspresentedin Table 6. This comparison shows little inﬂuence of the grid reﬁnement, and that the selected 200x200 manufacturing grid is suﬃcient to provide the required information on thicknessandlayupforaconvergedFEMsimulationforthe discrete thickness case. Discussion The FEM results of the discrete thickness converge with 50 elements in each direction, using a convergence criterion of 1%. On the other hand, it is diﬃcult to express parameters for a converged semi-analytical model using discrete thick- ness,asitdoesnotapproximatethediscontinuousthickness well, leading to an oscillation of the results even with 16th-order Legendre polynomials. Despite this oscillation, the global trend from the semi-analytical model approach the FEM results. Finally, it was veriﬁed that a 200x200 virtualmanufacturinggridisenoughtoconveytherequired thickness and layup information for the discrete thickness representation to the FEM simulation. 5.3. Validity of the Thickness Approximation Linear buckling analyses performed using the smeared thickness and semi-analytical modeling are further com- pared to the discrete thickness FEM simulation. Table 7 L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 11 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Table 7 Comparison of the semi-analytical smeared approximation model and the discrete FEM simulation eigenvalue. Discrete thickness FEM[-]Smeared thickness semi-analytical [-]Diﬀerence [%] Case 1 0.74 0.75 1.20 Case 2 8.52 8.12 -4.71 Case 3 24.66 25.79 4.59 Table 8 Eigenvalue comparison for the three veriﬁcation cases with ad- ditional layups. Case 1 Case 2 Case 3 Layup number 2 3 2 3 2 3 Smeared thickness semi-analytical [-] 0.71 0.77 5.68 13.23 25.89 25.07 Discrete thickness FEM [-] 0.70 0.77 5.84 13.54 25.92 25.03 Relative diﬀerence [%] 0.67 -0.25 -2.82 -2.27 -0.14 0.18 givesanoverviewoftheﬁrsteigenvalue obtainedforeach of thethree casesusing bothways ofsimulation. Tofurther assesswhetherdiﬀerencesarecaseorlayupspeciﬁc,twoad- ditionalrandomlayups(giveninAppendixB)aresimulated both ways for every case. These results can be found in Ta- ble 8. Firstly, considering all the data, the smeared modelling has a,5~relative error range with respect to its discrete thickness counterpart. This is a satisfactory agreement be- tweenthethicknessapproximationprovidedbythesmeared approach and the real thickness distribution of the discrete model. Secondly, the range of relative error is regarded as being intrinsic to a given layup sequence and not linked to the smeared approximation. It is assumed to stem from the inﬂuence of the local varying discrete thicknesses as op- posedtothesmoothsmearedapproximation,aspresentedin Table 5. The patterns of the ABD components distribution coincide,aswasthecasewiththetotalthicknessbuild-upof Section 4. The locally diﬀerent thickness can be of signif- icant impact, as the B and D matrices depend respectively onthesquareandcubicofthethickness,alongsidetheﬁbre orientation. Coupled to the complex and interlinked rela- tions forming the linear constitutive and geometric stiﬀness matrix, this altogether is believed to lead to the errors seen in Table 7 and 8. Finally, the eigenmode of the ﬁrst eigen- value of each test case is compared for both modelling of the smeared semi-analytical model, and the discrete thick- ness FEM analysis in Figures 24- 26, showing similar pat- terns. Therefore, it is concluded that the smeared thickness approximation suggested by [10] is well suited to represent the behaviour of its discrete thickness VAT laminate coun- terpart. 6. Tow Width Ratio Limit Theaforementionedresultsshowthatthesmearedthick- ness approach gives a good approximation for the linear buckling response of VAT laminates with discrete overlaps. However,thissimpliﬁedthicknessproﬁlehaslimitedrepre- Figure 24: Discrete (left) and smeared (right) thickness, 1st eigenmode of case 1. Figure 25: Discrete (top) and smeared (bottom) thickness, 1st eigenmode of case 2. Figure 26: Discrete (left) and smeared (right) thickness, 1st eigenmode of case 3. sentation capabilities. The approximation is well suited for globalresponses,whilelocalphenomenasuchasstresscon- centration,localstrengthandfailurecannotbepredicted,es- peciallyatandaroundoverlapregions. Theseaspectsshould be evaluated locally with a reﬁned model of the exact phys- ical representation, as the smeared thickness approximation does not take into account any information about tow inter- actions. Furthermore, the inﬂuence of the tow width is investi- gated on the validity of the smeared approximation, as it is expectedtoprogressivelybecomeaworseapproximationfor largertowwidths,sincetheoverlapsaremoredistinctasthey appear less frequently, and have a larger area. To this end, allthreecasesofSection3.4areevaluatedagainwithadis- crete thickness FEM analysis for diﬀerent tows widths. The resultsareshowninTable9,andtherelativediﬀerencewith their respective smeared eigenvalue in Table 10. In order to comparetheeﬀectonthediﬀerentgeometries,thetowwidth is normalised with the smallest dimension of the structure: thisdimensionwillremainthesmallestdimensioninanyro- tatedmanufacturingframe,whilstbeingthemostinﬂuenced by the larger tow widths, creating a more discrete and less L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 12 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Table 9 Buckling eigenvalue for diﬀerent tow widths. Tow width [mm] 3.175 6.35 25.4 40 65 100 175 350 Case 1[-] 0.73 0.74 0.74 0.74 0.75 0.76 0.76 0.77 Case 2[-] 8.57 8.52 8.40 8.50 8.58 8.72 8.79 8.14 Case 3[-] 24.65 24.66 25.03 24.48 25.63 25.61 25.63 25.23 Table 10 Buckling eigenvalue relative error for diﬀerent tow widths com- pared to the smeared approximation. Tow width [mm] 3.175 6.35 25.4 40 65 100 175 350 Case 1[-] -3.04 -1.19 -1.71 -1.81 -0.35 0.62 1.07 1.80 Case 2[-] 5.63 4.95 3.48 4.75 5.74 7.46 8.27 0.35 Case 3[-] -4.45 -4.39 -2.97 -5.08 -0.62 -0.72 -0.62 -2.20 Figure 27: First linear buckling eigenvalue error, of the smeared relative to the discrete thickness approach, for dif- ferent tow widths per minimum dimension ratios. smooth thickness distribution. The diﬀerence for the nor- malized dimensions is displayed in Figure 27. ThedatainFigure27doesnotshowtheexpecteddeteri- orationbetweenthediscretethicknessFEMsimulationsand smearedsemi-analyticalsolution. Despitethelackoftheap- parentinﬂuenceofthetowwidth,theauthorsdonotbelieve the smeared approach to hold valid for the complete simu- lated parameter range. As stated previously, the ﬁbre orien- tation is taken at exactly the grid points from the interpola- tion: with narrow tows, the grid points will always be close to a tow centerline. However, when large tows are present many grid points can be within one tow width, but the ori- entation at those grid points is then dictated by the center- line orientation and bending eﬀect through the width of the towduetosteering,ratherthantheinterpolationscheme,and would therefore induce errors. Therefore, the eﬀect of the tow width on the validity of the thickness approximation cannot be fully described by the current framework, since the exact ﬁbre orientation is badly predicted for large tow width ratios. However, by vi-sually comparing the discrete thickness given in Figure 28- 30 with their respective smeared proﬁle in Figure 18- 20, the authors propose that: the maximum tow width per min- imum dimension ratio should be 0.15. With this limit, it is guaranteed an accurate representation of the laminate stiﬀ- nessandofthethicknessdistribution,whenusingthesimple smearedapproach. Beyondthisthreshold,onestartsloosing global similarities and local particularities can no longer be smoothed out. Further investigations should be carried out with a detailed tow model taking the ﬁbre bending into ac- count, whereby sampling points can be assigned along the tow’s transverse direction, at a large distance from the tow centreline. Theguidelineshereinproposed,alongwiththegoodap- proximationofthesmearedapproach,canbehelpfulforde- signers to quickly assess VAT laminates with overlaps in a practical way for buckling as demonstrated in this paper, or just to obtain the pre-buckling or any general linear static state, as explained in Section 3.2. Moreover, the combina- tion with the semi-analytical modelling will serve as a ba- sis for future work on the design and optimisation of VAT structureswiththicknessvariationcoupledwiththesteering angles to increase buckling performance. 7. Conclusion In this paper, Variable Angle Tow (VAT) laminates are ﬁrstreplicatedvirtuallytoobtaintheoverlaplocationsalong the ﬁbre angle distribution. The manufacturing ﬁbre varia- tion is parameterized by means of a Lagrange interpolation scheme in the manufacturing coordinate system, which is considered alongside the complete rotation with respect to thelaminatecoordinatesystem. Subsequently,thetowpaths are plotted graphically with a given opacity, allowing to be superimposed and retrieve a discrete thickness proﬁle and overlap locations. This discrete thickness proﬁle is approx- imated with a continuous smeared approach, only dictated by the manufacturing angle, to simulate linear buckling. A procedure is also integrated to verify the manufacturability of each ply of the VAT laminates. Furthermore,theonesidedthicknessproﬁleoftheVAT laminates is incorporated in the simulation with an oﬀset in the ABD formulation, albeit limited to a manufacturing angle of60°. The semi-analytical formulation is based on the Kirchhoﬀ-Love plate kinematics and taking all ABD termsintoaccount,wherethedisplacementﬁeldsaresolved by means of the Ritz method, with Legendre polynomial shapefunctionsincorporatingdiﬀerentboundaryconditions through enrichment functions. However, before the thick- nessandlayupinformationcouldbeusedinthesimulations, both discrete and smeared thickness proﬁles are compared witheachothertoverifythevirtualmanufacturingmethod- ology, showing good agreement, with similar locations and shape of thickness build-up. With the pre-processing information veriﬁed, the buck- ling analysis and thickness approximation validity was in- vestigated. Inaﬁrstinstance,thesemi-analyticalmodelwas L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 13 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 28: Virtual manufacturing of case 1 with tow widths: 3, 25, 40, 65, 100, 175, 350 [ mm]; corresponding to a ratio of tow width to minimum plate dimension of: 0.002, 0.017, 0.027, 0.043, 0.067, 0.117, 0.233. After a ratio of 0.15, the local thickness features become more apparent, and the smeared approach is no longer recommended. veriﬁed against a Finite Element Method (FEM) model for diﬀerentgeometries,loadingandboundaryconditions. This was done based on the smeared approximation with a ran- dom variable stiﬀness layup, and showed good agreement between both simulation methods. Thereafter, the discrete thickness proﬁle was also modelled with both the FEM and semi-analytical modelling, with the ﬁrst method experienc- ing convergence issues due to the discontinuous thickness. The FEM method converged, whose results were in satis- factory agreement with the semi-analytical outcome using thesmearedapproximation,forboththebucklingeigenvalue and corresponding eigenmode. However, using the smeared approach for buckling sim- ulationsisnotvalidforanyarbitraryvalueoftowwidth,and thereforeaguidelineisproposedsuchthatthemaximumtow width should be smaller than 15% of the smallest plate di- mensionwhereathresholdof0.15onthetowwidthoverthe smallestplatedimensionissuggested, based on the limits of the visual correlation between the smeared and discrete thickness proﬁles. With this guideline, the smeared thick- ness approach can be used in the initial assessment of VAT laminateswithoverlapconﬁgurations,andonthedesignand optimisationofVATstructureswiththicknessvariationcou- pled with the steering angles. Appendix A. Smeared thickness FEM convergence The smeared thickness approach is also simulated by means of FEM, to verify the correct implementation of the semi-analytical model. Figure 31 shows the ﬁrst buckling eigenvalue for an increased number of elements, seeded equally in both x and y direction, with the thickness and layupinformationconstructedona200x200pointmanufac- turing grid. From this data, a 1% convergence is achieved with 15 elements for all cases, regardless of the dimensions of theTable 11 FEM reﬁned and converged eigenvalue comparison for the smeared thickness. Reﬁned model [-] Converged model [-] Diﬀerence [%] Case 1 0.75 0.75 0.6 Case 2 7.95 8.02 0.86 Case 3 25.35 25.58 0.92 Table 12 Converged smeared semi-analytical and FEM eigenvalue com- parison. Discrete thickness FEM[-]Smeared thickness semi-analytical [-]Diﬀerence [%] Case 1 0.75 0.75 0.19 Case 2 8.02 8.12 1.23 Case 3 25.58 25.79 0.84 plate. These 15 elements simulations are ﬁnally compared to a much more reﬁned virtual manufacturing model, to as- sesstheinﬂuenceofthegridtoconstructtheelements. This reﬁned model is created with a 1000x1000 point grid for the thickness and layup information, and 100 elements in each direction. The results of this comparison, given in Ta- ble 11, show that the reﬁnement of the FEM and manufac- turing model has little inﬂuence on the buckling outcome, meaning a 200x200 grid is enough for the virtual manufac- turing. The comparison between the converged ﬁrst eigenvalue of both the semi-analytical and FEM simulation for the smeared approximation are given in Table 12. This shows that the semi-analytical predictions are in good agreement with the FEM results, verifying the semi-analytical imple- mentation. L. Vertonghen and S. G. P. Castro: Preprint submitted to Elsevier Page 14 of 17 Modelling of ﬁbre steered plates with coupled thickness variation from overlapping continuous tows Figure 29: Virtualmanufacturingofcase2withtowwidths: 3, 25, 40, 65, 100, 175, 350 [ mm]; corresponding to a ratio of tow width to minimum plate dimension of: 0.015, 0.125, 0.200, 0.325, 0.500, 0.875, 1.750. After a ratio of 0.15, the local thickness features become more apparent, and the smeared approach is no longer recommended. B. Additional Layup Comparison Table 13 presents the two additional layups for each of the 3 cases to compare the validity of the smeared semi- analytical buckling simulation with the discrete FEM coun- terpart. CRediT authorship contribution statement L. Vertonghen: Conceptualization, methodology, soft- ware, formal analysis, investigation, validation, writing -originaldraftpreparation,writing-review editing. S.G.P. Castro:Conceptualization, supervision, software, data cu- ration, writing - original draft preparation, writing - review editing.