HL_SignedDistLevels.cc

#include "HL_SignedDistLevels.h"

#include "VelocityExtFromVelocityField.h"

void 
HL_SignedDistLevels::initializeBoundary() 
{
  FUNCNAME("HL_SignedDistLevels::initializeBoundary()");
  
  TraverseStack stack;
  FixVec<double, VERTEX> distVec(dim, NO_INIT);
  int elStatus;
  const int *elVertStatusVec;
  int NumVertIntPoints;
  int ElNum;
  clock_t time1;
  clock_t time2;
  double TIME;
  
  //=== transvering Mesh to SMI and add quantities
  time1 = clock();
  Mesh_to_SMI_and_quantity();
  time2 = clock();
  TIME = TIME_USED(time1, time2);
  
  cout <<"Zeit zum Transformieren nach SMI: "<<TIME<<" sec.\n";
  
  // ===== All non-boundary vertices are initialized with "infinity". =====
  sD_DOF->set(inftyValue);
  
  // ===== Traverse mesh and initialize boundary elements. =====
  const int nBasFcts = feSpace->getBasisFcts()->getNumber();
  DegreeOfFreedom *locInd = GET_MEMORY(DegreeOfFreedom, nBasFcts);

  ElInfo *elInfo;
  if (velExt && velExtType.isSet(VEL_EXT_FROM_VEL_FIELD)) {
    elInfo = stack.traverseFirst(feSpace->getMesh(),
				 -1, 
				 Mesh::CALL_LEAF_EL | 
				 Mesh::FILL_BOUND |
				 Mesh::FILL_COORDS |
				 Mesh::FILL_GRD_LAMBDA);
  } 
  else {
    elInfo = stack.traverseFirst(feSpace->getMesh(),
				 -1, 
				 Mesh::CALL_LEAF_EL | 
				 Mesh::FILL_BOUND |
				 Mesh::FILL_COORDS);
  }
  while(elInfo) {

    // Set elInfo in case velocity extension from velocity field is used.
    if (velExt && velExtType.isSet(VEL_EXT_FROM_VEL_FIELD)) {
      ((VelocityExtFromVelocityField *)velExt)->setElInfo(elInfo);
    }  
    
    // Get local indices of vertices.
    feSpace->getBasisFcts()->getLocalIndices(
			     const_cast<Element *>(elInfo->getElement()),
			     const_cast<DOFAdmin *>(feSpace->getAdmin()),
			     locInd); 
    
    // Get element status.
    elStatus = elLS->createElementLevelSet(elInfo);
    
    //Get some information for creating the first list
    
    elVertStatusVec = elLS->getElVertStatusVec();
    
    //Beginn creating the first list
    
    NumVertIntPoints = elLS->getNumVertIntPoints();
    ElNum = elInfo->getElement()->getIndex();
    
    createFirstList( elStatus, elVertStatusVec, ElNum, NumVertIntPoints);

    //end creating the first list
    
    // Is element cut by the interface ?
    if (elStatus == ElementLevelSet::LEVEL_SET_BOUNDARY) {
      
      // Reset element distance vector.
      for (int i=0; i<=dim; ++i) {
	distVec[i] = inftyValue;
      }
      
      // Mark all vertices as boundary vertices.
      for (int i=0; i<=dim; ++i) {
	(*bound_DOF)[locInd[i]] = 1.0;
      }
      
      // Calculate distance for all vertices.
      if (bndElDist->calcDistOnBoundaryElement(elInfo, distVec) !=
	  ElementLevelSet::LEVEL_SET_BOUNDARY) {
	ERROR_EXIT("error in distance calculation !\n");
	}
      else {
	
	// If distance is smaller, correct to new distance.
	for (int i=0; i<=dim; ++i) {
	  if ((*sD_DOF)[locInd[i]] > distVec[i]) {
	    (*sD_DOF)[locInd[i]] = distVec[i];
	    //If distance is corrected, calculate new velocity.
	    if(velExt != NULL)
	      {
		velExt->calcVelocityBoundary(locInd, i);
	      }
	  }
	}
      }
    }  // end of: elStatus == ElementLevelSet::LEVEL_SET_BOUNDARY
    
    else if (elLS->getNumVertIntPoints() != 0){
//       // Interface cuts element only in vertices.
//       elVertStatusVec = elLS->getElVertStatusVec();
      
//       for (int i=0; i<=dim; ++i) {
// 	if (elVertStatusVec[i] == ElementLevelSet::LEVEL_SET_BOUNDARY) {
// 	  (*bound_DOF)[locInd[i]] = 1.0;
// 	  (*sD_DOF)[locInd[i]] = 0.0;
// 	  }
//       }

      ERROR_EXIT("intersection point is element vertex. should never happen !\n");
    } 
    
    elInfo = stack.traverseNext(elInfo);
  }  // end of: mesh traverse
  
  FREE_MEMORY(locInd, DegreeOfFreedom, nBasFcts);
  
  return;
}

