update buoyancy

Author: bigfooted

SU2_CFD/src/iteration/CFluidIteration.cpp

@@ -238,7 +238,7 @@ bool CFluidIteration::Monitor(COutput* output, CIntegration**** integration, CGe
   output->SetHistoryOutput(geometry[val_iZone][val_iInst][MESH_0], solver[val_iZone][val_iInst][MESH_0],
                            config[val_iZone], config[val_iZone]->GetTimeIter(), config[val_iZone]->GetOuterIter(),
                            config[val_iZone]->GetInnerIter());
-  
+
   StopCalc = output->GetConvergence();
 
   /* --- Checking convergence of Fixed CL mode to target CL, and perform finite differencing if needed  --*/
@@ -345,12 +345,12 @@ void CFluidIteration::ComputeTurboPerformance(CSolver***** solver, CGeometry****
       bladesPrimitives.push_back(bladePrimitives);
     }
     TurbomachineryPerformance->ComputeTurbomachineryPerformance(bladesPrimitives);
-    
+
     auto nSpan = config_container[ZONE_0]->GetnSpanWiseSections();
     auto InState = TurbomachineryPerformance->GetBladesPerformances().at(ZONE_0).at(nSpan)->GetInletState();
     nSpan = config_container[nZone-1]->GetnSpanWiseSections();
     auto OutState =  TurbomachineryPerformance->GetBladesPerformances().at(nZone-1).at(nSpan)->GetOutletState();
-    
+
     TurbomachineryStagePerformance->ComputePerformanceStage(InState, OutState, config_container[nZone-1]);
   }
 }
@@ -658,7 +658,6 @@ void CFluidIteration::SetDualTime_Aeroelastic(CConfig* config) const {
         Marker_Tag = config->GetMarker_All_TagBound(iMarker);
         if (Marker_Tag == Monitoring_Tag) { owner = 1; break;
         }           owner = 0;
-       
 
       }
       plunge = config->GetAeroelastic_plunge(iMarker_Monitoring);

TestCases/py_wrapper/custom_source_buoyancy/lam_buoyancy_cavity.cfg

@@ -12,8 +12,7 @@
 %
 SOLVER= INC_NAVIER_STOKES
 KIND_TURB_MODEL= NONE
-MATH_PROBLEM= DIRECT
-RESTART_SOL= NO
+RESTART_SOL= YES
 
 % ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------%
 %
@@ -47,6 +46,7 @@ THERMAL_CONDUCTIVITY_CONSTANT= 0.0246295028571
 
 % ----------------------- BODY FORCE DEFINITION -------------------------------%
 %
+% switch off body force, we now do this with the python wrapper
 BODY_FORCE= NO
 BODY_FORCE_VECTOR= ( 0.0, -9.81, 0.0 )
 
@@ -67,12 +67,12 @@ MARKER_MONITORING= ( NONE )
 % ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------%
 %
 NUM_METHOD_GRAD= GREEN_GAUSS
-CFL_NUMBER= 50
+CFL_NUMBER= 100
 CFL_ADAPT= NO
-CFL_ADAPT_PARAM= ( 1.5, 0.5, 15.0, 1e10)
+CFL_ADAPT_PARAM= ( 0.9, 1.1, 10.0,10000, 1.0e-3)
 MAX_DELTA_TIME= 1E6
 RK_ALPHA_COEFF= ( 0.66667, 0.66667, 1.000000 )
-ITER=6
+ITER=1
 
 % ------------------------ LINEAR SOLVER DEFINITION ---------------------------%
 %
@@ -93,7 +93,7 @@ TIME_DISCRE_FLOW= EULER_IMPLICIT
 % --------------------------- CONVERGENCE PARAMETERS --------------------------%
 %
 CONV_RESIDUAL_MINVAL= -12
-CONV_STARTITER= 10
+CONV_STARTITER= 1
 CONV_CAUCHY_ELEMS= 100
 CONV_CAUCHY_EPS= 1E-6
 
@@ -113,6 +113,9 @@ VOLUME_ADJ_FILENAME= adjoint
 GRAD_OBJFUNC_FILENAME= of_grad.dat
 SURFACE_FILENAME= surface_flow
 SURFACE_ADJ_FILENAME= surface_adjoint
