### Installing Dependencies for Composite Optimization (Python)

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

This snippet demonstrates how to install the necessary Python packages (numpy, scipy, pymoo, composites) required to run the optimization tutorial. It uses the `pip` package installer via a shell command executed from within a Python environment (like a Jupyter notebook). The output is redirected to a temporary file.

```python
!python -m pip install numpy scipy pymoo composites > tmp.txt
```

--------------------------------

### Install Required Python Modules

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

This snippet installs the necessary Python libraries for running the optimization tutorial. It uses pip to install numpy, scipy, pymoo, and the composites library, redirecting output to a temporary file.

```Python
!python -m pip install numpy scipy pymoo composites > tmp.txt
```

--------------------------------

### Installing Required Python Modules

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Installs the `numpy`, `scipy`, `pymoo`, and `composites` libraries using pip from the command line, redirecting output to a temporary file. These modules are prerequisites for running the optimization code.

```python
!python -m pip install numpy scipy pymoo composites > tmp.txt
```

--------------------------------

### Commented Objective Function Example - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

This is a commented-out function definition (`objective_Workshop1`) that appears to be an example of an objective function for an optimization problem. It would likely calculate buckling and failure loads for a given design vector `x` and combine them into a single objective value to be minimized. It references `calc_buckling_FE` (not defined in this snippet) and `calc_failure_load_Haftka`.

```python
#def objective_Workshop1(x): # NOTE to be minimized
#    stack = discrete_stack_from_continuous_x(x)
#    prop = laminated_plate(stack=stack, plyt=ply_thickness, laminaprop=laminaprop)
#    lambda_cb = calc_buckling_FE(prop)
#    lambda_cs = calc_failure_load_Haftka(prop)
#    p = 0.08 # NOTE from the reference paper
#    obj = (1 - p)*min(lambda_cs, lambda_cb)
```

--------------------------------

### Testing Constraint Functions and Layup - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Provides a test case demonstrating how to use the `layup_thick_from_x`, `laminated_plate`, and `calc_buckling_analytical` functions with a sample design variable vector `test`. It prints the resulting layup stack, compares it to a reference, and calculates and prints the buckling load factor.

```python
test = [2, 1, 1, 1, 1, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
stack, thicknesses = layup_thick_from_x(test)
print(stack.tolist())
print(stack_ref_half)
print(len(stack))
prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
lambda_cb = calc_buckling_analytical(prop)
print(lambda_cb)
```

--------------------------------

### Initializing Laminated Plate Object - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Imports the `laminated_plate` class from the `composites` library and creates an instance named `prop_ref`. This object represents the composite laminate and is initialized with the stacking sequence (`stack_ref`), ply thickness (`ply_thickness`), and material properties (`laminaprop`).

```python
from composites import laminated_plate

prop_ref = laminated_plate(stack=stack_ref, plyt=ply_thickness, laminaprop=laminaprop)
```

--------------------------------

### Verifying Calculation Functions - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Calls the previously defined functions `calc_failure_load_Haftka` and `calc_buckling_analytical` with the reference laminate properties (`prop_ref`). It prints the calculated values and includes an assertion to check if the failure load factor is close to an expected value (13518.66). Requires `numpy`.

```python
lambda_cs_ref = calc_failure_load_Haftka(prop_ref)
print('lambda_cs_ref', lambda_cs_ref)
assert np.isclose(lambda_cs_ref, 13518.66, rtol=1e-3)

lambda_cb_ref_analytical = calc_buckling_analytical(prop_ref)
print('lambda_cb_ref_analytical', lambda_cb_ref_analytical)
```

--------------------------------

### Analyzing Optimal Solution (Python)

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Analyzes the optimal solution obtained from the optimization result ('res.X'). It converts the variable representation back into the laminate stack and thicknesses, calculates the material reserve factors for buckling ('ms_buckling_opt') and failure ('ms_failure_opt') relative to the target lambda, and prints various details of the optimal design including the number of layers, volume, reserve factors, and the optimal layup sequence.

```python
stack_opt, thicknesses_opt = layup_thick_from_x(res.X)
prop_opt = laminated_plate(stack=stack_opt, plyts=thicknesses_opt, laminaprop=laminaprop)
ms_buckling_opt = calc_buckling_analytical(prop_opt)/target_lambda - 1
ms_failure_opt = calc_failure_load_Haftka(prop_opt)/target_lambda - 1
print('num of layers', len(stack_opt))
print('volume', volume(res.X))
print('ms_buckling_opt', ms_buckling_opt)
print('ms_failure_opt', ms_failure_opt)
print('optimal layup', stack_opt)
```

--------------------------------

### Verifying Buckling and Failure Load Calculations - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Calls the previously defined functions `calc_failure_load_Haftka` and `calc_buckling_analytical` using the reference plate properties (`prop_ref`). It prints the calculated failure and buckling load factors and includes an assertion to verify the failure load against a known reference value (13518.66), confirming the calculation's correctness. Requires `numpy` for `isclose`.

