Add functionality for volumetric multiscale coupling to micro-mechanics simulations

Makefile

@@ -62,7 +62,7 @@ else
 	CC = mpicc
 endif
 
-FFLAGS = -Wall -O3 -fopenmp $(INCLUDES) ${ADDITIONAL_FFLAGS}
+FFLAGS = -Wall -O3 -fopenmp -fallow-argument-mismatch $(INCLUDES) ${ADDITIONAL_FFLAGS}
 # Note for GCC 10 or newer: add -fallow-argument-mismatch in the above flags
 FC = mpifort
 # FC = mpif90
@@ -73,8 +73,8 @@ include $(CCX)/Makefile.inc
 SCCXMAIN = ccx_$(CCX_VERSION).c
 
 # Append additional sources
-SCCXC += nonlingeo_precice.c dyna_precice.c CCXHelpers.c PreciceInterface.c
-SCCXF += getflux.f getkdeltatemp.f
+SCCXC += nonlingeo_precice.c dyna_precice.c CCXHelpers.c PreciceInterface.c linstatic_precice.c
+SCCXF += getflux.f getkdeltatemp.f getc3d8elementgausspointcoords.f getc3d4elementgausspointcoords.f
 
 
 

adapter/CCXHelpers.c

@@ -68,6 +68,17 @@ void getSurfaceElementsAndFaces(ITG setID, ITG *ialset, ITG *istartset, ITG *ien
   }
 }
 
+void getElementsIDs(ITG setID, ITG *ialset, ITG *istartset, ITG *iendset, ITG *elements)
+{
+
+  ITG i, k = 0;
+
+  for (i = istartset[setID] - 1; i < iendset[setID]; i++) {
+    elements[k] = ialset[i];
+    k++;
+  }
+}
+
 void getNodeCoordinates(ITG *nodes, ITG numNodes, int dim, double *co, double *v, int mt, double *coordinates)
 {
 
@@ -409,6 +420,18 @@ int getXloadIndexOffset(enum xloadVariable xloadVar)
   }
 }
 
+void getElementStrain(int strainIdx, int numIPTotal, double *eei, double *strainData)
+{
+  int i, idx;
+  // Loop through all element and respective gauss points
+  for (i = 0; i < numIPTotal; i++) {
+    idx                   = i * 6 + strainIdx; //TODO: Add explanation for 6
+    strainData[i * 3]     = eei[idx];
+    strainData[i * 3 + 1] = eei[idx + 1];
+    strainData[i * 3 + 2] = eei[idx + 2];
+  }
+}
+
 void setXload(double *xload, int *xloadIndices, double *values, int numValues, enum xloadVariable xloadVar)
 {
   ITG i;
@@ -470,6 +493,30 @@ void setNodeDisplacements(double *displacements, ITG numNodes, int dim, int *xbo
   }
 }
 
+void setElementStiffness(int stiffnessIdx, int numIPTotal, double *stiffnessData, double *xstiff)
+{
+  int i, idx;
+  // Loop through all element and respective gauss points
+  for (i = 0; i < numIPTotal; i++) {
+    idx             = i * 27 + stiffnessIdx; //TODO: Add explanation for 27
+    xstiff[idx]     = stiffnessData[i * 3];
+    xstiff[idx + 1] = stiffnessData[i * 3 + 1];
+    xstiff[idx + 2] = stiffnessData[i * 3 + 2];
+  }
+}
+
+void setElementStress(int stressIdx, int numIPTotal, double *stressData, double *stx)
+{
+  int i, idx;
+  // Loop through all element and respective gauss points
+  for (i = 0; i < numIPTotal; i++) {
+    idx          = i * 6 + stressIdx; //TODO: Add explanation for 6
+    stx[idx]     = stressData[i * 3];
+    stx[idx + 1] = stressData[i * 3 + 1];
+    stx[idx + 2] = stressData[i * 3 + 2];
+  }
+}
+
 bool isSteadyStateSimulation(ITG *nmethod)
 {
   return *nmethod == 1;

adapter/CCXHelpers.h

@@ -51,7 +51,18 @@ enum CouplingDataType { TEMPERATURE,
                         DISPLACEMENTS,
                         DISPLACEMENTDELTAS,
                         VELOCITIES,
-                        POSITIONS };
+                        POSITIONS,
+                        STRAIN1TO3,
+                        STRAIN4TO6,
+                        STRESS1TO3,
+                        STRESS4TO6,
+                        CMAT1,
+                        CMAT2,
+                        CMAT3,
+                        CMAT4,
+                        CMAT5,
+                        CMAT6,
+                        CMAT7 };
 
 /**
  * @brief Type of element used for faces mesh, where we assume only one type of element is used.
@@ -102,6 +113,16 @@ ITG getNumSetElements(ITG setID, ITG *istartset, ITG *iendset);
  */
 void getSurfaceElementsAndFaces(ITG setID, ITG *ialset, ITG *istartset, ITG *iendset, ITG *elements, ITG *faces);
 