-OUTPUT_WRT_FREQ= 100
-SCREEN_OUTPUT= (OUTER_ITER,INNER_ITER, RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y, RMS_TEMPERATURE)
+OUTPUT_WRT_FREQ= 10
+SCREEN_OUTPUT= (INNER_ITER, RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y, RMS_TEMPERATURE)
+
+WRT_PERFORMANCE= YES
+PYTHON_CUSTOM_SOURCE= YES
 

TestCases/py_wrapper/custom_source_buoyancy/run.py

@@ -26,17 +26,17 @@
 
 import sys
 import pysu2
-from mpi4py import MPI
+import numpy as np
 
+# with mpi:
+from mpi4py import MPI
+comm = MPI.COMM_WORLD
+rank = comm.Get_rank()
+# without mpi:
+#  comm = 0
 
 def main():
-  """
-  custom source to add buoyancy term.
-  """
-  # parallel
-  comm = MPI.COMM_WORLD
-  # serial
-  #comm = 0
+
 
   # Initialize the primal driver of SU2, this includes solver preprocessing.
   try:
@@ -45,83 +45,65 @@ def main():
     print('A TypeError occured in pysu2.CSinglezoneDriver : ', exception)
     raise
 
-  print("\n------------------------------ Begin Solver -----------------------------")
-  sys.stdout.flush()
+  if rank == 0:
+    print("\n------------------------------ Begin Solver -----------------------------")
+    sys.stdout.flush()
 
   # we need to add a source term to the energy equation. For this, we need to get the solver and the variable first.
   # we then loop over all points and for these points, we add the source term
   nDim = driver.GetNumberDimensions()
 
   # index to the flow solver
   iSOLVER = driver.GetSolverIndices()['INC.FLOW']
-  #print("index of flow solver = ",iSOLVER)
 
   # all the indices and the map to the names of the primitives
   primindex = driver.GetPrimitiveIndices()
-  #print("indices of primitives=",primindex)
-  #print("number of primitives:",len(primindex))
-
-  #print("number of elements:",driver.GetNumberElements())
 
   nVars = driver.GetNumberSolverVars(iSOLVER)
-  #print("number of solver variables:",nVars)
   varindex = primindex.copy()
-  for prim in varindex.copy():
-    if varindex[prim] >=nVars:
-      del varindex[prim]
-  varindex = dict(sorted(varindex.items(), key=lambda item: item[1]))
-
+  #for prim in varindex.copy():
+  #  if varindex[prim] >=nVars:
+  #    del varindex[prim]
+  #varindex = dict(sorted(varindex.items(), key=lambda item: item[1]))
 
-
-  print("solver variable names:",varindex)
   iDENSITY = primindex.get("DENSITY")
-  print("index of density = ",iDENSITY)
-
-  index_Vel = varindex.get("VELOCITY_X")
-  print("index of velocity = ",index_Vel)
-  custom_source_vector = [0.0 for i in range(nVars)]
-  print("custom source vector = ", custom_source_vector)
+  #index_Vel = varindex.get("VELOCITY_X")
 
-  #print("max. number of inner iterations: ",driver.GetNumberInnerIter());
-  #print("max nr of outer iterations: ",driver.GetNumberOuterIter());
-
-  # is in domain: isdomain = driver.GetNodeDomain(iPoint)
-  #for i_vertex in range(n_vertex)
-  #AllSolutions = driver.GetAllSolutions(iSOLVER)
   Body_Force_Vector = [0.0, -9.81, 0.0]
   DensityInc_0 = driver.GetDensity_FreeStreamND()
-  #print("rho freestream = ",DensityInc_0)
   Force_Ref = driver.GetForce_Ref()
-  #print("reference force = ",Force_Ref)
 
-  Iter = driver.GetNumberInnerIter()
-  print("1. inner iterations = ",Iter)
+  # super important to actually push the commands.
+  sys.stdout.flush()
+
+  # run N iterations
+  for inner_iter in range(11):
+    if (rank==0):
+      print("python iteration ", inner_iter)
 
-  for inner_iter in range(Iter):
     # set the source term, per point
-    print(driver.GetNumberNodes() - driver.GetNumberHaloNodes())
     for i_node in range(driver.GetNumberNodes() - driver.GetNumberHaloNodes()):
-      #SolutionVector =  driver.GetSolutionVector(iSOLVER,i_node)
-      PrimitiveVector =  driver.GetPrimitiveVector(iSOLVER,i_node)
-      DensityInc_i = PrimitiveVector[iDENSITY]
+      DensityInc_i = driver.Primitives()(i_node,iDENSITY)
 
       for iDim in range(nDim):