```python
lambda_cs_ref = calc_failure_load_Haftka(prop_ref)
print('lambda_cs_ref', lambda_cs_ref)
assert np.isclose(lambda_cs_ref, 13518.66, rtol=1e-3)

lambda_cb_ref_analytical = calc_buckling_analytical(prop_ref)
print('lambda_cb_ref_analytical', lambda_cb_ref_analytical)
```

--------------------------------

### Define Reference Case Geometry, Loading, and Stacking Sequence

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

This snippet sets up the parameters for the reference optimization case, including the plate dimensions, unit loads, and a reference 48-ply symmetric stacking sequence from Table 3 of the cited paper. It also calculates the initial number of independent design variables, doubling it for the ghost layer approach.

```Python
# Table 03, first row
a_value = 20. # [in]
b_value = 5. # [in]
load_Nxx_unit = 1. # [lb]
load_Nyy_unit = 0.125 # [lb]
stack_ref_half = [90]*2 + [+45,-45]*4 + [0]*4 + [+45, -45] + [0]*4 + [+45, -45] + [0]*2 
stack_ref = stack_ref_half + stack_ref_half[::-1] # symmetry condition

assert len(stack_ref) == 48
num_variables = len(stack_ref) // 4 # NOTE independent angles for a symmetric and balanced laminate

# NOTE make sure that you understand this
num_variables = 2*num_variables # NOTE twice as many variables needed for the ghost layer approach
```

--------------------------------

### Initializing Laminated Plate (Python)

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Creates a `laminated_plate` object using the `composites` library. It takes the stacking sequence (`stack_ref`), ply thickness (`ply_thickness`), and lamina properties (`laminaprop`) as input. This object likely contains the laminate's stiffness matrices (ABD).

```python
from composites import laminated_plate

prop_ref = laminated_plate(stack=stack_ref, plyt=ply_thickness, laminaprop=laminaprop)
```

--------------------------------

### Defining Reference Laminate Geometry, Loading, and Variables (Python)

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

This Python snippet sets up the reference data for the optimization task, based on the first design case from Table 3. It defines the plate dimensions (`a_value`, `b_value`), the unit biaxial loading (`load_Nxx_unit`, `load_Nyy_unit`), and the reference stacking sequence (`stack_ref`). It also calculates the initial number of independent variables, doubling it for the ghost layer approach to allow for ply addition/removal.

```python
# Table 03, first row
a_value = 20. # [in]
b_value = 5. # [in]
load_Nxx_unit = 1. # [lb]
load_Nyy_unit = 0.125 # [lb]
stack_ref_half = [90]*2 + [+45,-45]*4 + [0]*4 + [+45, -45] + [0]*4 + [+45, -45] + [0]*2
stack_ref = stack_ref_half + stack_ref_half[::-1] # symmetry condition

assert len(stack_ref) == 48
num_variables = len(stack_ref) // 4 # NOTE independent angles for a symmetric and balanced laminate

# NOTE make sure that you understand this
num_variables = 2*num_variables # NOTE twice as many variables needed for the ghost layer approach
```

--------------------------------

### Analyzing Optimization Results - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Analyzes the results obtained from the optimization (`res`). It extracts the optimal layup and thicknesses from the best solution (`res.X`), calculates the material reserve factors (margin of safety + 1) for buckling and failure based on the optimal layup, and prints the number of layers, volume, material reserve factors, and the optimal layup stack.

```python
stack_opt, thicknesses_opt = layup_thick_from_x(res.X)
prop_opt = laminated_plate(stack=stack_opt, plyts=thicknesses_opt, laminaprop=laminaprop)
ms_buckling_opt = calc_buckling_analytical(prop_opt)/target_lambda - 1
ms_failure_opt = calc_failure_load_Haftka(prop_opt)/target_lambda - 1
print('num of layers', len(stack_opt))
print('volume', volume(res.X))
print('ms_buckling_opt', ms_buckling_opt)
print('ms_failure_opt', ms_failure_opt)
print('optimal layup', stack_opt)
```

--------------------------------

### Running PyMoo Optimization - Python

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Executes the optimization process using the `pymoo.minimize` function. It takes the defined problem and method as input, sets the termination criteria (number of generations), a random seed for reproducibility, and enables saving the history. Finally, it prints the best solution found, its objective function value, and constraint violation.

```python
res = minimize(problem,
               method,
               termination=('n_gen', 40),
               seed=1,
               save_history=True
               )

print("Best solution found: %s" % res.X)
print("Function value: %s" % res.F)
print("Constraint violation: %s" % res.CV)
```

--------------------------------

### Initializing Laminated Plate Object - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Imports the `laminated_plate` class from the `composites` library and initializes a plate object. This object is created using a predefined stack sequence (`stack_ref`), the calculated ply thickness, and the defined lamina properties (`laminaprop`). The resulting object (`prop_ref`) contains the plate's stiffness matrices (ABD).