+/**
+ * @brief Gets the element IDs given a set ID
+ * @param setID: input set id
+ * @param ialset: CalculiX variable
+ * @param istartset: CalculiX variable
+ * @param iendset: CalculiX variable
+ * @param elements: output element IDs
+ */
+void getElementsIDs(ITG setID, ITG *ialset, ITG *istartset, ITG *iendset, ITG *elements);
+
 /**
  * @brief Gets the coordinates of a list of input node IDs
  * @param nodes: input node IDs
@@ -241,6 +262,15 @@ void getXbounIndices(ITG *nodes, ITG numNodes, int nboun, int *ikboun, int *ilbo
  */
 void getXforcIndices(ITG *nodes, ITG numNodes, int nforc, int *ikforc, int *ilforc, int *xforcIndices);
 
+/**
+ * @brief Gets the strain values at each Gauss point of each element
+ * @param strainIdx: CalculiX variable for the index of the strain values
+ * @param numIPTotal: CalculiX variable for the number of elements
+ * @param eei: CalculiX array for the element integration information
+ * @param strainData: CalculiX array for the strain values
+ */
+void getElementStrain(int strainIdx, int numIPTotal, double *eei, double *strainData);
+
 /**
  * @brief Modifies the values of a DFLUX or FILM boundary condition
  * @param xload: CalculiX array for the loads
@@ -317,6 +347,18 @@ void setNodeForces(double *forces, ITG numNodes, int dim, int *xforcIndices, dou
  */
 void setNodeDisplacements(double *displacements, ITG numNodes, int dim, int *xbounIndices, double *xboun);
 
+/**
+ * @brief Modifies the values of stress of the elements at the interface
+ * @param TODO
+ */
+void setElementStiffness(int stiffnessIdx, int numIPTotal, double *stiffnessData, double *xstiff);
+
+/**
+ * @brief Modifies the values of stress of the elements at the interface
+ * @param TODO
+ */
+void setElementStress(int stressIdx, int numIPTotal, double *stressData, double *stx);
+
 /**
  * @brief Returns whether it is a steady-state simulation based on the value of nmethod
  * @param nmethod: CalculiX variable with information regarding the type of analysis

adapter/ConfigReader.cpp

@@ -53,6 +53,10 @@ void ConfigReader_Read(char const *configFilename, char const *participantName,
       interface.facesMeshName = strdup(config["participants"][participantName]["interfaces"][i]["mesh"].as<std::string>().c_str());
     }
 
+    if (config["participants"][participantName]["interfaces"][i]["elements"]) {
+      interface.elementsMeshName = strdup(config["participants"][participantName]["interfaces"][i]["elements"].as<std::string>().c_str());
+    }
+
     std::string patchName = config["participants"][participantName]["interfaces"][i]["patch"].as<std::string>();
     std::transform(patchName.begin(), patchName.end(), patchName.begin(), toupper);
     interface.patchName = strdup(patchName.c_str());

adapter/ConfigReader.h

@@ -13,6 +13,7 @@
 typedef struct InterfaceConfig {
   char * facesMeshName;
   char * nodesMeshName;
+  char * elementsMeshName;
   char * patchName;
   int    map;
   int    numWriteData;