PDF
Optimization of a Variable-Stiffness Composite Plate with a Central Hole
Introduction
In traditional composite laminates, fibers within each ply are straight, and the global stiffness response is tailored through the selection of an appropriate stacking sequence. In contrast, variable-angle tow (VAT) laminates employ curvilinear fibers within each layer. By allowing the fiber orientation to vary spatially, VAT laminates significantly expand the design space and deliver improved structural performance. This enhanced capability can potentially reduce the required number of plies and the overall weight of the structure.
This example shows how to set up an optimization study to maximize the stiffness of a composite plate using the fiber orientation as the control variable.
Read more about the Composite Materials Module in the COMSOL blog, Introduction to the Composite Materials Module.
Model Definition
A one-layer square plate with a hole at its center is considered. The plate is subject to uniaxial traction. Only a quarter of the plate is modeled due to the symmetry of the boundary conditions and loads.
Figure 1: Geometry and applied loads.
To maximize the stiffness, the total elastic strain energy is chosen as the objective function to minimize, and the fiber angle with respect to the load direction is chosen as the control variable.
Material Properties
The composite lamina is made of carbon reinforced thermoplastic. The homogenized transversely isotropic material properties (Young’s modulus, shear modulus, and Poisson’s ratio) are given in Table 1.
Table 1: Material properties of a lamina.
Material property
Value
{E11, E22}
{138, 8.7} GPa
G12
5 GPa
{υ12, υ23}
{0.28, 0.45}
Results and Discussion
The optimized solution (Figure 2) is compared to a plate with straight fibers. Figure 3 and Figure 4 show the value of the total strain energy and the maximum displacement for different ply angles using a constant fiber angle. The VAT provides a stiffer solution (lower maximum displacement) with lower level of elastic strain energy. The plate with curvilinear fibers provide a lower maximum stress when comparing the stress fields for the optimal solutions for straight and curvilinear fibers (Figure 5, Figure 6).
Figure 2: Fiber orientation that maximize the plate stiffness for the given boundary conditions.
Figure 3: Strain energy values for composites with straight fibers. The horizontal dashed line indicates the optimum value using curvilinear fibers.
Figure 4: Maximum displacement for composites with straight fibers. The horizontal dashed line indicates the optimum value using curvilinear fibers.
Figure 5: Stress field for the composite plate with straight fibers oriented at 0 deg with respect to the x-axis. The deformation is amplified.
Figure 6: Stress field for the optimum solution with curvilinear fibers. The deformation is amplified.
Notes About the COMSOL Implementation
•
Modeling a composite laminate as a layered shell requires a surface geometry, referred to as a the base surface, and a Layered Material node, which adds an extra dimension (1D) in the normal direction. You can use the Layered Material functionality to model several layers stacked on top of each other, each having a different thickness, material properties, and fiber orientation. You can optionally specify the interface material between the layers and control the number of through-thickness mesh elements for each layer. To model a curvilinear fiber, add an expression as a function of the local coordinates in the Rotation field.
•
The third direction for the selected coordinate system in the Single Layer Material, Layered Material Link, or Layered Material Stack represents the normal direction. This is also the direction in which the layer stacking is interpreted from bottom to top, and therefore, it is crucial to visualize it during modeling. There are two ways to achieve this:
-
Using physics symbols: Go to the physics settings, find the Physics Symbols section, and select the Enable physics symbols checkbox. Then go to the material feature, for instance, Linear Elastic Material, to see the normal direction represented by green arrows.
-
Using result templates: When a solution dataset is available, use the result template Thickness and Orientation to plot the normal direction.
•
From a constitutive model point of view, you can either use the Layerwise (LW) theory available in the Layered Shell interface, or the Equivalent Single Layer (ESL) theory available in the Linear Elastic Material, Layered node in the Shell interface. The laminated composite presented in this example uses the Layered Shell interface.
•
The built-in Composites material library contains data for fiber and matrix constituents as well as for unidirectional and bidirectional laminae.
Application Library path: Composite_Materials_Module/Tutorials/fiber_angle_optimization
Modeling Instructions
From the File menu, choose New.
New
In the New window, click  Model Wizard.
Model Wizard
1
In the Model Wizard window, click  3D.
2
In the Select Physics tree, select Structural Mechanics > Layered Shell (lshell).
3
Click Add.
4
Click  Study.
5
In the Select Study tree, select General Studies > Stationary.
6
Click  Done.
Global Definitions
Parameters 1
Load the parameters from a file.
1
In the Model Builder window, under Global Definitions click Parameters 1.
2
In the Settings window for Parameters, locate the Parameters section.
3
Click  Load from File.
4
Browse to the model’s Application Libraries folder and double-click the file fiber_angle_optimization_parameters.txt.
Variables (Angle Ply for Straight Fibers)
1
In the Model Builder window, right-click Global Definitions and choose Variables.
2
In the Settings window for Variables, locate the Variables section.
3
In the table, enter the following settings:
Name
Expression
Unit
Description
plyAngle
theta
rad
 
