Skip to content

Commit

Permalink
Added more solvers and functions. (#734)
Browse files Browse the repository at this point in the history
* Removed the distributed arg for calcJac.

* Added DARhoSimpleCFoam.

* Added DAHeatTransferFoam and related functions.

* Added the missing files.

* Fixed the DATurb order.

* Used mpicc to auto-detect include and lib paths.

* Added DASimpleTFoam and related functions.

* Added the missing files.

* Added kOmegaSST and moved Coloring into DASolver.

* Fixed the code cov.

* Re-organized new turbModel structure and added Fv3 and laminar models.

* Added the moment function.

* Fixed an issue to avoid computing all dF partials.

* Removed un-used funcs in Python layer.

* Removed unused funcs in DASolver.

* Removed parts for DAFunction and re-organized codes.
  • Loading branch information
friedenhe authored Jan 3, 2025
1 parent d5680e8 commit e09e740
Show file tree
Hide file tree
Showing 117 changed files with 2,359 additions and 6,916 deletions.
37 changes: 32 additions & 5 deletions Allclean
Original file line number Diff line number Diff line change
Expand Up @@ -5,8 +5,35 @@ if [ -z "$WM_PROJECT" ]; then
exit 1
fi

cd src/adjoint && rm -rf lnInclude && cd -
cd src/adjoint && ./Allclean && cd -
cd src/pyUnitTests && ./Allclean && cd -
cd src/pyDASolvers && ./Allclean && cd -
cd src/utilities/coloring && ./Allclean && cd -
function cleanDAFoam()
{
cd src/adjoint && rm -rf lnInclude && cd -
cd src/adjoint && ./Allclean && cd -
cd src/newTurbModels && ./Allclean && cd -
cd src/pyUnitTests && ./Allclean && cd -
cd src/pyDASolvers && ./Allclean && cd -
}

# clean original mode
echo "***************** Cleaning original mode **************"
. $DAFOAM_ROOT_PATH/loadDAFoam.sh
cleanDAFoam
# clean ADR mode
echo "***************** Cleaning ADR mode **************"
. $DAFOAM_ROOT_PATH/loadDAFoam.sh
. $DAFOAM_ROOT_PATH/OpenFOAM/OpenFOAM-v1812-ADR/etc/bashrc
cleanDAFoam
# if COMPILE_DAFOAM_ADF is set, compile ADF mode
if [ -z "$COMPILE_DAFOAM_ADF" ]; then
echo "COMPILE_DAFOAM_ADF is not set. skip the ADF mode"
else
echo "***************** Cleaning ADF mode **************"
. $DAFOAM_ROOT_PATH/loadDAFoam.sh
. $DAFOAM_ROOT_PATH/OpenFOAM/OpenFOAM-v1812-ADF/etc/bashrc
cleanDAFoam
fi

# reset the OpenFOAM environment to the original mode
. $DAFOAM_ROOT_PATH/loadDAFoam.sh


6 changes: 4 additions & 2 deletions Allmake
Original file line number Diff line number Diff line change
Expand Up @@ -11,9 +11,9 @@ function makeDAFoam()
set -e
wmakeLnInclude src/adjoint
cd src/adjoint && ./Allmake && cd -
cd src/newTurbModels && ./Allmake && cd -
cd src/pyUnitTests && ./Allmake && cd -
cd src/pyDASolvers && ./Allmake && cd -
cd src/utilities/coloring && ./Allmake && cd -
# disable the -e flag
set +e
}
Expand All @@ -40,4 +40,6 @@ fi
# reset the OpenFOAM environment to the original mode
. $DAFOAM_ROOT_PATH/loadDAFoam.sh

pip install .
pip install .

ls -R dafoam/libs && echo " " && echo "*** Build Successful! ***" && echo " "
113 changes: 25 additions & 88 deletions dafoam/mphys/mphys_dafoam.py
Original file line number Diff line number Diff line change
Expand Up @@ -606,10 +606,6 @@ def setup(self):
else:
self.runColoring = True

# determine which function to compute the adjoint
self.evalFuncs = []
DASolver.setEvalFuncs(self.evalFuncs)