```python
from composites import laminated_plate

prop_ref = laminated_plate(stack=stack_ref, plyt=ply_thickness, laminaprop=laminaprop)
```

--------------------------------

### Analyzing Optimization Results - Python

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Analyzes the results obtained from the optimization. It extracts the optimal layup and thicknesses from the result, calculates the margin of safety for buckling and failure constraints based on the optimal properties, and prints the number of layers, volume, margins of safety, and the optimal layup.

```python
stack_opt, thicknesses_opt = layup_thick_from_x(res.X)
prop_opt = laminated_plate(stack=stack_opt, plyts=thicknesses_opt, laminaprop=laminaprop)
ms_buckling_opt = calc_buckling_analytical(prop_opt)/target_lambda - 1
ms_failure_opt = calc_failure_load_Haftka(prop_opt)/target_lambda - 1
print('num of layers', len(stack_opt))
print('volume', volume(res.X))
print('ms_buckling_opt', ms_buckling_opt)
print('ms_failure_opt', ms_failure_opt)
print('optimal layup', stack_opt)
```

--------------------------------

### Defining Geometry, Load, and Variables for Composite Optimization in Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Initializes geometric dimensions (`a_value`, `b_value`), unit loads (`load_Nxx_unit`, `load_Nyy_unit`), and a reference 48-ply stacking sequence (`stack_ref`) based on Table 3. It then calculates the number of independent variables for a symmetric and balanced laminate and doubles it for the ghost layer approach.

```python
# Table 03, first row
a_value = 20. # [in]
b_value = 5. # [in]
load_Nxx_unit = 1. # [lb]
load_Nyy_unit = 0.125 # [lb]
stack_ref_half = [90]*2 + [+45,-45]*4 + [0]*4 + [+45, -45] + [0]*4 + [+45, -45] + [0]*2
stack_ref = stack_ref_half + stack_ref_half[::-1] # symmetry condition

assert len(stack_ref) == 48
num_variables = len(stack_ref) // 4 # NOTE independent angles for a symmetric and balanced laminate

# NOTE make sure that you understand this
num_variables = 2*num_variables # NOTE twice as many variables needed for the ghost layer approach
```

--------------------------------

### Verifying Failure and Buckling Loads (Python)

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Calls the previously defined functions `calc_failure_load_Haftka` and `calc_buckling_analytical` with a reference plate property object (`prop_ref`). It prints the calculated failure load (`lambda_cs_ref`) and buckling load (`lambda_cb_ref_analytical`) and includes an assertion to check if the failure load matches an expected reference value.

```python
lambda_cs_ref = calc_failure_load_Haftka(prop_ref)
print('lambda_cs_ref', lambda_cs_ref)
assert np.isclose(lambda_cs_ref, 13518.66, rtol=1e-3)

lambda_cb_ref_analytical = calc_buckling_analytical(prop_ref)
print('lambda_cb_ref_analytical', lambda_cb_ref_analytical)

```

--------------------------------

### Executing Optimization with pymoo (Python)

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Executes the optimization process using the configured problem and method. It sets a termination criterion ('n_gen', 40 generations), a random seed for reproducibility, and saves the history of the optimization run. Finally, it prints the best solution found ('res.X'), the corresponding objective function value ('res.F'), and the constraint violation ('res.CV').

```python
res = minimize(problem,
               method,
               termination=('n_gen', 40),
               seed=1,
               save_history=True
               )

print("Best solution found: %s" % res.X)
print("Function value: %s" % res.F)
print("Constraint violation: %s" % res.CV)
```

--------------------------------

### Configuring Genetic Algorithm with pymoo (Python)

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Configures the optimization algorithm using pymoo's Genetic Algorithm (GA). It sets parameters like population size ('pop_size'), sampling method ('IntegerRandomSampling'), crossover operator ('SBX') with specific probability and eta values, mutation operator ('PM') with specific probability and eta values, and enables elimination of duplicate solutions.

```python
method = algorithm = GA(
    pop_size=20,
    sampling=IntegerRandomSampling(),
    crossover=SBX(prob=1.0, eta=3.0, vtype=float, repair=RoundingRepair()),
    mutation=PM(prob=1.0, eta=3.0, vtype=float, repair=RoundingRepair()),
    eliminate_duplicates=True)
```

--------------------------------

### Configuring Genetic Algorithm Method - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Configures the optimization algorithm using PyMoo's Genetic Algorithm (`GA`). It sets parameters such as population size (`pop_size`), sampling method (`IntegerRandomSampling`), crossover operator (`SBX`) with rounding repair, mutation operator (`PM`) with rounding repair, and duplicate elimination.