void 
HL_SignedDistLevels::Mesh_to_SMI_and_quantity()
{
  if(!smiAdapter)
    {
      smiAdapter = NEW SMIAdapter(1, 1,
				  const_cast<FiniteElemSpace*>(feSpace),
				  1, -1);
    }
  
  // cout << "\n\n\tSMI-Adapter angelegt !\n\n";
  
  //====== transfer Mesh to SMI ======================
  
  smiAdapter->addNeighbourInfo();
  
  // cout << "\n\n\tNachbarschafts-Infos in SMI ergaenzt !\n\n";
  
  smiAdapter-> transferMeshToSMI();
  
  //   cout << "\n\n\tGitter nach SMI geschrieben !\n\n";
  
  int smiError = SMI_Begin_write_transaction(1,1);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Begin_write_transaction() failed with error %d\n", smiError);
  
  int nul = 0;
  int *null = new int [dim+1];
  for (int i=0; i<=dim; i++)
    {
      null[i] = nul;
    }

  //which pair of element and node is saved in list
  smiError = 
    SMI_Add_quantity(1, 1, 1, SMI_LOC_ELEM, SMI_TYPE_INT, dim+1, null);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Add_quantity() failed with error %d\n", smiError);

  //which node is a boundary-node
  smiError = SMI_Add_quantity(1, 1, 2, SMI_LOC_NODE, SMI_TYPE_INT, 1, &nul);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Add_quantity() failed with error %d\n", smiError);

  //saves the number of tried updates
  smiError = SMI_Add_quantity(1, 1, 3, SMI_LOC_NODE, SMI_TYPE_INT, 1, &nul);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Add_quantity() failed with error %d\n", smiError);

  //saves the number of updates
  smiError = SMI_Add_quantity(1, 1, 4, SMI_LOC_NODE, SMI_TYPE_INT, 1, &nul);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Add_quantity() failed with error %d\n", smiError);

  //how often a node is saved in the second list
  smiError = SMI_Add_quantity(1, 1, 5, SMI_LOC_NODE, SMI_TYPE_INT, 1, &nul);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Add_quantity() failed with error %d\n", smiError);

  //saves the level of an element
  smiError = SMI_Add_quantity(1, 1, 6, SMI_LOC_ELEM, SMI_TYPE_INT, 1, 
			      &inftyValue);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Add_quantity() failed with error %d\n", smiError);

  // cout << "\n\n\tQuantities in SMI ergaenzt !\n\n";

  delete [] null;
}