# setup input and output for the solver
# we need to add states as outputs for all cases
local_state_size = DASolver.getNLocalAdjointStates()
Expand Down Expand Up @@ -667,17 +663,16 @@ def calcFFD2XvSeeds(self):
return xVDot

# calculate the residual
def apply_nonlinear(self, inputs, outputs, residuals):
DASolver = self.DASolver
# NOTE: we do not pass the states from inputs to the OF layer.
# this can cause potential convergence issue because the initial states
# in the inputs are set to all ones. So passing this all-ones states
# into the OF layer may diverge the primal solver. Here we can always
# use the states from the OF layer to compute the residuals.
# DASolver.setStates(outputs["%s_states" % self.discipline])

# get flow residuals from DASolver
residuals[self.stateName] = DASolver.getResiduals()
# def apply_nonlinear(self, inputs, outputs, residuals):
# DASolver = self.DASolver
# # NOTE: we do not pass the states from inputs to the OF layer.
# # this can cause potential convergence issue because the initial states
# # in the inputs are set to all ones. So passing this all-ones states
# # into the OF layer may diverge the primal solver. Here we can always
# # use the states from the OF layer to compute the residuals.
# # DASolver.setStates(outputs["%s_states" % self.discipline])
# # get flow residuals from DASolver
# residuals[self.stateName] = DASolver.getResiduals()

