From 71e1d84e5d436223f3d8d361a215da071b1a1133 Mon Sep 17 00:00:00 2001 From: sseraj Date: Fri, 4 Feb 2022 16:09:52 -0500 Subject: [PATCH 01/12] fixed cutCallBack signature for full update --- adflow/pyADflow.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/adflow/pyADflow.py b/adflow/pyADflow.py index 78a5bc9f5..c0fea49b9 100644 --- a/adflow/pyADflow.py +++ b/adflow/pyADflow.py @@ -3762,7 +3762,8 @@ def updateGeometryInfo(self, warpMesh=True): cutCallBack = self.getOption("cutCallBack") if cutCallBack is not None: xCen = self.adflow.utils.getcellcenters(1, n).T - cutCallBack(xCen, flag) + cellIDs = self.adflow.utils.getcellcgnsblockids(1, n) + cutCallBack(xCen, self.CGNSZoneNameIDs, cellIDs, flag) # Verify previous mesh failures self.adflow.killsignals.routinefailed = self.comm.allreduce( From da7fa5e8b409be002451d3f0c8225430bb59cc9e Mon Sep 17 00:00:00 2001 From: sseraj Date: Fri, 25 Feb 2022 11:16:18 -0500 Subject: [PATCH 02/12] allocated xCen for fast overset update in parallel --- src/overset/finalOversetCommStructures.F90 | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/overset/finalOversetCommStructures.F90 b/src/overset/finalOversetCommStructures.F90 index ac621d922..9b6b03af1 100644 --- a/src/overset/finalOversetCommStructures.F90 +++ b/src/overset/finalOversetCommStructures.F90 @@ -303,7 +303,8 @@ subroutine finalOversetCommStructures(level, sps) commPatternOverset(level, sps)%sendList(iSendProc)%block(n), & commPatternOverset(level, sps)%sendList(iSendProc)%indices(n, 3), & commPatternOverset(level, sps)%sendList(iSendProc)%interp(n, 8), & - commPatternOverset(level, sps)%sendList(iSendProc)%interpd(n, 8)) + commPatternOverset(level, sps)%sendList(iSendProc)%interpd(n, 8), & + commPatternOverset(level, sps)%sendList(iSendProc)%xCen(n, 3)) ! Now set the data do i=1, n @@ -323,6 +324,7 @@ subroutine finalOversetCommStructures(level, sps) call fracToWeights(realRecvBuf(jj+1:jj+3), & commPatternOverset(level, sps)%sendList(iSendProc)%interp(i, :)) commPatternOverset(level, sps)%sendList(iSendProc)%interpd(i, :) = zero + commPatternOverset(level, sps)%sendList(iSendProc)%xCen(i, :) = zero jj = jj + 3 end do end if From 32c2ce748af831b58ad9051be3cff937acc8e085 Mon Sep 17 00:00:00 2001 From: sseraj Date: Fri, 25 Feb 2022 11:53:04 -0500 Subject: [PATCH 03/12] re-enable fast overset update in AD routines --- src/adjoint/masterRoutines.F90 | 26 +++++++++++++------------- 1 file changed, 13 insertions(+), 13 deletions(-) diff --git a/src/adjoint/masterRoutines.F90 b/src/adjoint/masterRoutines.F90 index 639ca867f..8c1063585 100644 --- a/src/adjoint/masterRoutines.F90 +++ b/src/adjoint/masterRoutines.F90 @@ -83,9 +83,7 @@ subroutine master(useSpatial, famLists, funcValues, forces, & do sps=1, nTimeIntervalsSpectral ! Update overset connectivity if necessary - if (oversetPresent .and. & - (oversetUpdateMode == updateFast .or. & - oversetUpdateMode == updateFull)) then + if (oversetPresent .and. oversetUpdateMode == updateFast) then call updateOversetConnectivity(1_intType, sps) end if end do @@ -349,11 +347,12 @@ subroutine master_d(wdot, xdot, forcesDot, dwDot, famLists, funcValues, funcValu do sps=1, nTimeIntervalsSpectral ! Update overset connectivity if necessary - if (oversetPresent .and. & - (oversetUpdateMode == updateFast .or. & - oversetUpdateMode == updateFull)) then - print *,'Full overset update derivatives not implemented' - !call updateOversetConnectivity_d(1_intType, sps) + if (oversetPresent) then + if (oversetUpdateMode == updateFast) then + call updateOversetConnectivity_d(1_intType, sps) + else if (oversetUpdateMode == updateFull) then + print *,'Full overset update derivatives not implemented' + end if end if end do @@ -837,11 +836,12 @@ subroutine master_b(wbar, xbar, extraBar, forcesBar, dwBar, nState, famLists, & do sps=1, nTimeIntervalsSpectral ! Update overset connectivity if necessary - if (oversetPresent .and. & - (oversetUpdateMode == updateFast .or. & - oversetUpdateMode == updateFull)) then - print *,'Full overset update derivatives not implemented' - !call updateOversetConnectivity_b(1_intType, sps) + if (oversetPresent) then + if (oversetUpdateMode == updateFast) then + call updateOversetConnectivity_b(1_intType, sps) + else if (oversetUpdateMode == updateFull) then + print *,'Full overset update derivatives not implemented' + end if end if end do ! Now the adjoint of the coordinate exhcange From 54f6e24d8b98b23ec0c4bff39f96d1a6f0e94cdb Mon Sep 17 00:00:00 2001 From: sseraj Date: Fri, 25 Feb 2022 12:20:35 -0500 Subject: [PATCH 04/12] print on root only for full update --- src/adjoint/masterRoutines.F90 | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/src/adjoint/masterRoutines.F90 b/src/adjoint/masterRoutines.F90 index 8c1063585..2e68a668e 100644 --- a/src/adjoint/masterRoutines.F90 +++ b/src/adjoint/masterRoutines.F90 @@ -239,7 +239,7 @@ subroutine master_d(wdot, xdot, forcesDot, dwDot, famLists, funcValues, funcValu use constants use diffsizes, only : ISIZE1OFDrfbcdata, ISIZE1OFDrfviscsubface - use communication, only : adflow_comm_world + use communication, only : adflow_comm_world, myID use iteration, only : currentLevel use BCExtra_d, only : applyAllBC_Block_d use inputAdjoint, only : viscPC @@ -351,7 +351,9 @@ subroutine master_d(wdot, xdot, forcesDot, dwDot, famLists, funcValues, funcValu if (oversetUpdateMode == updateFast) then call updateOversetConnectivity_d(1_intType, sps) else if (oversetUpdateMode == updateFull) then - print *,'Full overset update derivatives not implemented' + if (myID == 0) then + print *,'Full overset update derivatives not implemented' + end if end if end if end do @@ -840,7 +842,9 @@ subroutine master_b(wbar, xbar, extraBar, forcesBar, dwBar, nState, famLists, & if (oversetUpdateMode == updateFast) then call updateOversetConnectivity_b(1_intType, sps) else if (oversetUpdateMode == updateFull) then - print *,'Full overset update derivatives not implemented' + if (myID == 0) then + print *,'Full overset update derivatives not implemented' + end if end if end if end do From 77509bf9d0ccf0e51a99c831dee6df1b4d1ee761 Mon Sep 17 00:00:00 2001 From: sseraj Date: Wed, 2 Mar 2022 11:50:39 -0500 Subject: [PATCH 05/12] update connectivity for fast mode only in blocketteRes --- src/NKSolver/blockette.F90 | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/NKSolver/blockette.F90 b/src/NKSolver/blockette.F90 index 25acea4e5..67c6af77b 100644 --- a/src/NKSolver/blockette.F90 +++ b/src/NKSolver/blockette.F90 @@ -189,9 +189,7 @@ subroutine blocketteRes(useDissApprox, useViscApprox, useUpdateIntermed, useFlow do sps=1, nTimeIntervalsSpectral ! Update overset connectivity if necessary - if (oversetPresent .and. & - (oversetUpdateMode == updateFast .or. & - oversetUpdateMode == updateFull)) then + if (oversetPresent .and. oversetUpdateMode == updateFast) then call updateOversetConnectivity(1_intType, sps) end if end do From 7ae6407cb0f199a34c11b607c1df5185c03315d0 Mon Sep 17 00:00:00 2001 From: sseraj Date: Mon, 27 Jun 2022 18:23:36 -0400 Subject: [PATCH 06/12] added overset update tests --- tests/reg_tests/test_zipper.py | 52 +++++++++++++++++++++++++++++++--- 1 file changed, 48 insertions(+), 4 deletions(-) diff --git a/tests/reg_tests/test_zipper.py b/tests/reg_tests/test_zipper.py index 8d7704898..705c31810 100644 --- a/tests/reg_tests/test_zipper.py +++ b/tests/reg_tests/test_zipper.py @@ -1,14 +1,16 @@ import unittest import os import copy +from parameterized import parameterized_class from adflow import ADFLOW -from reg_default_options import adflowDefOpts +from reg_default_options import adflowDefOpts, IDWarpDefOpts from reg_aeroproblems import ap_simple_cart_cube - +from idwarp import USMesh +from pygeo import DVGeometry baseDir = os.path.dirname(os.path.abspath(__file__)) -# Tests for zipper meshes using a mesh with two overlapping cubes +# Tests for overset and zipper meshes using a mesh with two overlapping cubes gridFile = os.path.join(baseDir, "../../input_files/cube_overset.cgns") commonTestOptions = { @@ -25,7 +27,29 @@ "nearWallDist": 1.0, } - +test_params = [ + { + "name": "frozen", + "options": { + "oversetUpdateMode": "frozen", + }, + }, + { + "name": "fast", + "options": { + "oversetUpdateMode": "fast", + }, + }, + { + "name": "full", + "options": { + "oversetUpdateMode": "full", + }, + }, +] + + +@parameterized_class(test_params) class TestCubeOverset(unittest.TestCase): N_PROCS = 2 @@ -42,12 +66,32 @@ def setUp(self): # These are the modified options common to these tests options.update(commonTestOptions) + # Add the parameterized options + options.update(self.options) + # Define the AeroProblem self.ap = copy.deepcopy(ap_simple_cart_cube) # Create the solver self.CFDSolver = ADFLOW(options=options, debug=False) + # Create mesh warping object + meshOptions = copy.copy(IDWarpDefOpts) + meshOptions.update({"gridFile": options["gridfile"]}) + self.CFDSolver.setMesh(USMesh(options=meshOptions)) + + # Create geometry object + ffdFile = os.path.join(baseDir, "../../input_files/cube_overset_ffd.xyz") + DVGeo = DVGeometry(ffdFile) + DVGeo.addLocalDV("shape", lower=-0.5, upper=0.5, axis="z") + self.CFDSolver.setDVGeo(DVGeo) + + # Apply shape change + dvDict = self.CFDSolver.DVGeo.getValues() + dvDict["shape"][4] = 0.1 + self.CFDSolver.setAeroProblem(self.ap) + self.CFDSolver.DVGeo.setDesignVars(dvDict) + def test_convergence(self): # Solve the flow From 36e4c7f8cb368293cf500b7115bb8753695451fb Mon Sep 17 00:00:00 2001 From: sseraj Date: Mon, 27 Jun 2022 18:24:08 -0400 Subject: [PATCH 07/12] renamed test_zipper --- tests/reg_tests/{test_zipper.py => test_overset.py} | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename tests/reg_tests/{test_zipper.py => test_overset.py} (100%) diff --git a/tests/reg_tests/test_zipper.py b/tests/reg_tests/test_overset.py similarity index 100% rename from tests/reg_tests/test_zipper.py rename to tests/reg_tests/test_overset.py From bf20ce30d9c65471a97868fb17486784d52db4f9 Mon Sep 17 00:00:00 2001 From: sseraj Date: Mon, 27 Jun 2022 18:35:14 -0400 Subject: [PATCH 08/12] updated docs --- doc/options.yaml | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/doc/options.yaml b/doc/options.yaml index 8a3481764..8e967937f 100644 --- a/doc/options.yaml +++ b/doc/options.yaml @@ -491,11 +491,14 @@ oversetUpdateMode: This is the cheapest option to use when the entire overset mesh is warped together, such as with USMesh in IDWarp. fast: > Only the weights are updated and the donors are unchanged. - This is usually a fairly small correction and is fast to run. This is used when the entire overset mesh is warped together, such as with USMesh in IDWarp. + This is fast to run but only works for small design changes. + Large design changes will cause the solver to stall because of nonphysical interpolation weights. full: > The overset connectivity is recomputed from stratch. This is usually only used when the component meshes are warped independently, such as with MultiUSMesh in IDWarp. + The hole cutting is likely to fail if there are large design changes. + Even if the hole cutting works, the aerodynamic derivatives will be inaccurate because the hole cutting routines are not differentiated. nRefine: desc: > From 48471a008d7bc938ee3b6582b2568077ee4aacb7 Mon Sep 17 00:00:00 2001 From: sseraj Date: Mon, 27 Jun 2022 18:41:35 -0400 Subject: [PATCH 09/12] add cutCallback to test --- tests/reg_tests/test_overset.py | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/tests/reg_tests/test_overset.py b/tests/reg_tests/test_overset.py index 705c31810..7419e11d1 100644 --- a/tests/reg_tests/test_overset.py +++ b/tests/reg_tests/test_overset.py @@ -12,6 +12,13 @@ # Tests for overset and zipper meshes using a mesh with two overlapping cubes +# This does not affect the hole cutting for this case but tests that cutCallback is working +def cutCallback(xCen, CGNSZoneNameIDs, cellIDs, flags): + + # Blank cells in the negative y-axis + flags[xCen[:, 1] < 0] = 1 + + gridFile = os.path.join(baseDir, "../../input_files/cube_overset.cgns") commonTestOptions = { "gridFile": gridFile, @@ -25,6 +32,7 @@ "useANKSolver": True, "useNKSolver": True, "nearWallDist": 1.0, + "cutCallback": cutCallback, } test_params = [ From 959be77520dae3cb2558cd23e54618fd55c0c273 Mon Sep 17 00:00:00 2001 From: sseraj Date: Tue, 19 Jul 2022 15:00:27 -0400 Subject: [PATCH 10/12] added fast update fraction check --- doc/options.yaml | 3 ++- src/overset/oversetAPI.F90 | 8 ++++---- src/overset/oversetCommUtilites.F90 | 17 ++++++++++++++++- 3 files changed, 22 insertions(+), 6 deletions(-) diff --git a/doc/options.yaml b/doc/options.yaml index 8e967937f..e025ccb85 100644 --- a/doc/options.yaml +++ b/doc/options.yaml @@ -493,7 +493,8 @@ oversetUpdateMode: Only the weights are updated and the donors are unchanged. This is used when the entire overset mesh is warped together, such as with USMesh in IDWarp. This is fast to run but only works for small design changes. - Large design changes will cause the solver to stall because of nonphysical interpolation weights. + Large design changes can shift the position of interpolate cells relative to their donors enough to cause nonphysical interpolation weights. + In this case, an error is raised. full: > The overset connectivity is recomputed from stratch. This is usually only used when the component meshes are warped independently, such as with MultiUSMesh in IDWarp. diff --git a/src/overset/oversetAPI.F90 b/src/overset/oversetAPI.F90 index fe667c1cf..485905a7b 100644 --- a/src/overset/oversetAPI.F90 +++ b/src/overset/oversetAPI.F90 @@ -2101,16 +2101,16 @@ end subroutine setExplicitHoleCut subroutine updateOverset(flag, n, closedFamList, nFam) ! This is the main gateway routine for updating the overset - ! connecitivty during a solution. It is wrapped and intended to be + ! connectivity during a solution. It is wrapped and intended to be ! called from Python. What this routine does depends on the value ! of oversetUpdateMode: ! updateFrozen: Nothing happens. The initial - ! connecitivty computed during initialzation is kept. + ! connectivity computed during initialzation is kept. - ! updteFast: Update just the weight, but leave the donors + ! updateFast: Update just the weight, but leave the donors ! unchanged. This is only applicable when the entire mesh is - ! warped at the same time with pyWarpuStruct. + ! warped at the same time like with USMesh in IDWarp. ! updateFull: Complete from scratch update. Run the full ! oversetComm routine. diff --git a/src/overset/oversetCommUtilites.F90 b/src/overset/oversetCommUtilites.F90 index ae565166a..2ca0cc882 100644 --- a/src/overset/oversetCommUtilites.F90 +++ b/src/overset/oversetCommUtilites.F90 @@ -1740,7 +1740,7 @@ subroutine updateOversetConnectivity(level, sps) ! mesh. It does *not* completely redo the connectivity. Rather, a ! newton search on the existing donors are performed using the ! updated coordinates. This type of update is only applicable if the - ! entire volume mesh is warped as one a-la pyWarpUstruct. This + ! entire volume mesh is warped as one like with USMesh in IDWarp. This ! actually ends up being a fairly small correction most of the time, ! however, desipite looks to the contrary is actually quite fast to ! run. @@ -1767,6 +1767,9 @@ subroutine updateOversetConnectivity(level, sps) real(kind=realType) :: frac(3), frac0(3), xCen(3) integer(kind=intType), dimension(8), parameter :: indices=(/1,2,4,3,5,6,8,7/) + ! Set a tolerance for checking whether fractions are between 0 and 1 + real(kind=realType) :: fracTol=1e-8 + ! Pointers to the overset comms to make it easier to read commPattern => commPatternOverset(level, sps) internal => internalOverset(level, sps) @@ -1913,6 +1916,12 @@ subroutine updateOversetConnectivity(level, sps) call newtonUpdate(xCen, & flowDoms(d1, level, sps)%x(i1-1:i1+1, j1-1:j1+1, k1-1:k1+1, :), frac0, frac) + ! Check if the fractions are between 0 and 1 + if (ANY(frac > one + fracTol) .or. ANY(frac < zero - fracTol)) then + print *, "Invalid overset connectivity update. Use 'frozen' or 'full' oversetUpdateMode instead." + error stop + end if + ! Set the new weights call fracToWeights(frac, internal%donorInterp(i, :)) enddo localInterp @@ -1951,6 +1960,12 @@ subroutine updateOversetConnectivity(level, sps) call newtonUpdate(xCen, & flowDoms(d2, level, sps)%x(i2-1:i2+1, j2-1:j2+1, k2-1:k2+1, :), frac0, frac) + ! Check if the fractions are between zero and one + if (ANY(frac > one + fracTol) .or. ANY(frac < zero - fracTol)) then + print *, "Invalid overset connectivity update. Use 'frozen' or 'full' oversetUpdateMode instead." + error stop + end if + ! Set the new weights call fracToWeights(frac, commPattern%sendList(ii)%interp(j, :)) enddo From 56a2abda1c733ec68155880cf1009e8f5de13fba Mon Sep 17 00:00:00 2001 From: sseraj Date: Tue, 19 Jul 2022 16:03:41 -0400 Subject: [PATCH 11/12] modified fraction check to work with complexify --- src/overset/oversetAPI.F90 | 2 +- src/overset/oversetCommUtilites.F90 | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/overset/oversetAPI.F90 b/src/overset/oversetAPI.F90 index 485905a7b..952838238 100644 --- a/src/overset/oversetAPI.F90 +++ b/src/overset/oversetAPI.F90 @@ -2106,7 +2106,7 @@ subroutine updateOverset(flag, n, closedFamList, nFam) ! of oversetUpdateMode: ! updateFrozen: Nothing happens. The initial - ! connectivity computed during initialzation is kept. + ! connectivity computed during initialization is kept. ! updateFast: Update just the weight, but leave the donors ! unchanged. This is only applicable when the entire mesh is diff --git a/src/overset/oversetCommUtilites.F90 b/src/overset/oversetCommUtilites.F90 index 2ca0cc882..6422d5cfd 100644 --- a/src/overset/oversetCommUtilites.F90 +++ b/src/overset/oversetCommUtilites.F90 @@ -1917,7 +1917,7 @@ subroutine updateOversetConnectivity(level, sps) flowDoms(d1, level, sps)%x(i1-1:i1+1, j1-1:j1+1, k1-1:k1+1, :), frac0, frac) ! Check if the fractions are between 0 and 1 - if (ANY(frac > one + fracTol) .or. ANY(frac < zero - fracTol)) then + if (MAXVAL(frac) > one + fracTol .or. MINVAL(frac) < zero - fracTol) then print *, "Invalid overset connectivity update. Use 'frozen' or 'full' oversetUpdateMode instead." error stop end if @@ -1961,7 +1961,7 @@ subroutine updateOversetConnectivity(level, sps) flowDoms(d2, level, sps)%x(i2-1:i2+1, j2-1:j2+1, k2-1:k2+1, :), frac0, frac) ! Check if the fractions are between zero and one - if (ANY(frac > one + fracTol) .or. ANY(frac < zero - fracTol)) then + if (MAXVAL(frac) > one + fracTol .or. MINVAL(frac) < zero - fracTol) then print *, "Invalid overset connectivity update. Use 'frozen' or 'full' oversetUpdateMode instead." error stop end if From 75df0bb39a810f4b80ca1954c3e3747e18576f38 Mon Sep 17 00:00:00 2001 From: sseraj Date: Tue, 19 Jul 2022 16:05:01 -0400 Subject: [PATCH 12/12] loosened fracTol --- src/overset/oversetCommUtilites.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/overset/oversetCommUtilites.F90 b/src/overset/oversetCommUtilites.F90 index 6422d5cfd..8c6dc4dc4 100644 --- a/src/overset/oversetCommUtilites.F90 +++ b/src/overset/oversetCommUtilites.F90 @@ -1768,7 +1768,7 @@ subroutine updateOversetConnectivity(level, sps) integer(kind=intType), dimension(8), parameter :: indices=(/1,2,4,3,5,6,8,7/) ! Set a tolerance for checking whether fractions are between 0 and 1 - real(kind=realType) :: fracTol=1e-8 + real(kind=realType) :: fracTol=1e-4 ! Pointers to the overset comms to make it easier to read commPattern => commPatternOverset(level, sps)