4
In the Label text field, type Variables (Angle Ply for Straight Fibers).
5
Locate the Variables section. In the table, enter the following settings:
Name
Expression
Unit
Description
plyAngle
theta
rad
Fiber orientation with respect X-axis
COMSOL Multiphysics is equipped with built-in material properties for a number of composite materials. Select the needed material from the Composites material folder in the built-in material library.
Add Material
1
In the Materials toolbar, click  Add Material to open the Add Material window.
2
Go to the Add Material window.
3
In the tree, select Composites > Laminae > Unidirectional fiber lamina: AS4/APC2 carbon/PEEK thermoplastic [fiber volume fraction 58%].
4
Right-click and choose Add to Global Materials.
5
In the Materials toolbar, click  Add Material to close the Add Material window.
Define a layered material with ply rotations as parameters to be optimized.
Global Definitions
Layered Material 1 (lmat1)
1
In the Model Builder window, under Global Definitions right-click Materials and choose Layered Material.
2
In the Settings window for Layered Material, locate the Layer Definition section.
3
In the table, enter the following settings:
Layer
Material
Rotation
Value
Thickness
Mesh elements
Layer 1
Unidirectional fiber lamina: AS4/APC2 carbon/PEEK thermoplastic [fiber volume fraction 58%] (mat1)
plyAngle
-
th
1
Geometry 1
Work Plane 1 (wp1)
In the Geometry toolbar, click  Work Plane.
Work Plane 1 (wp1) > Plane Geometry
In the Model Builder window, click Plane Geometry.
Work Plane 1 (wp1) > Square 1 (sq1)
1
In the Work Plane toolbar, click  Square.
2
In the Settings window for Square, locate the Size section.
3
In the Side length text field, type a.
Work Plane 1 (wp1) > Circle 1 (c1)
1
In the Work Plane toolbar, click  Circle.
2
In the Settings window for Circle, locate the Size and Shape section.
3
In the Radius text field, type r0.
Work Plane 1 (wp1) > Difference 1 (dif1)
1
In the Work Plane toolbar, click  Booleans and Partitions and choose Difference.
2
Select the object sq1 only.
3
In the Settings window for Difference, locate the Difference section.
4
Click to select the  Activate Selection toggle button for Objects to subtract.
5
Select the object c1 only.
6
Click the  Show Grid button in the Graphics toolbar.
7
In the Home toolbar, click  Build All.
Materials
Layered Material Link 1 (llmat1)
1
In the Model Builder window, under Component 1 (comp1) right-click Materials and choose Layers > Layered Material Link.
2
In the Settings window for Layered Material Link, locate the Orientation and Position section.
3
Click  Go to Source for Coordinate system.
Definitions (comp1)
Boundary System 1 (sys1)
1
In the Model Builder window, under Component 1 (comp1) > Definitions click Boundary System 1 (sys1).
2
In the Settings window for Boundary System, locate the Settings section.
3
Find the Coordinate names subsection. From the Axis list, choose x.
Layered Shell (lshell)
1
In the Model Builder window, under Component 1 (comp1) click Layered Shell (lshell).
2
In the Settings window for Layered Shell, click to expand the Discretization section.
3
From the Displacement field list, choose Quadratic serendipity.
Symmetry 1
1
In the Physics toolbar, click  Edges and choose Symmetry.
2
Select Edges 1 and 4 only.
Prescribed Displacement 1
1
In the Physics toolbar, click  Edges and choose Prescribed Displacement.
2
Select Edge 5 only.
3
In the Settings window for Prescribed Displacement, locate the Prescribed Displacement section.
4
From the Displacement in z direction list, choose Prescribed.
Boundary Load 1
1
In the Physics toolbar, click  Edges and choose Boundary Load.
2
Select Edge 5 only.
3
In the Settings window for Boundary Load, locate the Force section.
4
Specify the fA vector as
sigma0
x
Mesh 1
Free Triangular 1
1
In the Mesh toolbar, click  More Generators and choose Free Triangular.
2
Select Boundary 1 only.
3
In the Settings window for Free Triangular, click  Build All.
Size
1
In the Model Builder window, click Size.
2
In the Settings window for Size, locate the Element Size section.
3
From the Predefined list, choose Extra fine.
4
Click  Build All.
Definitions (comp1)
Maximum 1 (maxop1)
1
In the Definitions toolbar, click  Nonlocal Couplings and choose Maximum.
2
In the Settings window for Maximum, locate the Source Selection section.
3
From the Geometric entity level list, choose Boundary.
4
Select Boundary 1 only.
Study 1: Original Orientation
1
In the Model Builder window, click Study 1.
2
In the Settings window for Study, type Study 1: Original Orientation in the Label text field.
Parametric Sweep
1
In the Study toolbar, click  Parametric Sweep.
2
In the Settings window for Parametric Sweep, locate the Study Settings section.
3
Click  Add.
4
In the table, enter the following settings:
Parameter name
Parameter value list
Parameter unit
theta (Fiber orientation)
range(-90,5,90)
deg
5
In the Study toolbar, click  Compute.
Set default units for result presentation.
Results
Preferred Units 1
1
In the Results toolbar, click  Configurations and choose Preferred Units.
2
In the Settings window for Preferred Units, locate the Units section.
3
Click  Add Physical Quantity.
4
In the Physical Quantity dialog, select Solid Mechanics > Stress tensor (N/m^2) in the tree.
5
Click OK.
6
In the Settings window for Preferred Units, locate the Units section.
7
In the table, enter the following settings:
Quantity
Unit
Preferred unit
Stress tensor
N/m^2
MPa
8
Click  Add Physical Quantity.
9
In the Physical Quantity dialog, select General > Displacement (m) in the tree.
10
Click OK.
11
In the Settings window for Preferred Units, locate the Units section.
12
In the table, enter the following settings:
Quantity
Unit
Preferred unit
Displacement
m
mm
13
Select the Apply conversions to expressions with the same dimensions checkbox.
14
Click  Apply.
Layered Material (Original Orientation)
1
In the Model Builder window, expand the Results > Datasets node, then click Layered Material.
2
In the Settings window for Layered Material, type Layered Material (Original Orientation) in the Label text field.
Use the Mirror 3D datasets for full model visualization.
Mirror 3D 1
1
In the Results toolbar, click  More Datasets and choose Mirror 3D.
2
In the Settings window for Mirror 3D, locate the Data section.
3
From the Dataset list, choose Layered Material (Original Orientation).
Mirror 3D (Original Orientation)
1
In the Results toolbar, click  More Datasets and choose Mirror 3D.
2
In the Settings window for Mirror 3D, type Mirror 3D (Original Orientation) in the Label text field.
3
Locate the Data section. From the Dataset list, choose Mirror 3D 1.
4
Locate the Plane Data section. From the Plane list, choose xz-planes.
Stress (Original Orientation)
1
In the Model Builder window, under Results click Stress (lshell).
2
In the Settings window for 3D Plot Group, type Stress (Original Orientation) in the Label text field.
3
Locate the Data section. From the Dataset list, choose Mirror 3D (Original Orientation).
4
From the Parameter value (theta (deg)) list, choose 0.
5
Locate the Color Legend section. Select the Show maximum and minimum values checkbox.
6
Select the Show units checkbox.
7
In the Stress (Original Orientation) toolbar, click  Plot.
Now, setup the model to perform an optimization analysis using the fiber angle as control variable and the elastic strain energy as objective function.
Use a Control Function to optimize the fiber angle. Set the angle to be bounded between -90 and +90 degree and insert an initial guess.
Definitions (comp1)
Control Function 1 (cfunc1)
1
In the Definitions toolbar, click  Control Variables and choose Control Function.
2
In the Settings window for Control Function, locate the Output section.
3
In the fmin text field, type -pi/2.
4
In the fmax text field, type pi/2.
5
In the c0 text field, type sqrt(pi/4)*cos(cfunc1.x*2*pi).
6
Locate the Control Variable Discretization section. From the Control type list, choose Piecewise Bernstein polynomial.
7
Locate the Units section. Click  Select Input Quantity.
8
In the Physical Quantity dialog, type dim in the text field.
9
In the tree, select General > Dimensionless (1).
10
Click OK.
Control Function 2 (cfunc2)
1
Right-click Control Function 1 (cfunc1) and choose Duplicate.
2
In the Settings window for Control Function, locate the Output section.
3
In the c0 text field, type sqrt(pi/4)*sin(cfunc2.x*2*pi).
Global Definitions
Use the control functions to describe the distribution of the angle to be optimized by the solver.
Variables (Angle Ply for Curvilinear Fibers)
1
In the Model Builder window, right-click Global Definitions and choose Variables.
2
In the Settings window for Variables, locate the Variables section.
3
In the table, enter the following settings:
Name
Expression
Unit
Description
plyAngle
(comp1.cfunc1(X/a)*comp1.cfunc2(Y/a))[rad]
rad
Fiber orientation with respect X-axis
4
In the Label text field, type Variables (Angle Ply for Curvilinear Fibers).
Study 1: Original Orientation
Step 1: Stationary
Disable newly added variables node in the study step of Study 1: Original Orientation.
1
In the Model Builder window, under Study 1: Original Orientation click Step 1: Stationary.
2
In the Settings window for Stationary, locate the Physics and Variables Selection section.
3
Select the Modify model configuration for study step checkbox.
4
In the tree, select Global Definitions > Variables (Angle Ply for Curvilinear Fibers).
5
Click  Disable.
Add Study
1
In the Study toolbar, click  Add Study to open the Add Study window.
2
Go to the Add Study window.
3
Find the Studies subsection. In the Select Study tree, select General Studies > Stationary.
4
Click the Add Study button in the window toolbar.
5
In the Study toolbar, click  Add Study to close the Add Study window.
Study 2: Orientation Optimization
In the Settings window for Study, type Study 2: Orientation Optimization in the Label text field.
General Optimization
1
In the Study toolbar, click  Optimization and choose General Optimization.
2
In the Settings window for General Optimization, locate the Optimization Solver section.
3
From the Method list, choose GCMMA.
4
Locate the Objective Function section. In the table, enter the following settings:
Expression
Description
Evaluate for
comp1.lshell.Ws_tot
Total elastic strain energy
Stationary
Step 1: Stationary
1
In the Model Builder window, click Step 1: Stationary.
2
In the Settings window for Stationary, locate the Physics and Variables Selection section.
3
Select the Modify model configuration for study step checkbox.
4
In the tree, select Global Definitions > Variables (Angle Ply for Straight Fibers).
5
Click  Disable.
6
In the Study toolbar, click  Compute.
Results
Layered Material (Optimization)
1
In the Model Builder window, under Results > Datasets click Layered Material.
2
In the Settings window for Layered Material, type Layered Material (Optimization) in the Label text field.
Mirror 3D 3
1
In the Model Builder window, under Results > Datasets right-click Mirror 3D 1 and choose Duplicate.
2
In the Settings window for Mirror 3D, locate the Data section.
3
From the Dataset list, choose Layered Material (Optimization).
Mirror 3D (Optimized Orientation)
1
In the Model Builder window, right-click Mirror 3D (Original Orientation) and choose Duplicate.
2
In the Settings window for Mirror 3D, type Mirror 3D (Optimized Orientation) in the Label text field.
3
Locate the Data section. From the Dataset list, choose Mirror 3D 3.
3D Plot Group 3
In the Results toolbar, click  3D Plot Group.
Stress (Optimized Orientation)
1
In the Model Builder window, under Results click Stress (lshell).
2
In the Settings window for 3D Plot Group, type Stress (Optimized Orientation) in the Label text field.
3
Locate the Data section. From the Dataset list, choose Mirror 3D (Optimized Orientation).
4
Locate the Color Legend section. Select the Show maximum and minimum values checkbox.
5
Select the Show units checkbox.
6
In the Stress (Optimized Orientation) toolbar, click  Plot.
Optimized Fiber Orientation
1
In the Model Builder window, under Results click 3D Plot Group 3.
2
In the Settings window for 3D Plot Group, type Optimized Fiber Orientation in the Label text field.
3
Locate the Data section. From the Dataset list, choose Mirror 3D (Optimized Orientation).
4
Locate the Color Legend section. Select the Show maximum and minimum values checkbox.
5
Click to expand the Title section. From the Title type list, choose Manual.
6
In the Title text area, type Surface: Layer rotation (°).
7
Locate the Color Legend section. Select the Show units checkbox.
Surface 1
1
Right-click Optimized Fiber Orientation and choose Surface.
2
In the Settings window for Surface, locate the Expression section.
3
In the Expression text field, type llmat1.rot.
4
From the Unit list, choose °.
5
Click the  Go to XY View button in the Graphics toolbar.
Optimized Fiber Orientation
Right-click Surface 1 and choose Streamline.
Streamline 1
1
In the Settings window for Streamline, locate the Expression section.
2
In the x-component text field, type lshell.tm11.
3
In the y-component text field, type lshell.tm12.
4
In the z-component text field, type lshell.tm13.
5
Locate the Streamline Positioning section. From the Entry method list, choose Coordinates.
6
In the x text field, type r0.
7
In the y text field, type range(-a,a/50,a).
8
In the z text field, type 0.
9
Locate the Coloring and Style section. Find the Point style subsection. From the Color list, choose Black.
10
Find the Line style subsection. From the Type list, choose Tube.
11
Select the Radius scale factor checkbox.
12
In the Tube radius expression text field, type 1e-3.
Streamline 2
1
Right-click Results > Optimized Fiber Orientation > Streamline 1 and choose Duplicate.
2
In the Settings window for Streamline, locate the Streamline Positioning section.
3
In the x text field, type -r0.
Arrow Line 1
1
In the Model Builder window, right-click Optimized Fiber Orientation and choose Arrow Line.
2
In the Settings window for Arrow Line, locate the Expression section.
3
In the x-component text field, type lshell.bndl1.fax.
4
In the y-component text field, type lshell.bndl1.fay.
5
In the z-component text field, type lshell.bndl1.faz.
6
Locate the Arrow Positioning section. From the Placement list, choose Mesh vertices.
7
Locate the Coloring and Style section.
8
Select the Scale factor checkbox. In the associated text field, type 8e-11.
9
In the Optimized Fiber Orientation toolbar, click  Plot.
Total Elastic Strain Energy
1
In the Results toolbar, click  1D Plot Group.
2
In the Settings window for 1D Plot Group, type Total Elastic Strain Energy in the Label text field.
3
Locate the Legend section. From the Position list, choose Center.
Global 1
1
Right-click Total Elastic Strain Energy and choose Global.
2
In the Settings window for Global, locate the Data section.
3
From the Dataset list, choose Study 1: Original Orientation/Solution 1 (sol1).
4
Locate the y-Axis Data section. In the table, enter the following settings:
Expression
Unit
Description
lshell.Ws_tot
J
Total elastic strain energy
5
Click to expand the Legends section. From the Legends list, choose Manual.
6
In the table, enter the following settings:
Legends
Straight Fibers
Line Segments 1
1
In the Model Builder window, right-click Total Elastic Strain Energy and choose Line Segments.
2
In the Settings window for Line Segments, locate the Data section.
3
From the Dataset list, choose Study 2: Orientation Optimization/Solution 2 (sol2).
4
Locate the x-Coordinates section. In the table, enter the following settings:
Expression
Unit
Description
-90
1
 