def set_solver_input(self, inputs):
"""
Expand Down Expand Up @@ -726,7 +721,7 @@ def solve_nonlinear(self, inputs, outputs):

# get the objective functions
funcs = {}
DASolver.evalFunctions(funcs, evalFuncs=self.evalFuncs)
DASolver.evalFunctions(funcs)

# assign the computed flow states to outputs
states = DASolver.getStates()
Expand All @@ -740,6 +735,10 @@ def solve_nonlinear(self, inputs, outputs):
# set states
DASolver.setStates(states)

# We also need to just calculate the residual for the AD mode to initialize vars like URes
# We do not print the residual for AD, though
DASolver.solverAD.calcPrimalResidualStatistics("calc".encode())

def linearize(self, inputs, outputs, residuals):
# NOTE: we do not do any computation in this function, just print some information

Expand Down Expand Up @@ -767,7 +766,6 @@ def apply_linear(self, inputs, outputs, d_inputs, d_outputs, d_residuals, mode):
DASolver.setStates(outputs[self.stateName])

localAdjSize = DASolver.getNLocalAdjointStates()
localXvSize = DASolver.getNLocalPoints() * 3

if self.stateName in d_residuals:

Expand All @@ -782,12 +780,10 @@ def apply_linear(self, inputs, outputs, d_inputs, d_outputs, d_residuals, mode):
self.stateName,
"stateVar",
localAdjSize,
1,
jacInput,
self.residualName,
"residual",
localAdjSize,
1,
seed,
product,
)
Expand All @@ -798,62 +794,21 @@ def apply_linear(self, inputs, outputs, d_inputs, d_outputs, d_residuals, mode):
for inputName in list(inputDict.keys()):
inputType = inputDict[inputName]["type"]
inputSize = DASolver.solver.getInputSize(inputName, inputType)
inputDistributed = DASolver.solver.getInputDistributed(inputName, inputType)
product = np.zeros(inputSize)
jacInput = inputs[inputName]
DASolver.solverAD.calcJacTVecProduct(
inputName,
inputType,
inputSize,
inputDistributed,
jacInput,
self.residualName,
"residual",
localAdjSize,
1,
seed,
product,
)
d_inputs[inputName] += product

## this computes [dRdXv]^T*Psi using reverse mode AD
# if inputType == self.volCoordName:
# product = np.zeros(localXvSize)
# jacInput = inputs[self.volCoordName]
# DASolver.solverAD.calcJacTVecProduct(
# self.volCoordName,
# "volCoord",
# localXvSize,
# 1,
# jacInput,
# self.residualName,
# "residual",
# localAdjSize,
# 1,
# seed,
# product,
# )
# d_inputs[self.volCoordName] += product
# elif inputType == "patchVelocity":
# product = np.zeros(2)
# jacInput = inputs[inputName]
# DASolver.solverAD.calcJacTVecProduct(
# inputName,
# "patchVelocity",
# 2,
# 0,
# jacInput,
# self.residualName,
# "residual",
# localAdjSize,
# 1,
# seed,
# product,
# )
# d_inputs[inputName] += product
# else:
# raise AnalysisError("inputType %s not supported! " % inputType)

def solve_linear(self, d_outputs, d_residuals, mode):
# solve the adjoint equation [dRdW]^T * Psi = dFdW

Expand All @@ -880,7 +835,7 @@ def solve_linear(self, d_outputs, d_residuals, mode):

# run coloring
if self.DASolver.getOption("adjUseColoring") and self.runColoring:
self.DASolver.runColoring()
self.DASolver.solver.runColoring()
self.runColoring = False

if adjEqnSolMethod == "Krylov":
Expand Down Expand Up @@ -1149,20 +1104,11 @@ def setup(self):
inputDistributed = DASolver.solver.getInputDistributed(inputName, inputType)
self.add_input(inputName, distributed=inputDistributed, shape=inputSize, tags=["mphys_coupling"])

def mphys_add_funcs(self):
# add the function names to this component, called from runScript.py

# it is called function in DAOptions but it contains both objective and constraint functions
# add outputs
functions = self.DASolver.getOption("function")

self.funcs = []

for function in functions:
self.funcs.append(function)

# loop over the functions here and create the output
for f_name in self.funcs:
self.add_output(f_name, distributed=False, shape=1, units=None, tags=["mphys_result"])
for f_name in list(functions.keys()):
self.add_output(f_name, distributed=False, shape=1, units=None)

# get the objective function from DASolver
def compute(self, inputs, outputs):
Expand All @@ -1172,12 +1118,9 @@ def compute(self, inputs, outputs):
DASolver.setStates(inputs[self.stateName])

funcs = {}

if self.funcs is not None:
DASolver.evalFunctions(funcs, evalFuncs=self.funcs)
for f_name in self.funcs:
if f_name in funcs:
outputs[f_name] = funcs[f_name]
DASolver.evalFunctions(funcs)
for f_name in list(outputs.keys()):
outputs[f_name] = funcs[f_name]

# compute the partial derivatives of functions
def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode):
Expand Down Expand Up @@ -1207,7 +1150,7 @@ def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode):

# if the seed is zero, do not compute
if abs(seed) < 1e-14:
pass
continue

for inputName in list(d_inputs.keys()):
# compute dFdW * seed
Expand All @@ -1218,32 +1161,27 @@ def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode):
self.stateName,
"stateVar",
localAdjSize,
1,
jacInput,
functionName,
"function",
1,
0,
seed,
product,
)
d_inputs[self.stateName] += product
else:
inputType = inputDict[inputName]["type"]
inputSize = DASolver.solver.getInputSize(inputName, inputType)
inputDistributed = DASolver.solver.getInputDistributed(inputName, inputType)
product = np.zeros(inputSize)
jacInput = inputs[inputName]
DASolver.solverAD.calcJacTVecProduct(
inputName,
inputType,
inputSize,
inputDistributed,
jacInput,
functionName,
"function",
1,
0,
seed,
product,
)
Expand Down Expand Up @@ -1282,8 +1220,6 @@ def compute(self, inputs, outputs):
DASolver.setSurfaceCoordinates(x_a, DASolver.designSurfacesGroup)
DASolver.mesh.warpMesh()
solverGrid = DASolver.mesh.getSolverGrid()
# actually change the mesh in the C++ layer by setting xvVec
DASolver.xvFlatten2XvVec(solverGrid, DASolver.xvVec)
outputs["%s_vol_coords" % self.discipline] = solverGrid

# compute the mesh warping products in IDWarp
Expand Down Expand Up @@ -2111,7 +2047,8 @@ def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode):

class DAFoamActuator(ExplicitComponent):
"""
Component that updates actuator disk definition variables when actuator disks are displaced in an aerostructural case.
Component that updates actuator disk definition variables when actuator disks
are displaced in an aerostructural case.
"""

def initialize(self):
Expand Down
Loading

0 comments on commit e09e740

Please sign in to comment.