```python
method = algorithm = GA(
    pop_size=20,
    sampling=IntegerRandomSampling(),
    crossover=SBX(prob=1.0, eta=3.0, vtype=float, repair=RoundingRepair()),
    mutation=PM(prob=1.0, eta=3.0, vtype=float, repair=RoundingRepair()),
    eliminate_duplicates=True)
```

--------------------------------

### Executing PyMoo Optimization - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Executes the optimization process using the `minimize` function from PyMoo. It takes the defined `problem` and `method` as input, sets a termination criterion based on the number of generations (`n_gen`), specifies a random seed for reproducibility, and enables saving the optimization history. Finally, it prints the best solution found, its objective function value, and constraint violation.

```python
res = minimize(problem,
               method,
               termination=('n_gen', 40),
               seed=1,
               save_history=True
               )

print("Best solution found: %s" % res.X)
print("Function value: %s" % res.F)
print("Constraint violation: %s" % res.CV)
```

--------------------------------

### Importing Dependencies for Composite Analysis - Python

Source: https://github.com/saullocastro/composites/blob/master/doc/source/post-polar-plots-laminate-stiffness.ipynb

Imports the required libraries `numpy` for numerical operations, `matplotlib.pyplot` for plotting, and the `laminated_plate` class from the `composites` library to model composite laminates. These are essential prerequisites for the subsequent calculations and plotting.

```Python
import numpy as np
import matplotlib.pyplot as plt

from composites import laminated_plate
```

--------------------------------

### Generating Layup and Thickness from Vector - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Defines a function `layup_thick_from_x` that takes an integer vector `x` as input and generates a balanced and symmetric composite stacking sequence (`layup_bal`) and corresponding ply thicknesses (`plyts_bal`). The first half of the vector `x` determines the angles from `available_angles`, and the second half determines the thickness multipliers. Requires `numpy`.

```python
available_angles = [0, 45, 90]
ply_thickness = 0.005*25.4/1000 # [m]

def layup_thick_from_x(x):
    layup = []
    plyts = []
    for xi in x[0:len(x)//2]:
        layup.append(available_angles[xi])
    for xi in x[len(x)//2:]:
        plyts.append(ply_thickness*xi)
    layup_bal = []
    plyts_bal = []
    for angle_deg, plyt in zip(layup, plyts): # balancing
        layup_bal.append(angle_deg)
        plyts_bal.append(plyt)
        if angle_deg != 0 and angle_deg != 90:
            layup_bal.append(-angle_deg)
            plyts_bal.append(plyt)
        else:
            layup_bal.append(angle_deg)
            plyts_bal.append(plyt)
    layup_bal = layup_bal + layup_bal[::-1] # symmetry
    plyts_bal = plyts_bal + plyts_bal[::-1] # symmetry
    non_zero = np.logical_not(np.asarray(plyts_bal) == 0)
    return np.asarray(layup_bal)[non_zero], np.asarray(plyts_bal)[non_zero]
```

--------------------------------

### Defining and Testing Constraint Functions - Python

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Defines two constraint functions, `calc_constr_buckling` and `calc_constr_failure`, which return values less than or equal to 0 for feasible designs. It also includes test code to evaluate these functions and related properties for a sample layup.

```python
def calc_constr_buckling(x): # NOTE feasible when <= 0
    stack, thicknesses = layup_thick_from_x(x)
    prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
    lambda_cb = calc_buckling_analytical(prop)
    return 1 - lambda_cb/target_lambda

def calc_constr_failure(x): # NOTE feasible when <= 0
    stack, thicknesses = layup_thick_from_x(x)
    prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
    lambda_cs = calc_failure_load_Haftka(prop)
    return 1 - lambda_cs/target_lambda

test = [2, 1, 1, 1, 1, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
stack, thicknesses = layup_thick_from_x(test)
print(stack.tolist())
print(stack_ref_half)
print(len(stack))
prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
lambda_cb = calc_buckling_analytical(prop)
print(lambda_cb)
```

--------------------------------

### Calculating Buckling Load Analytically (Python)

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Defines a function `calc_buckling_analytical` that calculates the critical buckling load for a laminated plate using an analytical formula. It iterates through possible half-wave numbers (m, n) to find the minimum buckling load factor (`lambda_b_min`). Requires the plate properties (`prop`) object as input.

```python
from numpy import pi

def calc_buckling_analytical(prop):
    a = a_value*25.4/1000 # [m] along x
    b = b_value*25.4/1000 # [m] along y
    Nxx = load_Nxx_unit*4.448222/(25.4/1000) #[N/m]
    Nyy = load_Nyy_unit*4.448222/(25.4/1000) #[N/m]
    D11 = prop.ABD[3, 3]
    D12 = prop.ABD[3, 4]
    D22 = prop.ABD[4, 4]
    D66 = prop.ABD[5, 5]
    lambda_b_min = 1e30
    for m in range(1, 21):
        for n in range(1, 21):
            lambda_b = (pi**2*(D11*(m/a)**4
                             + 2*(D12 + 2*D66)*(m/a)**2*(n/b)**2
                             + D22*(n/b)**4)/((m/a)**2*Nxx + (n/b)**2*Nyy)
                       )
            lambda_b_min = min(lambda_b_min, lambda_b)
    return lambda_b_min
```

