From 21dc6b03ce30f41578733e19d3bbbde06d3b7634 Mon Sep 17 00:00:00 2001 From: Ping He Date: Mon, 29 Jul 2024 15:39:18 -0500 Subject: [PATCH] Fixed a potential issues for all-zero dF seeds (#667) * Fixed a potential issue for zero fBar. * Small change to the dF structure. --- dafoam/mphys/mphys_dafoam.py | 21 +++++++++++++++++---- 1 file changed, 17 insertions(+), 4 deletions(-) diff --git a/dafoam/mphys/mphys_dafoam.py b/dafoam/mphys/mphys_dafoam.py index cd204e3a..ddbcaf6d 100644 --- a/dafoam/mphys/mphys_dafoam.py +++ b/dafoam/mphys/mphys_dafoam.py @@ -893,7 +893,7 @@ def apply_linear(self, inputs, outputs, d_inputs, d_outputs, d_residuals, mode): d_inputs[inputName] += ACTDBarSub else: d_inputs[inputName] += ACTDBar - + # compute [dRdHSC]^T*Psi using reverse mode AD elif self.dvType[inputName] == "HSC": prodVec = PETSc.Vec().create(self.comm) @@ -1351,6 +1351,18 @@ def compute(self, inputs, outputs): # compute the partial derivatives of functions def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode): + # we first check if all the seeds in d_outputs are zeros. If yes, we return without calculation anything + n_non_zero_seeds = 0 + for func_name in d_outputs: + if d_outputs[func_name] != 0.0: + n_non_zero_seeds += 1 + if n_non_zero_seeds == 0: + return + elif n_non_zero_seeds > 1: + if self.comm.rank == 0: + print("************* Warning *************") + print("More than one non-zero seed found! ", d_outputs) + DASolver = self.DASolver # set the runStatus, this is useful when the actuator term is activated @@ -1381,7 +1393,8 @@ def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode): funcsBar = {} # assign value to funcsBar. NOTE: we only assign seed if d_outputs has - # non-zero values! + # non-zero values! We assume OM will pass only one non-zero seed here + # therefore, funcsBar should have only one key if self.funcs is None: raise AnalysisError("functions not set! Forgot to call mphys_add_funcs?") else: @@ -1394,7 +1407,7 @@ def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode): if self.comm.rank == 0: print("Computing partials for ", list(funcsBar.keys())) - + # update the obj func name for solve_linear later DASolver.setOption("solveLinearObjFuncName", list(funcsBar.keys())[0]) DASolver.updateDAOption() @@ -1478,7 +1491,7 @@ def compute_jacvec_product(self, inputs, d_inputs, d_outputs, mode): d_inputs[inputName] += ACTDBarSub * fBar else: d_inputs[inputName] += ACTDBar * fBar - + # compute dFdHSC elif self.dvType[inputName] == "HSC": dFdHSC = PETSc.Vec().create(self.comm)