Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Wrapper for AFLR4 AIM in Caps2tacs #197

Merged
merged 6 commits into from
Apr 10, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion examples/caps_wing/1_steady_analysis.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@

comm = MPI.COMM_WORLD
tacs_model = caps2tacs.TacsModel.build(csm_file="simple_naca_wing.csm", comm=comm)
tacs_model.egads_aim.set_mesh(
tacs_model.mesh_aim.set_mesh(
edge_pt_min=15,
edge_pt_max=20,
global_mesh_size=0.25,
Expand Down
2 changes: 1 addition & 1 deletion examples/caps_wing/2_unsteady_analysis.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@
# 1: build the tacs aim, egads aim wrapper classes
comm = MPI.COMM_WORLD
tacs_model = caps2tacs.TacsModel.build(csm_file="simple_naca_wing.csm", comm=comm)
tacs_model.egads_aim.set_mesh(
tacs_model.mesh_aim.set_mesh(
edge_pt_min=15,
edge_pt_max=20,
global_mesh_size=0.25,
Expand Down
2 changes: 1 addition & 1 deletion examples/caps_wing/3_sizing_optimization.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@
comm = MPI.COMM_WORLD
# can also switch to large_naca_wing.csm file here if you want and it will automatically update the DVs
tacs_model = caps2tacs.TacsModel.build(csm_file="large_naca_wing.csm", comm=comm)
tacs_model.egads_aim.set_mesh( # need a refined-enough mesh for the derivative test to pass
tacs_model.mesh_aim.set_mesh( # need a refined-enough mesh for the derivative test to pass
edge_pt_min=15,
edge_pt_max=20,
global_mesh_size=0.01,
Expand Down
2 changes: 1 addition & 1 deletion examples/caps_wing/4_sizing_and_shape.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@
# --------------------------------------------------------------#
comm = MPI.COMM_WORLD
tacs_model = caps2tacs.TacsModel.build(csm_file="large_naca_wing.csm", comm=comm)
tacs_model.egads_aim.set_mesh( # need a refined-enough mesh for the derivative test to pass
tacs_model.mesh_aim.set_mesh( # need a refined-enough mesh for the derivative test to pass
edge_pt_min=15,
edge_pt_max=20,
global_mesh_size=0.01,
Expand Down
106 changes: 106 additions & 0 deletions examples/caps_wing/5_steady_aflr_analysis.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,106 @@
"""
Sean Engelstad, March 2023
GT SMDO Lab, Dr. Graeme Kennedy
Caps to TACS example
"""

from tacs import caps2tacs
from mpi4py import MPI

# run a steady elastic structural analysis in TACS using the tacsAIM wrapper caps2tacs submodule
# -------------------------------------------------------------------------------------------------
# 1: build the tacs aim, egads aim wrapper classes

comm = MPI.COMM_WORLD
tacs_model = caps2tacs.TacsModel.build(
csm_file="simple_naca_wing.csm", comm=comm, mesh="aflr"
)
tacs_model.mesh_aim.set_mesh(
ff_growth=1.4, min_scale=0.05, max_scale=0.5, use_quad=True
).register_to(tacs_model)

aluminum = caps2tacs.Isotropic.aluminum().register_to(tacs_model)

# setup the thickness design variables + automatic shell properties
nribs = int(tacs_model.get_config_parameter("nribs"))
nspars = int(tacs_model.get_config_parameter("nspars"))
for irib in range(1, nribs + 1):
caps2tacs.ThicknessVariable(
caps_group=f"rib{irib}", value=0.05, material=aluminum
).register_to(tacs_model)
for ispar in range(1, nspars + 1):
caps2tacs.ThicknessVariable(
caps_group=f"spar{ispar}", value=0.05, material=aluminum
).register_to(tacs_model)
caps2tacs.ThicknessVariable(
caps_group="OML", value=0.03, material=aluminum
).register_to(tacs_model)

# add constraints and loads
caps2tacs.PinConstraint("root").register_to(tacs_model)
caps2tacs.GridForce("OML", direction=[0, 0, 1.0], magnitude=100).register_to(tacs_model)

# add analysis functions to the model
caps2tacs.AnalysisFunction.ksfailure(ksWeight=50.0, safetyFactor=1.5).register_to(
tacs_model
)
caps2tacs.AnalysisFunction.mass().register_to(tacs_model)

# run the pre analysis to build tacs input files
# alternative is to call tacs_aim.setup_aim().pre_analysis() with tacs_aim = tacs_model.tacs_aim
tacs_model.setup(include_aim=True)

# ----------------------------------------------------------------------------------
# 2. Run the TACS steady elastic structural analysis, forward + adjoint

# choose method 1 or 2 to demonstrate the analysis : 1 - fewer lines, 2 - more lines
method = 2

# show both ways of performing the structural analysis in more or less detail
if method == 1: # less detail version
tacs_model.pre_analysis()
tacs_model.run_analysis()
tacs_model.post_analysis()

print("Tacs model static analysis outputs...\n")
for func in tacs_model.analysis_functions:
print(f"func {func.name} = {func.value}")

for var in tacs_model.variables:
derivative = func.get_derivative(var)
print(f"\td{func.name}/d{var.name} = {derivative}")
print("\n")

elif method == 2: # longer way that directly uses pyTACS & static problem routines
SPs = tacs_model.createTACSProbs(addFunctions=True)

# solve each structural analysis problem (in this case 1)
tacs_funcs = {}
tacs_sens = {}
function_names = tacs_model.function_names
for caseID in SPs:
SPs[caseID].solve()
SPs[caseID].evalFunctions(tacs_funcs, evalFuncs=function_names)
SPs[caseID].evalFunctionsSens(tacs_sens, evalFuncs=function_names)
SPs[caseID].writeSolution(
baseName="tacs_output", outputDir=tacs_model.analysis_dir
)

# print the output analysis functions and sensitivities
print("\nTacs Analysis Outputs...")
for func_name in function_names:
# find the associated tacs key (tacs key = (loadset) + (func_name))
for tacs_key in tacs_funcs:
if func_name in tacs_key:
break

func_value = tacs_funcs[tacs_key].real
print(f"\tfunctional {func_name} = {func_value}")

struct_sens = tacs_sens[tacs_key]["struct"]
print(f"\t{func_name} structDV gradient = {struct_sens}")

xpts_sens = tacs_sens[tacs_key]["Xpts"]
print(f"\t{func_name} coordinate derivatives = {xpts_sens}")

print("\n")
18 changes: 16 additions & 2 deletions examples/caps_wing/simple_naca_wing.csm
Original file line number Diff line number Diff line change
Expand Up @@ -101,6 +101,7 @@ patbeg ispar nspars
# add caps attributes
select face
attribute capsGroup !$spar+ispar
ATTRIBUTE AFLR_GBC $TRANSP_UG3_GBC
attribute _color $green
patend

Expand All @@ -120,6 +121,7 @@ patbeg index ninnerRibs
# add caps attributes
select face
attribute capsGroup !$rib+irib
ATTRIBUTE AFLR_GBC $TRANSP_UG3_GBC
attribute _color $green

# union with previous ribs/spars
Expand All @@ -132,10 +134,12 @@ patbeg index ninnerRibs
endif

patend
store internalStructure

### 3. Complete the wing geometry ###

# intersect the internal structure with the solid wing
restore internalStructure
restore solidWing
intersect

Expand Down Expand Up @@ -165,11 +169,21 @@ udprim editAttr filename <<

# add load attribute to OML
select face $capsGroup $OML
#ATTRIBUTE AFLR_GBC $TRANSP_UG3_GBC
attribute capsLoad $OML

# use this to refine certain edges
#ATTRIBUTE AFLR4_Edge_Refinement_Weight 1.0

select face
#ATTRIBUTE AFLR_GBC $TRANSP_UG3_GBC
attribute capsMesh $OML

# add AIM attribute to specify the analyses to use
select body
attribute capsAIM $egadsTessAIM;tacsAIM
ATTRIBUTE AFLR4_Cmp_ID 1
attribute capsMeshLength 1.0
attribute capsAIM $aflr4AIM;tacsAIM;egadsTessAIM

end
|||||
||||||||||||
112 changes: 112 additions & 0 deletions examples/caps_wing/unsteady_analysis_test.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,112 @@
"""
Sean Engelstad, Febuary 2023
GT SMDO Lab, Dr. Graeme Kennedy
Caps to TACS example
"""

from tacs import functions, caps2tacs
import numpy as np
from mpi4py import MPI

# -------------------------------------------------------------------------------------------------
# 1: build the tacs aim, egads aim wrapper classes
comm = MPI.COMM_WORLD
tacs_model = caps2tacs.TacsModel.build(csm_file="simple_naca_wing.csm", comm=comm)
tacs_model.egads_aim.set_mesh(
edge_pt_min=15,
edge_pt_max=20,
global_mesh_size=0.25,
max_surf_offset=0.01,
max_dihedral_angle=15,
).register_to(tacs_model)

aluminum = caps2tacs.Isotropic.aluminum().register_to(tacs_model)

# setup the thickness design variables + automatic shell properties
nribs = int(tacs_model.get_config_parameter("nribs"))
nspars = int(tacs_model.get_config_parameter("nspars"))
for irib in range(1, nribs + 1):
caps2tacs.ThicknessVariable(
caps_group=f"rib{irib}", value=0.05, material=aluminum
).register_to(tacs_model)
for ispar in range(1, nspars + 1):
caps2tacs.ThicknessVariable(
caps_group=f"spar{ispar}", value=0.05, material=aluminum
).register_to(tacs_model)
caps2tacs.ThicknessVariable(
caps_group="OML", value=0.03, material=aluminum
).register_to(tacs_model)

# add constraints and loads
caps2tacs.PinConstraint("root").register_to(tacs_model)
caps2tacs.GridForce("OML", direction=[0, 0, 1.0], magnitude=100).register_to(tacs_model)

# run the pre analysis to build tacs input files
tacs_model.setup(include_aim=True)

# ----------------------------------------------------------------------------------
# 2. Run the TACS unsteady elastic structural analysis, forward + adjoint

# create a transient problem in pyTACS
fea_solver = tacs_model.fea_solver
fea_solver.initialize()
TP = fea_solver.createTransientProblem(
"sinusoidalWing", tInit=0.0, tFinal=1.0, numSteps=10
)
timeSteps = TP.getTimeSteps()

# get the subcase of the problem and its load ID
for subcase in fea_solver.bdfInfo.subcases.values():
if "LOAD" in subcase:
loadsID = subcase["LOAD"][0]

# loop over each load
for itime, time in enumerate(timeSteps):
# compute the load scale at this time step
load_scale = np.sin(0.5 * time)

# add the loads from a static problem
TP.addLoadFromBDF(itime, loadsID, scale=load_scale)

# Add functions to the transient problem
function_names = ["mass", "ks_vmfailure"]
TP.addFunction("mass", functions.StructuralMass)
TP.addFunction("ks_vmfailure", functions.KSFailure, safetyFactor=1.5, ksWeight=50.0)

# run the tacs unsteady analysis, forward and adjoint
print("Starting unsteady tacs wing analysis, takes about 5 minutes")
tacs_funcs = {}
tacs_sens = {}
TP.solve()
TP.evalFunctions(tacs_funcs, evalFuncs=function_names)
TP.evalFunctionsSens(tacs_sens, evalFuncs=function_names)
TP.writeSolution(baseName="tacs_output", outputDir=tacs_model.analysis_dir)

struct_derivs = tacs_sens["ks_vmfailure"]["struct"]
adjoint_TD = struct_derivs[0]

# complex step derivative test
assembler = fea_solver.assembler
xvec = assembler.createDesignVec()
assembler.getDesignVars(xvec)
xarray = xvec.getArray()

# This assumes that the TACS variables are not distributed and are set
# only on the tacs_comm root processor.
h = 1e-30
if comm.rank == 0:
xarray[0] += 1j * h

assembler.setDesignVars(xvec)


tacs_funcs = {}
TP.evalFunctions(tacs_funcs)
ks_value = tacs_funcs["ks_vmfailure"]
complex_step_TD = ks_value.imag / h

relative_error = (adjoint_TD - complex_step_TD) / complex_step_TD

print(f"approximate derivative = {adjoint_TD}")
print(f"complex step derivative = {complex_step_TD}")
print(f"relative error = {relative_error}")
60 changes: 60 additions & 0 deletions tacs/caps2tacs/aflr_aim.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,60 @@
__all__ = ["AflrAim"]


class AflrAim:
def __init__(self, caps_problem, comm, root=0):
"""
MPI wrapper class for AflrAIM from ESP/CAPS
Requires use of ESPbeta at the moment
"""

self.caps_problem = caps_problem
self.comm = comm
self.root = root

# holds aflr4 aim
self._aim = None

self._build_aim()
return

@property
def root_proc(self) -> bool:
return self.comm.rank == self.root

@property
def aim(self):
"""surface mesher aim aka aflr4 aim"""
return self._aim

@property
def analysis_dir(self):
_analysis_dir = None
if self.comm.rank == self.root:
_analysis_dir = self.aim.analysisDir
_analysis_dir = self.comm.bcast(_analysis_dir, root=self.root)
return _analysis_dir

def _build_aim(self):
if self.root_proc:
self._aim = self.caps_problem.analysis.create(aim="aflr4AIM", name="aflr4")
return

def set_mesh(
self, ff_growth=1.4, min_scale=0.05, max_scale=0.5, use_quad=False, no_prox=True
):
# set surface mesh properties
if self.root_proc:
self.aim.input.ff_cdfr = ff_growth
self.aim.input.min_scale = min_scale
self.aim.input.max_scale = max_scale
self.aim.input.AFLR4_Quad = use_quad
self._aim.input.no_prox = no_prox
return self

def register_to(self, tacs_aim):
"""
cascade method to register the egads aim to the tacs aim wrapper class
"""
tacs_aim.register(self)
return self
Loading