--------------------------------

### Defining Optimization Problem with PyMoo - Python

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Defines a custom `ElementwiseProblem` class for the optimization using PyMoo. It sets the number of variables, objectives, constraints, variable types, and bounds. The `_evaluate` method calculates the objective function (volume) and constraint violations (buckling and failure) for a given set of design variables `x`.

```python
from pymoo.core.problem import ElementwiseProblem
from pymoo.optimize import minimize
from pymoo.algorithms.soo.nonconvex.ga import GA
from pymoo.operators.sampling.rnd import IntegerRandomSampling
from pymoo.operators.crossover.sbx import SBX
from pymoo.operators.mutation.pm import PM
from pymoo.operators.repair.rounding import RoundingRepair


n_plies = 48
n_var = 2*(n_plies//4)

# define permutation problem
class Problem(ElementwiseProblem):
    def __init__(self):
        super().__init__(n_var = n_var,
                         n_obj = 1,
                         xl    = [0]*n_var, # NOTE variables always changing from 0 to xu
                         xu    = [len(available_angles)-1]*(n_var//2) + [1]*(n_var//2),
                         n_ieq_constr = 2,
                         vtype = int)
    
    def _evaluate(self, x, out, *args, **kwargs):
        x = np.asarray(x, dtype=int)
        layup, plyts = layup_thick_from_x(x)
        out['F'] = volume(x) # objective function
        out['G'] = [calc_constr_buckling(x), calc_constr_failure(x)]
    

problem = Problem()
```

--------------------------------

### Configuring Genetic Algorithm Method - PyMoo Python

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Configures the Genetic Algorithm (GA) method to be used for the optimization. It specifies parameters such as population size, sampling strategy (IntegerRandomSampling), crossover operator (SBX with RoundingRepair), mutation operator (PM with RoundingRepair), and duplicate elimination.

```python
method = algorithm = GA(
    pop_size=20,
    sampling=IntegerRandomSampling(),
    crossover=SBX(prob=1.0, eta=3.0, vtype=float, repair=RoundingRepair()),
    mutation=PM(prob=1.0, eta=3.0, vtype=float, repair=RoundingRepair()),
    eliminate_duplicates=True)
```

--------------------------------

### Defining Buckling and Failure Constraints - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Defines two constraint functions, `calc_constr_buckling` and `calc_constr_failure`. These functions calculate the margin of safety for buckling and failure, respectively, based on the laminate properties derived from the design variables `x`. The constraints are formulated such that they are feasible when the return value is less than or equal to 0.

```python
def calc_constr_buckling(x): # NOTE feasible when <= 0
    stack, thicknesses = layup_thick_from_x(x)
    prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
    lambda_cb = calc_buckling_analytical(prop)
    return 1 - lambda_cb/target_lambda

def calc_constr_failure(x): # NOTE feasible when <= 0
    stack, thicknesses = layup_thick_from_x(x)
    prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
    lambda_cs = calc_failure_load_Haftka(prop)
    return 1 - lambda_cs/target_lambda
```

--------------------------------

### Defining and Testing Constraint Functions (Python)

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Defines two inequality constraint functions: one for buckling ('calc_constr_buckling') and one for failure ('calc_constr_failure'). These functions return a value <= 0 for a feasible solution. They calculate the buckling and failure loads based on the laminate properties derived from the layup variables (x) and compare them against a 'target_lambda'. A test case demonstrates how to use these functions and related utility functions ('layup_thick_from_x', 'laminated_plate').

```python
def calc_constr_buckling(x): # NOTE feasible when <= 0
    stack, thicknesses = layup_thick_from_x(x)
    prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
    lambda_cb = calc_buckling_analytical(prop)
    return 1 - lambda_cb/target_lambda

def calc_constr_failure(x): # NOTE feasible when <= 0
    stack, thicknesses = layup_thick_from_x(x)
    prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
    lambda_cs = calc_failure_load_Haftka(prop)
    return 1 - lambda_cs/target_lambda

test = [2, 1, 1, 1, 1, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
stack, thicknesses = layup_thick_from_x(test)
print(stack.tolist())
print(stack_ref_half)
print(len(stack))
prop = laminated_plate(stack=stack, plyts=thicknesses, laminaprop=laminaprop)
lambda_cb = calc_buckling_analytical(prop)
print(lambda_cb)
```

--------------------------------