void 
HL_SignedDistLevels::createFirstList(const int elStatus, 
				     const int *elVertStatusVec, 
				     int ElNum, 
				     const int NumVertIntPoints)
{
  int *sv = new int[dim+1];
  int *nodeIndices;
  int nul = 0;
  int smiError;
  
  if (elStatus == ElementLevelSet::LEVEL_SET_BOUNDARY)
    {
      List_Element.push(ElNum);

      smiError = 
	SMI_Set_quantity_values(1, 1, 6, SMI_TYPE_INT, 1, 1, &ElNum, &nul);
      TEST_EXIT(smiError == SMI_OK)
	("SMI_Set_quantity_values() failed with error %d\n", smiError);

      //saving which nodes are boundary-nodes
      for (int i=0; i<=dim; i++)
	{
	  sv[i] = 1;
	}

      smiError =
	SMI_Get_elems(1, 1, 1, const_cast<DegreeOfFreedom*>(&ElNum),
		      NULL, &nodeIndices, NULL, NULL);
      TEST_EXIT(smiError == SMI_OK)
	("SMI_Get_elems() failed with error %d\n", smiError);
      
      smiError = 
	SMI_Set_quantity_values(1, 1, 2, SMI_TYPE_INT, 1, 3, nodeIndices, sv);
      TEST_EXIT(smiError == SMI_OK)
	("SMI_Set_quantity_values() failed with error %d\n", smiError);
    }
  //if the elemet isn't a boundary-element, but the interface cuts the FE in two nodes
//   else if (NumVertIntPoints == dim)
//     {
//       List_Element.push( ElNum );           

//       smiError = 
// 	SMI_Set_quantity_values(1, 1, 6, SMI_TYPE_INT, 1, 1, &ElNum, &nul);
//       TEST_EXIT(smiError == SMI_OK)
// 	("SMI_Set_quantity_values() failed with error %d\n", smiError);
      
//       //saving which nodes are boundary-nodes
//       for (int i=0; i<=dim; i++)
// 	{
// 	  if(elVertStatusVec[i] == ElementLevelSet::LEVEL_SET_BOUNDARY)
// 	    {
// 	      sv[i] = 1;
// 	    }
// 	  else
// 	    {
// 	      sv[i] = 0;
// 	    }

// 	  smiError = 
// 	    SMI_Get_elems(1, 1, 1, const_cast<DegreeOfFreedom*>(&ElNum),
// 			  NULL, &nodeIndices, NULL, NULL);
// 	  TEST_EXIT(smiError == SMI_OK)
// 	    ("SMI_Get_elems() failed with error %d\n", smiError);

// 	  smiError = 
// 	    SMI_Set_quantity_values(1, 1, 2, SMI_TYPE_INT, 1, 3,
// 				    nodeIndices, sv);
// 	  TEST_EXIT(smiError == SMI_OK)
// 	    ("SMI_Set_quantity_values() failed with error %d\n", smiError);
// 	}
//     }
  
  delete [] sv;
}
  