90
1
 
5
Locate the y-Coordinates section. In the table, enter the following settings:
Expression
Unit
Description
lshell.Ws_tot
J
Total elastic strain energy
lshell.Ws_tot
J
Total elastic strain energy
6
Click to expand the Coloring and Style section. Find the Line style subsection. From the Line list, choose Dashed.
7
From the Color list, choose From theme.
8
Click to expand the Legends section. Select the Show legends checkbox.
9
From the Legends list, choose Manual.
10
In the table, enter the following settings:
Legends
Variable Angle Tow
11
In the Total Elastic Strain Energy toolbar, click  Plot.
Maximum Displacement
1
Right-click Total Elastic Strain Energy and choose Duplicate.
2
In the Settings window for 1D Plot Group, type Maximum Displacement in the Label text field.
3
Click to expand the Title section. From the Title type list, choose Label.
4
Locate the Plot Settings section.
5
Select the y-axis label checkbox. In the associated text field, type Maximum Displacement (mm).
Global 1
1
In the Model Builder window, expand the Maximum Displacement node, then click Global 1.
2
In the Settings window for Global, locate the y-Axis Data section.
3
In the table, enter the following settings:
Expression
Unit
Description
maxop1(lshell.atxd1(0,lshell.disp))
mm
Maximum 1
Line Segments 1
1
In the Model Builder window, click Line Segments 1.
2
In the Settings window for Line Segments, locate the y-Coordinates section.
3
In the table, enter the following settings:
Expression
Unit
Description
maxop1(lshell.atxd1(0,lshell.disp))
mm
Maximum 1
maxop1(lshell.atxd1(0,lshell.disp))
mm
Maximum 1
4
In the Maximum Displacement toolbar, click  Plot.