### Defining Optimization Problem with pymoo (Python)

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Defines the optimization problem by subclassing pymoo's 'ElementwiseProblem'. It sets the number of variables ('n_var'), objectives ('n_obj'), lower ('xl') and upper ('xu') bounds for variables, number of inequality constraints ('n_ieq_constr'), and variable type ('vtype'). The '_evaluate' method calculates the objective function ('out['F']') using the 'volume' function and the constraint violations ('out['G']') using the previously defined constraint functions.

```python
from pymoo.core.problem import ElementwiseProblem
from pymoo.optimize import minimize
from pymoo.algorithms.soo.nonconvex.ga import GA
from pymoo.operators.sampling.rnd import IntegerRandomSampling
from pymoo.operators.crossover.sbx import SBX
from pymoo.operators.mutation.pm import PM
from pymoo.operators.repair.rounding import RoundingRepair


n_plies = 48
n_var = 2*(n_plies//4)

# define permutation problem
class Problem(ElementwiseProblem):
    def __init__(self):
        super().__init__(n_var = n_var,
                         n_obj = 1,
                         xl    = [0]*n_var, # NOTE variables always changing from 0 to xu
                         xu    = [len(available_angles)-1]*(n_var//2) + [1]*(n_var//2),
                         n_ieq_constr = 2,
                         vtype = int)
    
    def _evaluate(self, x, out, *args, **kwargs):
        x = np.asarray(x, dtype=int)
        layup, plyts = layup_thick_from_x(x)
        out['F'] = volume(x) # objective function
        out['G'] = [calc_constr_buckling(x), calc_constr_failure(x)]
    

problem = Problem()
```

--------------------------------

### Calculating Analytical Buckling Load - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Defines a function `calc_buckling_analytical` that calculates the critical buckling load factor for a laminated plate based on analytical formulas. It takes a plate property object (`prop`), plate dimensions (`a_value`, `b_value`), and applied loads (`load_Nxx_unit`, `load_Nyy_unit`) as input. It iterates through potential buckling modes (m, n) to find the minimum load factor. Requires `numpy` for pi.

```python
from numpy import pi

def calc_buckling_analytical(prop):
    a = a_value*25.4/1000 # [m] along x
    b = b_value*25.4/1000 # [m] along y
    Nxx = load_Nxx_unit*4.448222/(25.4/1000) #[N/m]
    Nyy = load_Nyy_unit*4.448222/(25.4/1000) #[N/m]
    D11 = prop.ABD[3, 3]
    D12 = prop.ABD[3, 4]
    D22 = prop.ABD[4, 4]
    D66 = prop.ABD[5, 5]
    lambda_b_min = 1e30
    for m in range(1, 21):
        for n in range(1, 21):
            lambda_b = (pi**2*(D11*(m/a)**4
                             + 2*(D12 + 2*D66)*(m/a)**2*(n/b)**2
                             + D22*(n/b)**4)/((m/a)**2*Nxx + (n/b)**2*Nyy)
                       )
            lambda_b_min = min(lambda_b_min, lambda_b)
    return lambda_b_min
```

--------------------------------

### Calculating Strain Failure Load (Python)

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Defines a function `calc_failure_load_Haftka` that calculates the load factor leading to strain failure using an optimization approach. It defines an inner function `strain_MS` to calculate the margin of safety based on applied loads and laminate stiffness (A matrix). It uses `scipy.optimize.minimize` to find the load factor (`lbd`) where the margin of safety is zero. Requires the plate properties (`prop`) object as input.

```python
import numpy as np
import scipy.optimize as opt

def calc_failure_load_Haftka(prop):
    def strain_MS(lbd):
        Nxx = -load_Nxx_unit*4.448222/(25.4/1000) #[N/m]
        Nyy = -load_Nyy_unit*4.448222/(25.4/1000) #[N/m]
        Nxx = lbd*Nxx*1.5
        Nyy = lbd*Nyy*1.5
        vecN = np.asarray([Nxx, Nyy])
        A11 = prop.ABD[0, 0]
        A12 = prop.ABD[0, 1]
        A22 = prop.ABD[1, 1]
        exx, eyy = np.linalg.inv(np.array([[A11, A12], [A12, A22]])) @ vecN
        margin_of_safety = 1e15
        for thetadeg in prop.stack:
            cost = np.cos(np.deg2rad(thetadeg))
            sint = np.sin(np.deg2rad(thetadeg))
            epsilon_i_1 = cost**2*exx + sint**2*eyy
            epsilon_i_2 = sint**2*exx + cost**2*eyy
            gamma_i_12 = sint**2*(eyy - exx)
            ms_new = min(
                    epsilon_1_allowable/abs(epsilon_i_1) - 1,
                    epsilon_2_allowable/abs(epsilon_i_2) - 1,
                    gamma_12_allowable/abs(gamma_i_12) - 1
                    )
            margin_of_safety = min(margin_of_safety, ms_new)
        return margin_of_safety
    lbd_init = 100
    positiveMS = opt.NonlinearConstraint(strain_MS, 0., np.inf, jac='2-point')
    res = opt.minimize(strain_MS, lbd_init, tol=1e-6, bounds=((100, None),),
            constraints=[positiveMS], jac='2-point')
    assert res.success
    return res.x[0]
```