void 
HL_SignedDistLevels::HL_updateIteration()
{
  int numNodes;
  int *nodeIndices;
  int max_q3 = 0;
  int max_q4 = 0;
  int sum_q3 = 0;
  int sum_q4 = 0;
  clock_t time1;
  clock_t time2;
  double TIME;
  int control;
  
  int smiError = 
    SMI_Get_all_nodes(1, 1, const_cast<DegreeOfFreedom*>(&numNodes),
		      &nodeIndices);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Get_all_nodes() failed with error %d\n", smiError);

  int *values_q3 = new int [numNodes];
  int *values_q4 = new int [numNodes];
  
  GET_PARAMETER(0,"SignedDist->count updates", "%d", &count_updates);
  
  time1 = clock();
  control = traverseListElement();
  time2 = clock();
  TIME = TIME_USED(time1, time2);
  
  cout<<"Zeit zum durchlaufen der ersten (Elementen-) Liste: "<<TIME<<" sec.\n\n";
  if(control == 1)
    {
      createLevels ();
    }
  
  /////////////////////////////////
  
  
  //print chosen level in a file
  if (print_level == 1)
    {
      DOFVector<double> *level_DOF = NEW DOFVector<double>(sD_DOF->getFeSpace(), 
							   "level_DOF");
      
      int numElems;
      int *elemIndices;
      int *NodeIndices;
      
      smiError = SMI_Get_all_elems(1, 1, &numElems, &elemIndices);
      TEST_EXIT(smiError == SMI_OK)
	("SMI_Get_all_elems() failed with error %d\n", smiError);
      
      int *value_q6 = new int [numElems];
      
      smiError =
	SMI_Get_quantity_values(1, 1, 6, SMI_TYPE_INT, 1, numElems,
				elemIndices, value_q6);
      TEST_EXIT(smiError == SMI_OK)
	("SMI_Get_quantity_values() failed with error %d\n", smiError);
      
      for(int i=0; i<numElems; i++)
	{	  
	  smiError = 
	    SMI_Get_elems(1, 1, 1, &elemIndices[i], NULL, &NodeIndices, 
			  NULL, NULL);
	  TEST_EXIT(smiError == SMI_OK)
	    ("SMI_Get_elems() failed with error %d\n", smiError);

	  for(int j=0; j<=dim; j++)
	    {
	      (*level_DOF)[NodeIndices[j]] = -1;
	    }
	}
      for(int i=0; i<numElems; i++)
	{	  
	  smiError = 
	    SMI_Get_elems(1, 1, 1, &elemIndices[i], NULL, &NodeIndices, 
			  NULL, NULL);
	  TEST_EXIT(smiError == SMI_OK)
	    ("SMI_Get_elems() failed with error %d\n", smiError);

	  if (value_q6[i] == chosen_level)
	    {
	      for(int j=0; j<dim+1; j++)
		{
		  (*level_DOF)[NodeIndices[j]] = 0;
		}
	    }
	}
      
      FileWriter *levelFileWriter = 
	NEW FileWriter("level->output",
		       level_DOF->getFeSpace()->getMesh(),
		       level_DOF);
      levelFileWriter->writeFiles (adaptInfo, false);
      
      DELETE levelFileWriter;
    }
  
  //end of the part for printing chosen levels
  
  
  time1 = clock();
  traversingListLevel_ ( sD_DOF );
  time2 = clock();
  TIME = TIME_USED(time1,time2);
  
  cout<<"Zeit zum Durchlaufen der zweiten  Liste (Liste von Listen): "<<TIME<<" sec.\n\n";
  
  std::string smiOutFile;
  GET_PARAMETER(0, "SignedDist->count updates->output->filename", &smiOutFile);
  cout << "count updates Ausgabe-Datei:  " << smiOutFile.c_str() << "\n\n";
  
  ofstream out (smiOutFile.c_str());
  
  smiError = 
    SMI_Get_quantity_values(1, 1, 3, SMI_TYPE_INT, 1, numNodes, nodeIndices,
			    values_q3);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Get_quantity_values() failed with error %d\n", smiError);

  smiError =
    SMI_Get_quantity_values(1, 1, 4, SMI_TYPE_INT, 1, numNodes, nodeIndices,
			    values_q4);  
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Get_quantity_values() failed with error %d\n", smiError);
  
  for (int i=0; i<numNodes; i++)
    {
      out<<nodeIndices[i]<<" "<<values_q3[i]<<" "<<values_q4[i]<<"\n";
      
      if(max_q3 < values_q3[i])
	{
	  max_q3 = values_q3[i];
	}
      if(max_q4 < values_q4[i])
	{
	  max_q4 = values_q4[i];
	}
      sum_q3 = sum_q3 + values_q3[i];
      sum_q4 = sum_q4 + values_q4[i];
    }
  
  out<<"\n\n maximale Anzahl an versuchten Updates auf einem Knoten: "<<max_q3<<"\n maximale Anzahl an durchgefuehrten Updates auf einem Knoten:  "<<max_q4<<"\n";
  out<<"\n Summe aller versuchten Updates: "<<sum_q3<<"\n Summe aller durchgefuehrten updates: "<<sum_q4<<"\n";
  
  out<<"Anzahl an Knoten: "<<numNodes<<"\n";
  
  out.close();
  
  smiError = SMI_End_write_transaction(1, 1);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_End_write_transaction() failed with error %d\n", smiError);
  
  return;
}
 
int 
HL_SignedDistLevels::traverseListElement()
{
  int Element;
  
  bool neighboursExist;
  bool neighbourInNewListSet = false;
  
  while (List_Element.size() != 0)
    {
      Element = List_Element.front();
      
      neighboursExist =  collectNeighbours_setLevels (Element, 0, &neighbourInNewListSet);
      
      List_Element.pop();
    }
  
  //within the function "collectNeighbours_setLevels" the new pairs are saved in the list "level"
  //this list will be one entry of the list "Level_"
  if (neighbourInNewListSet)
    {
      Level_.push_back (Level);
    }
  
  while (Level.size() != 0)
    {
      Level.pop();
    }
  if (!neighbourInNewListSet)
    {
      return 0;
    }
  else
    {
      return 1;
    }
  
}

