su2code · cvencro · Aug 5, 2019 · Aug 7, 2019 · Aug 7, 2019 · Aug 11, 2019
diff --git a/Common/include/geometry/dual_grid/CPoint.hpp b/Common/include/geometry/dual_grid/CPoint.hpp
@@ -60,6 +60,9 @@ class CPoint {
   su2activevector Volume;                 /*!< \brief Volume or Area of the control volume in 3D and 2D. */
   su2activevector Volume_n;               /*!< \brief Volume at time n. */
   su2activevector Volume_nM1;             /*!< \brief Volume at time n-1. */
+  su2activevector Volume_Old;             /*!< \brief Old containers for Volume. */
+  su2activevector Volume_n_Old;           /*!< \brief Old containers for Volume at time n. */
+  su2activevector Volume_nM1_Old;         /*!< \brief Old containers for Volume at time n-1. */
   su2activevector Periodic_Volume;        /*!< \brief Missing component of volume or area of a control volume on a periodic marker in 3D and 2D. */
 
   su2vector<bool> Domain;                 /*!< \brief Indicates if a point must be computed or belong to another boundary */
@@ -524,9 +527,9 @@ class CPoint {
   inline su2double GetVolume_n(unsigned long iPoint) const { return Volume_n(iPoint); }
 
   /*!
-   * \brief Get the volume of the control volume at time n+1.
+   * \brief Get the volume of the control volume at time n-1.
    * \param[in] iPoint - Index of the point.
-   * \return Volume of the control volume at time n+1
+   * \return Volume of the control volume at time n-1
    */
   inline su2double GetVolume_nM1(unsigned long iPoint) const { return Volume_nM1(iPoint); }
 
@@ -536,10 +539,35 @@ class CPoint {
   void SetVolume_n();
 
   /*!
-   * \brief Set the volume of the control volume at time n+1.
+   * \brief Set the volume of the control volume at time n-1.
    */
   void SetVolume_nM1();
 
+  /*!
+   * \brief Set the volume of the control volume at time n using n-1.
+   */
+  void SetVolume_n_from_OldnM1();
+
+  /*!
+   * \brief Set the volume of the control volume at current time using time n.
+   */
+  void SetVolume_from_Oldn();
+
+  /*!
+   * \brief Set the Volume to Volume_Old.
+   */
+  void SetVolume_Old();
+
+  /*!
+   * \brief Set the Volume_n to Volume_n_Old.
+   */
+  void SetVolume_n_Old();
+
+  /*!
+   * \brief Set the Volume_nM1 to Volume_nM1_Old.
+   */
+  void SetVolume_nM1_Old();
+
   /*!
    * \brief Set the parent control volume of an agglomerated control volume.
    * \param[in] iPoint - Index of the point.

diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp
@@ -2214,7 +2214,8 @@ enum MPI_QUANTITIES {
   MESH_DISPLACEMENTS   = 27,  /*!< \brief Mesh displacements at the interface. */
   SOLUTION_TIME_N      = 28,  /*!< \brief Solution at time n. */
   SOLUTION_TIME_N1     = 29,  /*!< \brief Solution at time n-1. */
-  PRIMITIVE            = 30   /*!< \brief Primitive solution communication. */
+  PRIMITIVE            = 30,  /*!< \brief Primitive solution communication. */
+  SOLUTION_VEL_PRED    = 31
 };
 
 /*!

diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp
@@ -4834,20 +4834,15 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_
 
       Restart_Flow = false;
 
-      if (GetKind_GridMovement() != RIGID_MOTION &&
-          GetKind_GridMovement() != NO_MOVEMENT) {
-        SU2_MPI::Error(string("Dynamic mesh movement currently only supported for the discrete adjoint solver for\n") +
-                       string("GRID_MOVEMENT = RIGID_MOTION."), CURRENT_FUNCTION);
-      }
-
       if (Unst_AdjointIter- long(nTimeIter) < 0){
         SU2_MPI::Error(string("Invalid iteration number requested for unsteady adjoint.\n" ) +
                        string("Make sure EXT_ITER is larger or equal than UNST_ADJOINT_ITER."),
                        CURRENT_FUNCTION);
       }
 
       /*--- If the averaging interval is not set, we average over all time-steps ---*/
-
+    }
+    if (Time_Domain) {
       if (Iter_Avg_Objective == 0.0) {
         Iter_Avg_Objective = nTimeIter;
       }
@@ -9811,10 +9806,12 @@ void CConfig::SetMultizone(CConfig *driver_config, CConfig **config_container){
     switch (config_container[iZone]->GetKind_Solver()) {
     case EULER: case NAVIER_STOKES: case RANS:
     case INC_EULER: case INC_NAVIER_STOKES: case INC_RANS:
-    case NEMO_EULER: case NEMO_NAVIER_STOKES:
+    case NEMO_EULER: case NEMO_NAVIER_STOKES:     
+    case DISC_ADJ_EULER: case DISC_ADJ_NAVIER_STOKES: case DISC_ADJ_RANS:
+    case DISC_ADJ_INC_EULER: case DISC_ADJ_INC_NAVIER_STOKES: case DISC_ADJ_INC_RANS:
       fluid_zone = true;
       break;
-    case FEM_ELASTICITY:
+    case FEM_ELASTICITY: case DISC_ADJ_FEM:
       structural_zone = true;
       Relaxation = true;
       break;

diff --git a/Common/src/geometry/dual_grid/CPoint.cpp b/Common/src/geometry/dual_grid/CPoint.cpp
@@ -68,6 +68,12 @@ void CPoint::FullAllocation(unsigned short imesh, const CConfig *config) {
     Volume_nM1.resize(npoint) = su2double(0.0);
   }
 
+  if (config->GetDynamic_Grid() && config->GetDiscrete_Adjoint()) {
+    Volume_Old.resize(npoint) = su2double(0.0);
+    Volume_n_Old.resize(npoint) = su2double(0.0);
+    Volume_nM1_Old.resize(npoint) = su2double(0.0);
+  }
+
   if(config->GetAD_Mode() && config->GetMultizone_Problem()) {
     AD_InputIndex.resize(npoint,nDim) = 0;
     AD_OutputIndex.resize(npoint,nDim) = 0;
@@ -163,6 +169,31 @@ void CPoint::SetVolume_nM1() {
   parallelCopy(Volume_n.size(), Volume_n.data(), Volume_nM1.data());
 }
 
+void CPoint::SetVolume_Old() {
+  assert(Volume_Old.size() == Volume.size());
+  parallelCopy(Volume.size(), Volume.data(), Volume_Old.data());
+}
+
+void CPoint::SetVolume_n_Old() {
+  assert(Volume_n_Old.size() == Volume_n.size());
+  parallelCopy(Volume_n.size(), Volume_n.data(), Volume_n_Old.data());
+}
+
+void CPoint::SetVolume_nM1_Old() {
+  assert(Volume_nM1_Old.size() == Volume_nM1.size());
+  parallelCopy(Volume_nM1.size(), Volume_nM1.data(), Volume_nM1_Old.data());
+}
+
+void CPoint::SetVolume_n_from_OldnM1() {
+  assert(Volume_n.size() == Volume_nM1_Old.size());
+  parallelCopy(Volume_nM1_Old.size(), Volume_nM1_Old.data(), Volume_n.data());
+}
+
+void CPoint::SetVolume_from_Oldn() {
+  assert(Volume.size() == Volume_n_Old.size());
+  parallelCopy(Volume_n_Old.size(), Volume_n_Old.data(), Volume.data());
+}
+
 void CPoint::SetCoord_n() {
   assert(Coord_n.size() == Coord.size());
   parallelCopy(Coord.size(), Coord.data(), Coord_n.data());

diff --git a/Common/src/grid_movement/CVolumetricMovement.cpp b/Common/src/grid_movement/CVolumetricMovement.cpp
@@ -438,8 +438,10 @@ void CVolumetricMovement::ComputeSolid_Wall_Distance(CGeometry *geometry, CConfi
   for(iMarker=0; iMarker<config->GetnMarker_All(); ++iMarker) {
     if( (config->GetMarker_All_KindBC(iMarker) == EULER_WALL ||
          config->GetMarker_All_KindBC(iMarker) == HEAT_FLUX)  ||
-       (config->GetMarker_All_KindBC(iMarker) == ISOTHERMAL) ||
-        (config->GetMarker_All_KindBC(iMarker) == CHT_WALL_INTERFACE)) {
+        (config->GetMarker_All_KindBC(iMarker) == ISOTHERMAL) ||
+        (config->GetMarker_All_KindBC(iMarker) == CHT_WALL_INTERFACE) ||
+        (config->GetMarker_All_KindBC(iMarker) == CLAMPED_BOUNDARY ) ||
+        (config->GetMarker_All_KindBC(iMarker) == LOAD_BOUNDARY )) {
       nVertex_SolidWall += geometry->GetnVertex(iMarker);
     }
   }
@@ -457,8 +459,10 @@ void CVolumetricMovement::ComputeSolid_Wall_Distance(CGeometry *geometry, CConfi
   for (iMarker=0; iMarker<config->GetnMarker_All(); ++iMarker) {
     if ( (config->GetMarker_All_KindBC(iMarker) == EULER_WALL ||
          config->GetMarker_All_KindBC(iMarker) == HEAT_FLUX)  ||
-       (config->GetMarker_All_KindBC(iMarker) == ISOTHERMAL)  ||
-         (config->GetMarker_All_KindBC(iMarker) == CHT_WALL_INTERFACE)) {
+         (config->GetMarker_All_KindBC(iMarker) == ISOTHERMAL)  ||
+         (config->GetMarker_All_KindBC(iMarker) == CHT_WALL_INTERFACE) ||
+         (config->GetMarker_All_KindBC(iMarker) == CLAMPED_BOUNDARY ) ||
+         (config->GetMarker_All_KindBC(iMarker) == LOAD_BOUNDARY )) {
       for (iVertex=0; iVertex<geometry->GetnVertex(iMarker); ++iVertex) {
         iPoint = geometry->vertex[iMarker][iVertex]->GetNode();
         PointIDs[jj++] = iPoint;

diff --git a/SU2_CFD/include/drivers/CDriver.hpp b/SU2_CFD/include/drivers/CDriver.hpp
@@ -818,7 +818,8 @@ class CDriver {
    */
   void SetSourceTerm_DispAdjoint(unsigned short iMarker, unsigned long iVertex, passivedouble val_AdjointX,
                                  passivedouble val_AdjointY, passivedouble val_AdjointZ);
-
+  void SetSourceTerm_VelAdjoint(unsigned short iMarker, unsigned long iVertex, passivedouble val_AdjointX,
+                                 passivedouble val_AdjointY, passivedouble val_AdjointZ);
   /*!
    * \brief Get the undeformed mesh coordinates
    * \param[in] iMarker - Marker identifier.

diff --git a/SU2_CFD/include/iteration/CAdjFluidIteration.hpp b/SU2_CFD/include/iteration/CAdjFluidIteration.hpp
@@ -67,7 +67,7 @@ class CAdjFluidIteration : public CFluidIteration {
   void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-               unsigned short val_iInst) override;
+               unsigned short val_iInst, bool CrossTerm) override;
 
   /*!
    * \brief Updates the containers for the adjoint fluid system.

diff --git a/SU2_CFD/include/iteration/CDiscAdjFEAIteration.hpp b/SU2_CFD/include/iteration/CDiscAdjFEAIteration.hpp
@@ -90,7 +90,7 @@ class CDiscAdjFEAIteration : public CIteration {
   void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-               unsigned short val_iInst) override;
+               unsigned short val_iInst, bool CrossTerm) override;
 
   /*!
    * \brief Updates the containers for the discrete adjoint mean flow system.

diff --git a/SU2_CFD/include/iteration/CDiscAdjFluidIteration.hpp b/SU2_CFD/include/iteration/CDiscAdjFluidIteration.hpp
@@ -86,7 +86,7 @@ class CDiscAdjFluidIteration : public CIteration {
   void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-               unsigned short val_iInst) override;
+               unsigned short val_iInst, bool CrossTerm) override;
 
   /*!
    * \brief Updates the containers for the discrete adjoint fluid system.

diff --git a/SU2_CFD/include/iteration/CDiscAdjHeatIteration.hpp b/SU2_CFD/include/iteration/CDiscAdjHeatIteration.hpp
@@ -79,7 +79,7 @@ class CDiscAdjHeatIteration : public CIteration {
   void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-               unsigned short val_iInst) override;
+               unsigned short val_iInst, bool CrossTerm) override;
 
   /*!
    * \brief Perform a single iteration of the adjoint fluid system.

diff --git a/SU2_CFD/include/iteration/CFEAIteration.hpp b/SU2_CFD/include/iteration/CFEAIteration.hpp
@@ -60,7 +60,7 @@ class CFEAIteration : public CIteration {
   void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-               unsigned short val_iInst) override;
+               unsigned short val_iInst, bool CrossTerm) override;
 
   /*!
    * \brief Iterate the structural system for a number of Inner_Iter iterations.

diff --git a/SU2_CFD/include/iteration/CFEMFluidIteration.hpp b/SU2_CFD/include/iteration/CFEMFluidIteration.hpp
@@ -75,7 +75,7 @@ class CFEMFluidIteration : public CFluidIteration {
   void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-               unsigned short val_iInst) override;
+               unsigned short val_iInst, bool CrossTerm) override;
 
   /*!
    * \brief Updates the containers for the finite element flow system.

diff --git a/SU2_CFD/include/iteration/CFluidIteration.hpp b/SU2_CFD/include/iteration/CFluidIteration.hpp
@@ -67,7 +67,7 @@ class CFluidIteration : public CIteration {
   void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-               unsigned short val_iInst) override;
+               unsigned short val_iInst, bool CrossTerm) override;
 
   /*!
    * \brief Iterate the fluid system for a number of Inner_Iter iterations.

diff --git a/SU2_CFD/include/iteration/CHeatIteration.hpp b/SU2_CFD/include/iteration/CHeatIteration.hpp
@@ -58,7 +58,7 @@ class CHeatIteration : public CFluidIteration {
   void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-               unsigned short val_iInst) override;
+               unsigned short val_iInst, bool CrossTerm) override;
 
   /*!
    * \brief Updates the containers for the heat system.

diff --git a/SU2_CFD/include/iteration/CIteration.hpp b/SU2_CFD/include/iteration/CIteration.hpp
@@ -136,7 +136,7 @@ class CIteration {
   virtual void Iterate(COutput* output, CIntegration**** integration, CGeometry**** geometry, CSolver***** solver,
                        CNumerics****** numerics, CConfig** config, CSurfaceMovement** surface_movement,
                        CVolumetricMovement*** grid_movement, CFreeFormDefBox*** FFDBox, unsigned short val_iZone,
-                       unsigned short val_iInst){}
+                       unsigned short val_iInst, bool CrossTerm){}
 
   /*!
    * \brief A virtual member.

diff --git a/SU2_CFD/include/solvers/CDiscAdjSolver.hpp b/SU2_CFD/include/solvers/CDiscAdjSolver.hpp
@@ -129,8 +129,9 @@ class CDiscAdjSolver final : public CSolver {
    *        after tape has been evaluated.
    * \param[in] geometry - The geometrical definition of the problem.
    * \param[in] config - The particular config.
+   * \param[in] CrossTerm - Boolean to determine if this is a cross term extraction.
    */
-  void ExtractAdjoint_Solution(CGeometry *geometry, CConfig *config) override;
+  void ExtractAdjoint_Solution(CGeometry *geometry, CConfig *config, bool CrossTerm) override;
 
   /*!
    * \brief A virtual member.

diff --git a/SU2_CFD/include/solvers/CFEASolver.hpp b/SU2_CFD/include/solvers/CFEASolver.hpp
@@ -449,7 +449,10 @@ class CFEASolver : public CSolver {
                     CNumerics *numerics,
                     const CConfig *config,
                     unsigned short val_marker) final;
-
+  void BC_Velocity(CGeometry *geometry,
+                    CNumerics *numerics,
+                    const CConfig *config,
+                    unsigned short val_marker) final;
   /*!
    * \brief Iterate using an implicit Newmark solver.
    * \param[in] geometry - Geometrical definition of the problem.
@@ -638,6 +641,7 @@ class CFEASolver : public CSolver {
    * \param[in] config - Configuration of the problem.
    */
   void PredictStruct_Displacement(CGeometry *geometry, CConfig *config) final;
+  void PredictStruct_Velocity(CGeometry *geometry, CConfig *config) final;
 
   /*!
    * \brief Computation of Aitken's coefficient.

diff --git a/SU2_CFD/include/solvers/CFVMFlowSolverBase.inl b/SU2_CFD/include/solvers/CFVMFlowSolverBase.inl
@@ -1369,7 +1369,7 @@ void CFVMFlowSolverBase<V, FlowRegime>::Pressure_Forces(const CGeometry* geometr
 
     RefTemp = Temperature_Inf;
     RefDensity = Density_Inf;
-    if (dynamic_grid) {
+    if (dynamic_grid && !config->GetMultizone_Problem()) {
       Mach2Vel = sqrt(Gamma * Gas_Constant * RefTemp);
       Mach_Motion = config->GetMach_Motion();
       RefVel2 = (Mach_Motion * Mach2Vel) * (Mach_Motion * Mach2Vel);

diff --git a/SU2_CFD/include/solvers/CMeshSolver.hpp b/SU2_CFD/include/solvers/CMeshSolver.hpp
@@ -87,6 +87,7 @@ class CMeshSolver final : public CFEASolver {
    * \param[in] config - Definition of the particular problem.
    */
   void ComputeGridVelocity(CGeometry *geometry, CConfig *config);
+  void ComputeGridVelocity_FromBoundary(CGeometry **geometry, CNumerics **numerics, CConfig *config);
 
   /*!
    * \brief Update the coarse multigrid levels after the grid movement.
@@ -101,6 +102,7 @@ class CMeshSolver final : public CFEASolver {
    * \param[in] config - Definition of the particular problem.
    */
   void SetBoundaryDisplacements(CGeometry *geometry, CNumerics *numerics, CConfig *config);
+  void SetBoundaryVelocities(CGeometry *geometry, CNumerics *numerics, CConfig *config);
 
 public:
   /*!

diff --git a/SU2_CFD/include/solvers/CSolver.hpp b/SU2_CFD/include/solvers/CSolver.hpp
@@ -985,6 +985,10 @@ class CSolver {
                                    CNumerics *numerics,
                                    const CConfig *config,
                                    unsigned short val_marker) { }
+  inline virtual void BC_Velocity(CGeometry *geometry,
+                                   CNumerics *numerics,
+                                   const CConfig *config,
+                                   unsigned short val_marker) { }
 
   /*!
    * \brief A virtual member.
@@ -3514,6 +3518,8 @@ class CSolver {
    */
   inline virtual void PredictStruct_Displacement(CGeometry *geometry,
                                                  CConfig *config) { }
+  inline virtual void PredictStruct_Velocity(CGeometry *geometry,
+                                                 CConfig *config) { }
 
   /*!
    * \brief A virtual member.
@@ -3775,6 +3781,15 @@ class CSolver {
    * \param[in] config - The particular config.
    */
   inline virtual void ExtractAdjoint_Solution(CGeometry *geometry, CConfig *config){}
+
+  /*!
+   * \brief A virtual member.
+   * \param[in] geometry - The geometrical definition of the problem.
+   * \param[in] solver_container - The solver container holding all solutions.
+   * \param[in] config - The particular config.
+   * \param[in] CrossTerm - Boolean to determine if this is a cross term extraction.
+   */
+  inline virtual void ExtractAdjoint_Solution(CGeometry *geometry, CConfig *config, bool CrossTerm){}
 
   /*!
    * \brief A virtual member.

diff --git a/SU2_CFD/include/variables/CDiscAdjFEABoundVariable.hpp b/SU2_CFD/include/variables/CDiscAdjFEABoundVariable.hpp
@@ -43,6 +43,7 @@ class CDiscAdjFEABoundVariable final : public CDiscAdjFEAVariable {
   MatrixType FlowTraction_Sens;        /*!< \brief Adjoint of the flow tractions. */
   MatrixType SourceTerm_DispAdjoint;   /*!< \brief Source term applied into the displacement
                                                    adjoint coming from external solvers. */
+  MatrixType SourceTerm_VelAdjoint;
 
   CVertexMap<unsigned> VertexMap;      /*!< \brief Object that controls accesses to the variables of this class. */
 
@@ -101,6 +102,10 @@ class CDiscAdjFEABoundVariable final : public CDiscAdjFEAVariable {
     if (!VertexMap.GetVertexIndex(iPoint)) return;
     SourceTerm_DispAdjoint(iPoint,iDim) = val;
   }
+  inline void SetSourceTerm_VelAdjoint(unsigned long iPoint, unsigned long iDim, su2double val) override {
+    if (!VertexMap.GetVertexIndex(iPoint)) return;
+    SourceTerm_VelAdjoint(iPoint,iDim) = val;
+  }
 
   /*!
    * \brief Get the source term applied into the displacement adjoint coming from external solvers
@@ -111,6 +116,10 @@ class CDiscAdjFEABoundVariable final : public CDiscAdjFEAVariable {
     if (!VertexMap.GetVertexIndex(iPoint)) return 0.0;
     return SourceTerm_DispAdjoint(iPoint,iDim);
   }
+  inline su2double GetSourceTerm_VelAdjoint(unsigned long iPoint, unsigned long iDim) const override {
+    if (!VertexMap.GetVertexIndex(iPoint)) return 0.0;
+    return SourceTerm_VelAdjoint(iPoint,iDim);
+  }
 
   /*!
    * \brief Get whether a node is on the boundary