--------------------------------

### Calculating Strain Failure Load Factor - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Defines a function `calc_failure_load_Haftka` that calculates the load factor at which the composite laminate fails based on a strain criterion. It uses `scipy.optimize.minimize` to find the load factor (`lbd`) where the minimum margin of safety (`strain_MS`) is zero. It requires `numpy` and `scipy.optimize`.

```python
import numpy as np
import scipy.optimize as opt

def calc_failure_load_Haftka(prop):
    def strain_MS(lbd):
        Nxx = -load_Nxx_unit*4.448222/(25.4/1000) #[N/m]
        Nyy = -load_Nyy_unit*4.448222/(25.4/1000) #[N/m]
        Nxx = lbd*Nxx*1.5
        Nyy = lbd*Nyy*1.5
        vecN = np.asarray([Nxx, Nyy])
        A11 = prop.ABD[0, 0]
        A12 = prop.ABD[0, 1]
        A22 = prop.ABD[1, 1]
        exx, eyy = np.linalg.inv(np.array([[A11, A12], [A12, A22]])) @ vecN
        margin_of_safety = 1e15
        for thetadeg in prop.stack:
            cost = np.cos(np.deg2rad(thetadeg))
            sint = np.sin(np.deg2rad(thetadeg))
            epsilon_i_1 = cost**2*exx + sint**2*eyy
            epsilon_i_2 = sint**2*exx + cost**2*eyy
            gamma_i_12 = sint**2*(eyy - exx)
            ms_new = min(
                    epsilon_1_allowable/abs(epsilon_i_1) - 1,
                    epsilon_2_allowable/abs(epsilon_i_2) - 1,
                    gamma_12_allowable/abs(gamma_i_12) - 1
                    )
            margin_of_safety = min(margin_of_safety, ms_new)
        return margin_of_safety
    lbd_init = 100
    positiveMS = opt.NonlinearConstraint(strain_MS, 0., np.inf, jac='2-point')
    res = opt.minimize(strain_MS, lbd_init, tol=1e-6, bounds=((100, None),),
            constraints=[positiveMS], jac='2-point')
    assert res.success
    return res.x[0]
```

--------------------------------

### (Commented) Objective Function for Optimization - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

This commented-out function `objective_Workshop1` outlines a potential objective function for a composite optimization problem. It would take a design vector `x`, convert it to a stack, calculate buckling and failure loads (using `calc_buckling_FE` and `calc_failure_load_Haftka`), and combine them into a single value to be minimized, potentially weighted by a factor `p`. It depends on `laminated_plate` and the calculation functions.

```python
#def objective_Workshop1(x): # NOTE to be minimized
#    stack = discrete_stack_from_continuous_x(x)
#    prop = laminated_plate(stack=stack, plyt=ply_thickness, laminaprop=laminaprop)
#    lambda_cb = calc_buckling_FE(prop)
#    lambda_cs = calc_failure_load_Haftka(prop)
#    p = 0.08 # NOTE from the reference paper
#    obj = (1 - p)*min(lambda_cs, lambda_cb)
```

--------------------------------

### Calculating Analytical Buckling Load - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Defines a function `calc_buckling_analytical` that calculates the critical buckling load factor for a laminated plate using an analytical formula based on plate dimensions (`a`, `b`), applied loads (`Nxx`, `Nyy`), and the plate's bending stiffness matrix (`D`). It iterates through mode shapes (m, n) to find the minimum buckling load factor. Requires the `numpy` library.

```python
from numpy import pi

def calc_buckling_analytical(prop):
    a = a_value*25.4/1000 # [m] along x
    b = b_value*25.4/1000 # [m] along y
    Nxx = load_Nxx_unit*4.448222/(25.4/1000) #[N/m]
    Nyy = load_Nyy_unit*4.448222/(25.4/1000) #[N/m]
    D11 = prop.ABD[3, 3]
    D12 = prop.ABD[3, 4]
    D22 = prop.ABD[4, 4]
    D66 = prop.ABD[5, 5]
    lambda_b_min = 1e30
    for m in range(1, 21):
        for n in range(1, 21):
            lambda_b = (pi**2*(D11*(m/a)**4
                             + 2*(D12 + 2*D66)*(m/a)**2*(n/b)**2
                             + D22*(n/b)**4)/((m/a)**2*Nxx + (n/b)**2*Nyy)
                       )
            lambda_b_min = min(lambda_b_min, lambda_b)
    return lambda_b_min
```

--------------------------------

### Defining Volume Objective Function (Python)

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Defines the objective function to be minimized. It calculates the volume of the composite laminate based on the layup variables (x), using predefined dimensions 'a_value' and 'b_value' and the calculated thicknesses from 'layup_thick_from_x'. The volume is returned in cubic meters.