void 
HL_SignedDistLevels::createLevels ()
{
  bool neighbourExists = true;
  bool stop = false;
  Vert_Struct vertStruct;
  int i = 0;
  bool elementInNextListSet = true;
  int value_q2;
  int smiError;
  
  while(elementInNextListSet)
    {
      stop = false;
      neighbourExists = false;
      elementInNextListSet = false;
      
      while(Level_[i].size() != 0) 
	{
	  vertStruct =  Level_[i].front();
	  
	  neighbourExists = collectNeighbours_setLevels (vertStruct.ElNum, i+1, &elementInNextListSet);
	  
	  smiError = 
	    SMI_Get_quantity_values(1, 1, 2, SMI_TYPE_INT, 1, 1, 
				    &vertStruct.VertNum, &value_q2);
	  TEST_EXIT(smiError == SMI_OK)
	    ("SMI_Get_quantity_values() failed with error %d\n", smiError);

	  if(value_q2 == 0)
	    {
	      helpLevel.push (vertStruct);
	    }
	  Level_[i].pop();
	}
      
      Level_[i] = helpLevel;
      while (helpLevel.size() != 0)
	{
	  helpLevel.pop();
	}
      
      if(elementInNextListSet)
	{
	  Level_.push_back (Level);
	}
      
      while (Level.size() != 0)
	{
	  Level.pop();
	}
      i++;
    }
}

bool 
HL_SignedDistLevels::collectNeighbours_setLevels (const int Element,
						  const int currentIndex,
						  bool *elementInNewListSet)
{
  Vert_Struct vertStruct;
  
  int *neighbour=new int[dim+1];
  int *oppVertices=new int[dim+1];
  int *value = new int [dim+1];
  int saved_level;
  int value_q5;
  int *nodeIndices;
  int *nodeIndicesOfElem;
  int smiError;

  bool neighbourExists = false;
  
  smiAdapter->getNeighbourInfo(Element, neighbour, oppVertices);
  
  for ( int i=0; i<=dim; i++)
    {
      //if this neighbour exists it will be saved in a list according to it's level
      if(neighbour[i] != -1)
	{
	  neighbourExists = true;

	  //calling old level and saving new level if necessary
	  smiError =
	    SMI_Get_quantity_values(1, 1, 6, SMI_TYPE_INT, 1, 1, 
				    &neighbour[i], &saved_level);
	  TEST_EXIT(smiError == SMI_OK)
	    ("SMI_Get_quantity_values() failed with error %d\n", smiError);

	  if(saved_level > currentIndex)
	    {
	      saved_level = currentIndex+1;

	      smiError =
		SMI_Set_quantity_values (1, 1, 6, SMI_TYPE_INT, 1, 1, 
					 &neighbour[i], &saved_level);
	      TEST_EXIT(smiError == SMI_OK)
		("SMI_Set_quantity_values() failed with error %d\n", 
		 smiError);
	      
	      smiError =
		SMI_Get_elems(1, 1, 1, &neighbour[i], NULL, &nodeIndices, 
			      NULL,NULL);
	      TEST_EXIT(smiError == SMI_OK)
		("SMI_Get_elems() failed with error %d\n", 
		 smiError);
	      
	      smiError = 
		SMI_Get_quantity_values(1, 1, 1, SMI_TYPE_INT, dim+1, 1,
					&neighbour[i], value);
	      TEST_EXIT(smiError == SMI_OK)
		("SMI_Get_quantity_values() failed with error %d\n", 
		 smiError);

	      if (value[oppVertices[i]] == 0)
		{
		  vertStruct.ElNum = neighbour[i] ;
		  vertStruct.VertNum = oppVertices[i];                           
		  Level.push( vertStruct );
		  value[oppVertices[i]] = 1;
		  (*elementInNewListSet) = true;
		  
		  //counts how often a node is saved in the list of lists
		  if(count_in_list == 1)
		    {
		      smiError = 
			SMI_Get_elems(1, 1, 1,
				      const_cast<DegreeOfFreedom*>(&(vertStruct.ElNum)), 
				      NULL, &nodeIndicesOfElem, NULL, NULL);
		      TEST_EXIT(smiError == SMI_OK)
			("SMI_Get_elems() failed with error %d\n", smiError);

		      smiError = 
			SMI_Get_quantity_values(1, 1, 5, SMI_TYPE_INT, 1, 1,
						const_cast<DegreeOfFreedom*>(& nodeIndicesOfElem[vertStruct.VertNum]),
						&value_q5);
		      TEST_EXIT(smiError == SMI_OK)
			("SMI_Get_quantity_values() failed with error %d\n", 
			 smiError);

		      value_q5 = value_q5 + 1;

		      smiError =
			SMI_Set_quantity_values(1, 1, 5, SMI_TYPE_INT, 1, 1,
						const_cast<DegreeOfFreedom*>(& nodeIndicesOfElem[vertStruct.VertNum]),
						&value_q5);
		      TEST_EXIT(smiError == SMI_OK)
			("SMI_Set_quantity_values() failed with error %d\n", 
			 smiError);
		    }
		}
	      //save which pairs of elements and nodes is saved in a list
	      smiError = 
		SMI_Set_quantity_values(1, 1, 1, SMI_TYPE_INT, dim+1, 1,
// 					&neighbour[i], value);
					neighbour+i, value);
	      TEST_EXIT(smiError == SMI_OK)
		("SMI_Set_quantity_values() failed with error %d\n", 
		 smiError);
	    }
	}
    }
  
  delete [] neighbour;
  delete [] oppVertices;
  delete [] value;
  
  return neighbourExists;
  
  
}

