Skip to content

Commit

Permalink
Added the actuator direction as design vars. (#525)
Browse files Browse the repository at this point in the history
  • Loading branch information
friedenhe committed Nov 22, 2023
1 parent 17d9460 commit 234fbd3
Show file tree
Hide file tree
Showing 21 changed files with 268 additions and 132 deletions.
14 changes: 7 additions & 7 deletions dafoam/mphys/mphys_dafoam.py
Original file line number Diff line number Diff line change
Expand Up @@ -653,7 +653,7 @@ def setup(self):
elif dvType == "BC": # add boundary conditions
self.add_input(dvName, distributed=False, shape_by_conn=True, tags=["mphys_coupling"])
elif dvType == "ACTD": # add actuator parameter variables
nACTDVars = 10
nACTDVars = 13
if "comps" in list(designVariables[dvName].keys()):
nACTDVars = len(designVariables[dvName]["comps"])
self.add_input(dvName, distributed=False, shape=nACTDVars, tags=["mphys_coupling"])
Expand Down Expand Up @@ -859,7 +859,7 @@ def apply_linear(self, inputs, outputs, d_inputs, d_outputs, d_residuals, mode):
# compute [dRdActD]^T*Psi using reverse mode AD
elif self.dvType[inputName] == "ACTD":
prodVec = PETSc.Vec().create(self.comm)
prodVec.setSizes((PETSc.DECIDE, 10), bsize=1)
prodVec.setSizes((PETSc.DECIDE, 13), bsize=1)
prodVec.setFromOptions()
DASolver.solverAD.calcdRdActTPsiAD(
DASolver.xvVec, DASolver.wVec, resBarVec, inputName.encode(), prodVec
Expand Down Expand Up @@ -1195,7 +1195,7 @@ def setup(self):
elif dvType == "BC": # add boundary conditions
self.add_input(dvName, distributed=False, shape_by_conn=True, tags=["mphys_coupling"])
elif dvType == "ACTD": # add actuator parameter variables
nACTDVars = 10
nACTDVars = 13
if "comps" in list(designVariables[dvName].keys()):
nACTDVars = len(designVariables[dvName]["comps"])
self.add_input(dvName, distributed=False, shape=nACTDVars, tags=["mphys_coupling"])
Expand Down Expand Up @@ -1367,7 +1367,7 @@ def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode):
# compute dFdActD
elif self.dvType[inputName] == "ACTD":
dFdACTD = PETSc.Vec().create(self.comm)
dFdACTD.setSizes((PETSc.DECIDE, 10), bsize=1)
dFdACTD.setSizes((PETSc.DECIDE, 13), bsize=1)
dFdACTD.setFromOptions()
DASolver.solverAD.calcdFdACTAD(
DASolver.xvVec, DASolver.wVec, objFuncName.encode(), inputName.encode(), dFdACTD
Expand Down Expand Up @@ -2249,15 +2249,15 @@ def setup(self):
self.fvSourceDict = self.DASolver.getOption("fvSource")

for fvSource, _ in self.aerostructDict["fvSource"].items():
self.add_input("dv_actuator_%s" % fvSource, shape=(7), distributed=False, tags=["mphys_coupling"])
self.add_input("dv_actuator_%s" % fvSource, shape=(10), distributed=False, tags=["mphys_coupling"])
self.add_input("x_prop_%s" % fvSource, shape_by_conn=True, distributed=True, tags=["mphys_coupling"])

self.add_output("actuator_%s" % fvSource, shape_by_conn=(10), distributed=False, tags=["mphys_coupling"])
self.add_output("actuator_%s" % fvSource, shape_by_conn=(13), distributed=False, tags=["mphys_coupling"])

def compute(self, inputs, outputs):
# Loop over all actuator disks
for fvSource, _ in self.aerostructDict["fvSource"].items():
actuator = np.zeros(10)
actuator = np.zeros(13)
# Update variables on root proc
if self.comm.rank == 0:
actuator[3:] = inputs["dv_actuator_%s" % fvSource][:]
Expand Down
4 changes: 2 additions & 2 deletions dafoam/pyDAFoam.py
Original file line number Diff line number Diff line change
Expand Up @@ -2227,7 +2227,7 @@ def calcTotalDerivsAOA(self, objFuncName, designVarName, dFScaling=1.0, accumula

def calcTotalDerivsACT(self, objFuncName, designVarName, designVarType, dFScaling=1.0, accumulateTotal=False):

nDVTable = {"ACTP": 9, "ACTD": 10, "ACTL": 11}
nDVTable = {"ACTP": 9, "ACTD": 13, "ACTL": 11}
nDVs = nDVTable[designVarType]

# calculate dFdACT
Expand Down Expand Up @@ -2692,7 +2692,7 @@ def solveAdjoint(self):
elif designVarDict[designVarName]["designVarType"] in ["ACTL", "ACTP", "ACTD"]:
if self.getOption("useAD")["mode"] == "fd":
designVarType = designVarDict[designVarName]["designVarType"]
nDVTable = {"ACTP": 9, "ACTD": 10, "ACTL": 11}
nDVTable = {"ACTP": 9, "ACTD": 13, "ACTL": 11}
nDVs = nDVTable[designVarType]
# calculate dRdACT
dRdACT = PETSc.Mat().create(PETSc.COMM_WORLD)
Expand Down
24 changes: 15 additions & 9 deletions src/adjoint/DAFvSource/DAFvSource.C
Original file line number Diff line number Diff line change
Expand Up @@ -143,18 +143,24 @@ void DAFvSource::syncDAOptionToActuatorDVs()
scalarList centerList;
diskSubDict.readEntry<scalarList>("center", centerList);

// we have 9 design variables for each disk
scalarList dvList(10);
scalarList dirList;
diskSubDict.readEntry<scalarList>("direction", dirList);

// we have 13 design variables for each disk
scalarList dvList(13);
dvList[0] = centerList[0];
dvList[1] = centerList[1];
dvList[2] = centerList[2];
dvList[3] = diskSubDict.getScalar("innerRadius");
dvList[4] = diskSubDict.getScalar("outerRadius");
dvList[5] = diskSubDict.getScalar("scale");
dvList[6] = diskSubDict.getScalar("POD");
dvList[7] = diskSubDict.getScalar("expM");
dvList[8] = diskSubDict.getScalar("expN");
dvList[9] = diskSubDict.getScalar("targetThrust");
dvList[3] = dirList[0];
dvList[4] = dirList[1];
dvList[5] = dirList[2];
dvList[6] = diskSubDict.getScalar("innerRadius");
dvList[7] = diskSubDict.getScalar("outerRadius");
dvList[8] = diskSubDict.getScalar("scale");
dvList[9] = diskSubDict.getScalar("POD");
dvList[10] = diskSubDict.getScalar("expM");
dvList[11] = diskSubDict.getScalar("expN");
dvList[12] = diskSubDict.getScalar("targetThrust");

// set actuatorDiskDVs_
actuatorDiskDVs_.set(diskName, dvList);
Expand Down
21 changes: 9 additions & 12 deletions src/adjoint/DAFvSource/DAFvSourceActuatorDisk.C
Original file line number Diff line number Diff line change
Expand Up @@ -222,21 +222,18 @@ void DAFvSourceActuatorDisk::calcFvSource(volVectorField& fvSource)
word diskName = fvSourceSubDict.toc()[idxI];
dictionary diskSubDict = fvSourceSubDict.subDict(diskName);

scalarList direction;
diskSubDict.readEntry<scalarList>("direction", direction);
vector dirNorm = {direction[0], direction[1], direction[2]};
vector center = {actuatorDiskDVs_[diskName][0], actuatorDiskDVs_[diskName][1], actuatorDiskDVs_[diskName][2]};
vector dirNorm = {actuatorDiskDVs_[diskName][3], actuatorDiskDVs_[diskName][4], actuatorDiskDVs_[diskName][5]};
dirNorm = dirNorm / mag(dirNorm);
vector center = {
actuatorDiskDVs_[diskName][0], actuatorDiskDVs_[diskName][1], actuatorDiskDVs_[diskName][2]};
scalar innerRadius = actuatorDiskDVs_[diskName][3];
scalar outerRadius = actuatorDiskDVs_[diskName][4];
scalar innerRadius = actuatorDiskDVs_[diskName][6];
scalar outerRadius = actuatorDiskDVs_[diskName][7];
word rotDir = diskSubDict.getWord("rotDir");
// we will calculate or read scale later
scalar scale;
scalar POD = actuatorDiskDVs_[diskName][6];
scalar POD = actuatorDiskDVs_[diskName][9];
scalar eps = diskSubDict.getScalar("eps");
scalar expM = actuatorDiskDVs_[diskName][7];
scalar expN = actuatorDiskDVs_[diskName][8];
scalar expM = actuatorDiskDVs_[diskName][10];
scalar expN = actuatorDiskDVs_[diskName][11];
// Now we need to compute normalized eps in the radial direction, i.e. epsRStar this is because
// we need to smooth the radial distribution of the thrust, here the radial location is
// normalized as rStar = (r - rInner) / (rOuter - rInner), so to make epsRStar consistent with this
Expand Down Expand Up @@ -310,12 +307,12 @@ void DAFvSourceActuatorDisk::calcFvSource(volVectorField& fvSource)
tmpThrustSumAll += fAxial * mesh_.V()[cellI];
}
reduce(tmpThrustSumAll, sumOp<scalar>());
scalar targetThrust = actuatorDiskDVs_[diskName][9];
scalar targetThrust = actuatorDiskDVs_[diskName][12];
scale = targetThrust / tmpThrustSumAll;
}
else
{
scale = actuatorDiskDVs_[diskName][5];
scale = actuatorDiskDVs_[diskName][8];
}

// now we have the correct scale, repeat the loop to assign fvSource
Expand Down
62 changes: 50 additions & 12 deletions src/adjoint/DAPartDeriv/DAPartDerivdFdACTD.C
Original file line number Diff line number Diff line change
Expand Up @@ -128,6 +128,8 @@ void DAPartDerivdFdACTD::calcPartDerivMat(

scalarList center;
diskModelSubDict.readEntry<scalarList>("center", center);
scalarList direction;
diskModelSubDict.readEntry<scalarList>("direction", direction);
scalar innerRadius = diskModelSubDict.getScalar("innerRadius");
scalar outerRadius = diskModelSubDict.getScalar("outerRadius");
scalar POD = diskModelSubDict.getScalar("POD");
Expand Down Expand Up @@ -161,42 +163,60 @@ void DAPartDerivdFdACTD::calcPartDerivMat(
diskModelSubDict.set("center", center);
}
else if (i == 3)
{
// perturb dir_x
direction[0] += delta;
diskModelSubDict.set("direction", direction);
}
else if (i == 4)
{
// perturb dir_y
direction[1] += delta;
diskModelSubDict.set("direction", direction);
}
else if (i == 5)
{
// perturb dir_z
direction[2] += delta;
diskModelSubDict.set("direction", direction);
}
else if (i == 6)
{
// perturb innerRadius
innerRadius += delta;
diskModelSubDict.set("innerRadius", innerRadius);
}
else if (i == 4)
else if (i == 7)
{
// perturb outerRadius
outerRadius += delta;
diskModelSubDict.set("outerRadius", outerRadius);
}
else if (i == 5)
else if (i == 8)
{
// perturb scale
scale += delta;
diskModelSubDict.set("scale", scale);
}
else if (i == 6)
else if (i == 9)
{
// perturb POD
POD += delta;
diskModelSubDict.set("POD", POD);
}
else if (i == 7)
else if (i == 10)
{
// perturb expM
expM += delta;
diskModelSubDict.set("expM", expM);
}
else if (i == 8)
else if (i == 11)
{
// perturb expN
expN += delta;
diskModelSubDict.set("expN", expN);
}
else if (i == 9)
else if (i == 12)
{
targetThrust += delta;
diskModelSubDict.set("targetThrust", targetThrust);
Expand Down Expand Up @@ -229,42 +249,60 @@ void DAPartDerivdFdACTD::calcPartDerivMat(
diskModelSubDict.set("center", center);
}
else if (i == 3)
{
// reset dir_x
direction[0] -= delta;
diskModelSubDict.set("direction", direction);
}
else if (i == 4)
{
// reset dir_y
direction[1] -= delta;
diskModelSubDict.set("direction", direction);
}
else if (i == 5)
{
// reset dir_z
direction[2] -= delta;
diskModelSubDict.set("direction", direction);
}
else if (i == 6)
{
// reset innerRadius
innerRadius -= delta;
diskModelSubDict.set("innerRadius", innerRadius);
}
else if (i == 4)
else if (i == 7)
{
// reset outerRadius
outerRadius -= delta;
diskModelSubDict.set("outerRadius", outerRadius);
}
else if (i == 5)
else if (i == 8)
{
// reset scale
scale -= delta;
diskModelSubDict.set("scale", scale);
}
else if (i == 6)
else if (i == 9)
{
// reset POD
POD -= delta;
diskModelSubDict.set("POD", POD);
}
else if (i == 7)
else if (i == 10)
{
// reset expM
expM -= delta;
diskModelSubDict.set("expM", expM);
}
else if (i == 8)
else if (i == 11)
{
// reset expN
expN -= delta;
diskModelSubDict.set("expN", expN);
}
else if (i == 9)
else if (i == 12)
{
targetThrust -= delta;
diskModelSubDict.set("targetThrust", targetThrust);
Expand Down
2 changes: 1 addition & 1 deletion src/adjoint/DAPartDeriv/DAPartDerivdFdACTD.H
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@ class DAPartDerivdFdACTD
protected:

/// number of design variables for actuator disk
label nActDVs_ = 10;
label nActDVs_ = 13;

public:
TypeName("dFdACTD");
Expand Down
Loading

0 comments on commit 234fbd3

Please sign in to comment.