```python
def volume(x): # NOTE to be minimized
    a = a_value*25.4/1000 # [m] along x
    b = b_value*25.4/1000 # [m] along y
    stack, thicknesses = layup_thick_from_x(x)
    return sum(thicknesses)*a*b # [m^3]
```

--------------------------------

### Defining Allowable Strain Limits in Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt2.ipynb

Sets the maximum allowable strain values for the principal strains (epsilon_1, epsilon_2) and shear strain (gamma_12) based on values provided in the reference paper. These values are used for the strain failure constraint in the optimization problem.

```python
# NOTE allowable strains from the reference paper
epsilon_1_allowable = 0.008
epsilon_2_allowable = 0.029
gamma_12_allowable = 0.015
```

--------------------------------

### Defining PyMoo Optimization Problem - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Defines a custom optimization problem class `Problem` inheriting from `pymoo.core.problem.ElementwiseProblem`. It sets the number of variables (`n_var`), objectives (`n_obj`), lower/upper bounds (`xl`, `xu`), number of inequality constraints (`n_ieq_constr`), and variable type (`vtype`). The `_evaluate` method calculates the objective function (`volume`) and constraint violations (`calc_constr_buckling`, `calc_constr_failure`) for a given design vector `x`.

```python
from pymoo.core.problem import ElementwiseProblem
from pymoo.optimize import minimize
from pymoo.algorithms.soo.nonconvex.ga import GA
from pymoo.operators.sampling.rnd import IntegerRandomSampling
from pymoo.operators.crossover.sbx import SBX
from pymoo.operators.mutation.pm import PM
from pymoo.operators.repair.rounding import RoundingRepair


n_plies = 48
n_var = 2*(n_plies//4)

# define permutation problem
class Problem(ElementwiseProblem):
    def __init__(self):
        super().__init__(n_var = n_var,
                         n_obj = 1,
                         xl    = [0]*n_var, # NOTE variables always changing from 0 to xu
                         xu    = [len(available_angles)-1]*(n_var//2) + [1]*(n_var//2),
                         n_ieq_constr = 2,
                         vtype = int)
    
    def _evaluate(self, x, out, *args, **kwargs):
        x = np.asarray(x, dtype=int)
        layup, plyts = layup_thick_from_x(x)
        out['F'] = volume(x) # objective function
        out['G'] = [calc_constr_buckling(x), calc_constr_failure(x)]
    

problem = Problem()
```

--------------------------------

### Defining Volume Objective Function - Python

Source: https://github.com/saullocastro/composites/blob/master/tutorials/opt1.ipynb

Defines the objective function `volume` to be minimized during optimization. It calculates the volume of the composite laminate based on its dimensions (`a_value`, `b_value`) and the total thickness derived from the design variables `x`.

```python
def volume(x): # NOTE to be minimized
    a = a_value*25.4/1000 # [m] along x
    b = b_value*25.4/1000 # [m] along y
    stack, thicknesses = layup_thick_from_x(x)
    return sum(thicknesses)*a*b # [m^3]
```

--------------------------------

### Defining Target Load Factor (Python)

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Sets a variable `target_lambda` to a specific value (10000). This variable likely represents a design load factor used in subsequent calculations or optimization processes.

```python
target_lambda = 10000
```

--------------------------------

### Generating Layup and Thickness from Vector (Python)

Source: https://github.com/saullocastro/composites/blob/master/doc/source/opt-lightweight-discrete-ghost-layer.ipynb

Defines a function `layup_thick_from_x` that takes an integer vector `x` as input and generates a balanced and symmetric composite layup and corresponding ply thicknesses. The first half of the vector defines angles from `available_angles`, and the second half defines thickness multipliers. It enforces balancing and symmetry rules and removes zero-thickness plies.

```python
available_angles = [0, 45, 90]
ply_thickness = 0.005*25.4/1000 # [m]

def layup_thick_from_x(x):
    layup = []
    plyts = []
    for xi in x[0:len(x)//2]:
        layup.append(available_angles[xi])
    for xi in x[len(x)//2:]:
        plyts.append(ply_thickness*xi)
    layup_bal = []
    plyts_bal = []
    for angle_deg, plyt in zip(layup, plyts): # balancing
        layup_bal.append(angle_deg)
        plyts_bal.append(plyt)
        if angle_deg != 0 and angle_deg != 90:
            layup_bal.append(-angle_deg)
            plyts_bal.append(plyt)
        else:
            layup_bal.append(angle_deg)
            plyts_bal.append(plyt)
    layup_bal = layup_bal + layup_bal[::-1] # symmetry
    plyts_bal = plyts_bal + plyts_bal[::-1] # symmetry
    non_zero = np.logical_not(np.asarray(plyts_bal) == 0)
    return np.asarray(layup_bal)[non_zero], np.asarray(plyts_bal)[non_zero]
```