void 
HL_SignedDistLevels::traversingListLevel_ ( DOFVector<double> *boundVal_DOF )
{
  bool Continue = true;
  
  Vert_Struct El_Vert;
  
  int Vert;
  
  int locVert;
  int value_q3;
  int value_q4;
  int value_q5;
  int counter = 0;
  int *nodeIndices;
  int *nodeIndicesOfElem;
  double *coords = new double [(dim+1)*dim];
  int *valuesINT=new int[dim+1];
  double update;
  int counter_list = 0;
  int counter_lists_of_lists = 0;
  FixVec<WorldVector<double> *, VERTEX> coordsOfNodes(dim);
  for ( int i=0; i<=dim; i++)
    {
      coordsOfNodes[i] = new WorldVector<double>;
    }
  FixVec<double,VERTEX> boundVal(dim,NO_INIT);
  WorldVector<double> helpVec;
  int smiError;
  
  //for repeating the traverse of the list "Level_"
  while (Continue == true)
    {
      Continue = false;
      
      //for traversing of the list "Level_"
      for (int j=0; j < (signed int) Level_.size();j++)
	{ 
	  //for traversing the entries of the List "Level_"                    

	  while (Level_[j].size() != 0)
	    {
	      El_Vert = Level_[j].front();
	      Vert = El_Vert.VertNum;
	      
	      smiError = 
		SMI_Get_elems (1, 1, 1, &(El_Vert.ElNum), NULL, 
			       &(nodeIndices), NULL, NULL);
	      TEST_EXIT(smiError == SMI_OK)
		("SMI_Get_elems() failed with error %d\n", smiError);

	      smiError =
		SMI_Get_nodes (1, 1, 3, dim, nodeIndices, coords);
	      TEST_EXIT(smiError == SMI_OK)
		("SMI_Get_nodes() failed with error %d\n", smiError);
	      
	      //sorting of the nodes for calculating the update
	      counter = 0;
	      for (int i=0; i<=dim; i++)
		{
		  if (i == El_Vert.VertNum)      
		    {
		      (*(coordsOfNodes[dim]))[0]=coords[i*2];
		      (*(coordsOfNodes[dim]))[1]=coords[i*2+1];
		      boundVal[dim] = (*boundVal_DOF)[nodeIndices[i]];
		      locVert = i;
		    }
		  else
		    {
		      (*(coordsOfNodes[counter]))[0]=coords[i*2];
		      (*(coordsOfNodes[counter]))[1]=coords[i*2+1];
		      boundVal[counter] = (*boundVal_DOF)[nodeIndices[i]];
		      counter = counter+1;
		    }
		}
	      
	      //save permutation of vertexes for calculation of the velocity
	      if(velExt != NULL)
		{
		  velExt->setPermutation(locVert, 0);
		}
	      update = elUpdate->calcElementUpdate( coordsOfNodes, boundVal );
	      
	      //for counting tried updates
	      if(count_updates == 1)
		{
		  smiError =
		    SMI_Get_quantity_values(1, 1, 3, SMI_TYPE_INT, 1, 1,
					    &nodeIndices[locVert], &value_q3);
		  TEST_EXIT(smiError == SMI_OK)
		    ("SMI_Get_quantity_values() failed with error %d\n", 
		     smiError);

		  value_q3 = value_q3 + 1;

		  smiError = 
		    SMI_Set_quantity_values(1, 1, 3, SMI_TYPE_INT, 1, 1,
					    const_cast<DegreeOfFreedom*>(&nodeIndices[locVert]),
					    &value_q3);
		  TEST_EXIT(smiError == SMI_OK)
		    ("SMI_Set_quantity_values() failed with error %d\n", 
		     smiError);
		}
	      
	      //checking whether the calculated update will be set as new distance
	      if (update < (*boundVal_DOF)[nodeIndices[locVert]] && ((*boundVal_DOF)[nodeIndices[locVert]]-update) > tol)
		{
		  (*boundVal_DOF)[nodeIndices[locVert]] = update;
		  //If distance is corrected, calculate new velocity.
		  if(velExt != NULL)
		    {
		      velExt->calcVelocity(nodeIndices, locVert);
		    }

		  //for counting successful updates
		  if(count_updates == 1)
		    {
		      smiError =
			SMI_Get_quantity_values(1, 1, 4, SMI_TYPE_INT, 1, 1,
						&nodeIndices[locVert],
						&value_q4); 
		      TEST_EXIT(smiError == SMI_OK)
			("SMI_Get_quantity_values() failed with error %d\n", 
			 smiError);

		      value_q4 = value_q4 + 1;

		      smiError =
			SMI_Set_quantity_values(1, 1, 4, SMI_TYPE_INT, 1, 1,
						const_cast<DegreeOfFreedom*>(&nodeIndices[locVert]),
						&value_q4);
		      TEST_EXIT(smiError == SMI_OK)
			("SMI_Set_quantity_values() failed with error %d\n", 
			 smiError);
		    }
		  //??
		  Continue = search_and_include_comb(El_Vert.ElNum, 
						     Vert, nodeIndices, j);
		} 
	      
	      smiError =
		SMI_Get_quantity_values(1, 1, 1, SMI_TYPE_INT, dim+1, 1,
					&El_Vert.ElNum, valuesINT);
	      TEST_EXIT(smiError == SMI_OK)
		("SMI_Get_quantity_values() failed with error %d\n", 
		 smiError);

	      Level_[j].pop();
	      valuesINT[locVert] = 0;

	      smiError =
		SMI_Set_quantity_values(1, 1, 1, SMI_TYPE_INT, dim+1, 1,
					&El_Vert.ElNum, valuesINT);
	      TEST_EXIT(smiError == SMI_OK)
		("SMI_Set_quantity_values() failed with error %d\n", 
		 smiError);
	      
	      //counts how often a node is saved in the list of lists
	      if(count_in_list == 1)
		{
		  smiError =
		    SMI_Get_elems(1, 1, 1, 
				  const_cast<DegreeOfFreedom*>(&(El_Vert.ElNum)), 
				  NULL, &nodeIndicesOfElem, NULL, NULL);
		  TEST_EXIT(smiError == SMI_OK)
		    ("SMI_Get_elems() failed with error %d\n", smiError);
		  
		  smiError =
		    SMI_Get_quantity_values(1, 1, 5, SMI_TYPE_INT, 1, 1,
					    const_cast<DegreeOfFreedom*>(& nodeIndicesOfElem[El_Vert.VertNum]),
					    &value_q5);
		  TEST_EXIT(smiError == SMI_OK)
		    ("SMI_Get_quantity_values() failed with error %d\n", 
		     smiError);

		  value_q5 = value_q5 - 1;

		  smiError =
		    SMI_Set_quantity_values(1, 1, 5, SMI_TYPE_INT, 1, 1,
					    const_cast<DegreeOfFreedom*>(& nodeIndicesOfElem[El_Vert.VertNum]),
					    &value_q5);
		  TEST_EXIT(smiError == SMI_OK)
		    ("SMI_Set_quantity_values() failed with error %d\n", 
		     smiError);
		}
	      
	      //for printing which node is how often saved in a list after a special number of tried updates
	      if(count_in_list == 1)
		{
		  counter_list++;
		  if(counter_list % print_in_list == 0)
		    {
		      print_quantity_5(counter_list);
		    }
		}	 
	    }
	}
      //for printing which node is how often saved in a list at the end of the ..th traverse of the list "Level_"
      if(count_in_list == 1)
	{
	  counter_lists_of_lists++;
	  if(counter_lists_of_lists == print_in_list_2)
	    {
	      print_quantity_5(counter_lists_of_lists);
	    }
	}
      
    }
  
  delete [] coords;
  delete [] valuesINT;
  for ( int i=0; i<Global::getGeo(VERTEX,dim); i++)
    {
      delete coordsOfNodes[i];
    }
}