-        custom_source_vector[iDim+1] = -(DensityInc_i - DensityInc_0) * Body_Force_Vector[iDim] / Force_Ref
-
-      #driver.SetPointCustomSource(iSOLVER, i_node, custom_source_vector)
+        custom_source_vector = (DensityInc_i - DensityInc_0) * Body_Force_Vector[iDim] / Force_Ref
+        driver.UserDefinedSource(iSOLVER).Set(i_node,iDim+1,custom_source_vector)
 
-    print("   ***   inner iteration:",inner_iter)
     driver.Preprocess(inner_iter)
-
-    # Run one time iteration.
     driver.Run()
+
     driver.Postprocess()
     driver.Update()
 
     # Monitor the solver and output solution to file if required.
-    #driver.Monitor(inner_iter)
+    #stopcalc = driver.Monitor(inner_iter)
     driver.Output(inner_iter)
 
+    #if (stopcalc):
+    #    if (rank==0):
+    #        "Max iterations or convergence criteria reached, stopping."
+    #    break;
+
   # Finalize the solver and exit cleanly.
   driver.Finalize()
 

TestCases/py_wrapper/turbulent_premixed_psi/psi.cfg

@@ -67,7 +67,6 @@ AXISYMMETRIC= YES
 SPECIFIED_INLET_PROFILE= YES
 INLET_MATCHING_TOLERANCE=1e-4
 INLET_FILENAME= inlet.dat
-%INLET_INTERPOLATION_FUNCTION= LINEAR_1D
 INC_INLET_TYPE=  VELOCITY_INLET
 INC_INLET_DAMPING= 0.1
 MARKER_INLET= ( inlet, 673, 40.0, 1.0,  0.0, 0.0)
@@ -82,9 +81,10 @@ MARKER_OUTLET= ( outlet, 0.0 )
 %
 NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES
 %
-CFL_NUMBER= 10.0
+CFL_NUMBER= 2.00
 CFL_REDUCTION_SPECIES= 1.0
 CFL_REDUCTION_TURB= 1.0
+%CFL_ADAPT= YES
 CFL_ADAPT= YES
 CFL_ADAPT_PARAM= ( 0.95, 1.01, 1.0, 250, 1.0e-4, 0)
 
@@ -95,7 +95,7 @@ ITER= 1
 %
 LINEAR_SOLVER= FGMRES
 LINEAR_SOLVER_PREC= ILU
-LINEAR_SOLVER_ERROR= 1E-8
+LINEAR_SOLVER_ERROR= 1E-12
 LINEAR_SOLVER_ITER= 10
 
 % -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------%
@@ -131,7 +131,6 @@ SPECIES_CLIPPING_MAX= 1.0
 %
 CONV_NUM_METHOD_TURB= BOUNDED_SCALAR
 MUSCL_TURB= NO
-
 %
 % --------------------------- CONVERGENCE PARAMETERS --------------------------%
 %
@@ -141,9 +140,7 @@ CONV_STARTITER= 10
 %
 % ------------------------- INPUT/OUTPUT INFORMATION --------------------------%
 %
-%MESH_FILENAME= psi_coarse.su2
 MESH_FILENAME= psi.su2
-%MESH_FILENAME= psi_fine.su2
 %
 SCREEN_OUTPUT= INNER_ITER WALL_TIME \
                RMS_PRESSURE RMS_VELOCITY-X RMS_VELOCITY-Y RMS_TKE RMS_DISSIPATION RMS_SPECIES_0 \
@@ -158,11 +155,12 @@ MARKER_ANALYZE= gas_inlet, air_axial_inlet, outlet
 MARKER_ANALYZE_AVERAGE= AREA
 %
 OUTPUT_FILES= RESTART, PARAVIEW_MULTIBLOCK
-VOLUME_OUTPUT= RESIDUAL, PRIMITIVE
+VOLUME_OUTPUT= RESIDUAL, PRIMITIVE, SPECIES_UDS_0
 OUTPUT_WRT_FREQ= 100
 %
 READ_BINARY_RESTART= YES
 RESTART_FILENAME= restart
 SOLUTION_FILENAME= solution
 %
 WRT_PERFORMANCE= YES
+PYTHON_CUSTOM_SOURCE= YES