bool 
HL_SignedDistLevels::search_and_include_comb(int ElNum, 
					     int Vert, 
					     int *nodeIndices, 
					     const int Index)
{
  bool stop_l = false;
  bool stop_r = false;
  bool Continue = false;
  bool continue_out = false;
  int neighbour_l;
  int neighbour_r;
  int Vert_1_l;
  int Vert_2_l;
  int Vert_1_r;
  int Vert_2_r;
  int counter;
  int elem_l;
  int elem_r;
  int node_l;
  int node_r;
  int VertGlobal;
  int *nodeIndicesOfElem;
  
  int smiError = 
    SMI_Get_elems(1, 1, 1, const_cast<DegreeOfFreedom*>(&(ElNum)), NULL, 
		  &nodeIndicesOfElem, NULL, NULL);
  TEST_EXIT(smiError == SMI_OK)
    ("SMI_Get_elems() failed with error %d\n", smiError);

  VertGlobal = nodeIndicesOfElem[Vert];
  
  //here the process of including elements/nodes into the second list
  stop_l = false;
  stop_r = false;	      
  
  elem_l = ElNum;
  if (elem_l == -1)
    {
      stop_l = true;
    }
  elem_r = ElNum;
  if (elem_r == -1)
    {
      stop_r = true;
    }
  
  counter = 0;
  for (int i=0; i<=dim; i++)
    {
      if(i != Vert && counter == 1)  
	{
	  node_r = i;                
	  counter = counter + 1;
	}
      if(i != Vert && counter == 0)  
	{
	  node_l = i;                
	  counter = counter + 1;
	}
      
    }
  
  while (!stop_l || !stop_r)
    {
      //get next neighbour
      if (!stop_l)
	{
	  stop_l = getNextNeighbour (elem_l, VertGlobal, node_l, neighbour_l, Vert_1_l, Vert_2_l);
	}
      if(!stop_r)
	{
	  stop_r = getNextNeighbour (elem_r, VertGlobal, node_r, neighbour_r, Vert_1_r, Vert_2_r);
	}
      
      if(neighbour_l == elem_r)
	{
	  break;
	}
      
      //indclude into the second list (only if possible)
      if (!stop_l)
	{
	  Continue = includeIntoList (neighbour_l, Vert_1_l, Vert_2_l, Index);
	  if(Continue == true)
	    {
	      continue_out = Continue;
	    }
	  if(neighbour_l == neighbour_r)
	    {
	      break;
	    }
	  //"elem_l", "node_l", "elem_r", "node_r" have to be set on the next elements
	  node_l = getNext_node_l_r(elem_l,neighbour_l, node_l,VertGlobal);
	  elem_l = neighbour_l;
	  
	}
      if (!stop_r && neighbour_l != neighbour_r)
	{
	  Continue = includeIntoList (neighbour_r, Vert_1_r, Vert_2_r, Index);
	  if(Continue == true)
	    {
	      continue_out = Continue;
	    }
	  //"elem_l", "node_l", "elem_r", "node_r" have to be set on the next elements
	  node_r = getNext_node_l_r(elem_r,neighbour_r, node_r, VertGlobal);
	  elem_r = neighbour_r;
	}
    }
  return continue_out;                  
}

bool 
HL_SignedDistLevels::includeIntoList (int ElNum_in, 
				      int VertNum_1_in, 
				      int VertNum_2_in, 
				      const int Index)
{
  Vert_Struct vertStruct_0;
  Vert_Struct vertStruct_1;
  
  int *nodeIndices;
  int *valuesINT = new int[dim+1];
  int *locVert = new int[dim];
  int value_q2;
  int value_q5;
  int value_q6;
  int smiError;
  
  bool Continue = false;
  
  if(ElNum_in != -1)
    {