From 6c0542745423c1b8eb6be8a01f6aa477c262eeb1 Mon Sep 17 00:00:00 2001
From: secher <secher>
Date: Fri, 16 Feb 2007 09:39:20 +0000
Subject: [PATCH] gradient calculate and fields on cells

---
 src/FILTERGUI/SelectParams.cxx | 103 +++++++++++++++++++++++----------
 1 file changed, 73 insertions(+), 30 deletions(-)

diff --git a/src/FILTERGUI/SelectParams.cxx b/src/FILTERGUI/SelectParams.cxx
index 08a0beb..08d1f42 100644
--- a/src/FILTERGUI/SelectParams.cxx
+++ b/src/FILTERGUI/SelectParams.cxx
@@ -651,7 +651,7 @@ void SelectParams::process()
     QMessageBox::information( this,
 			      "Filtering",
 			      "Unable to process filtering.\n"
-			      "You must select a reference field on nodes." );
+			      "You must select a reference field on nodes or on cells." );
     reject();
     return;
   }
@@ -661,7 +661,6 @@ void SelectParams::process()
 
   // call ensight driver MED to generate input ensight mesh ,
   // input ensight boundary mesh and input criteria ensight field for filtoo
-
   int id;
   ENSIGHT_MESH_DRIVER myMeshDriver("/tmp/input.geom",_mesh);
   myMeshDriver.addSupport(supB);
@@ -699,37 +698,20 @@ void SelectParams::process()
 void SelectParams::generateCriteria() throw(SALOME_Exception)
 {
   double val, min, max;
+  double sigmaV;
   bool isGVal;
   MED_EN::medEntityMesh typ;
+  const int *revC;
+  const int *indC;
+  FIELD<double> *volume;
+  set <int> listElements;
+  set <int>::iterator elemIt ;
 
-  // if reference field is on elements create reference field on nodes
   if(_myDField)
     typ = _myDField->getSupport()->getEntity();
   else
     typ = _myIField->getSupport()->getEntity();
 
-  switch(typ){
-  case MED_CELL:
-    MESSAGE("REFERENCE FIELD ON CELL");
-    throw SALOME_Exception("Filter run only on reference field on nodes");
-    break;
-  case MED_FACE:
-    MESSAGE("REFERENCE FIELD ON FACE");
-    throw SALOME_Exception("Filter run only on reference field on nodes");
-    break;
-  case MED_EDGE:
-    MESSAGE("REFERENCE FIELD ON EDGE");
-    throw SALOME_Exception("Filter run only on reference field on nodes");
-    break;
-  case MED_NODE:
-    MESSAGE("REFERENCE FIELD ON NODES");
-    break;
-  case MED_ALL_ENTITIES:
-    MESSAGE("REFERENCE FIELD ON ALL ELEMENTS");
-    throw SALOME_Exception("Filter run only on reference field on nodes");
-    break;
-  }
-
   // create support on nodes
   SUPPORT *sup = new SUPPORT(_mesh,"Support",MED_NODE);
 
@@ -741,13 +723,74 @@ void SelectParams::generateCriteria() throw(SALOME_Exception)
   // read number of nodes
   int NumberOf = sup->getNumberOfElements(MED_ALL_ELEMENTS);
 
+  switch(typ){
+  case MED_CELL:
+    if(_myDField){
+      // calculate reverse connectivity to have the list of elements which contains node i
+      revC = _myDField->getSupport()->getMesh()->getReverseConnectivity(MED_NODAL,typ);
+      indC = _myDField->getSupport()->getMesh()->getReverseConnectivityIndex(MED_NODAL,typ);
+      // calculate volume field on mesh
+      volume = _myDField->getSupport()->getMesh()->getVolume(_myDField->getSupport());
+    }
+    else{
+      // calculate reverse connectivity to have the list of elements which contains node i
+      revC = _myIField->getSupport()->getMesh()->getReverseConnectivity(MED_NODAL,typ);
+      indC = _myIField->getSupport()->getMesh()->getReverseConnectivityIndex(MED_NODAL,typ);
+      // calculate volume field on mesh
+      volume = _myIField->getSupport()->getMesh()->getVolume(_myIField->getSupport());
+    }
+    break;
+  default:
+    break;
+  }
+
   for (int i=1; i<NumberOf+1; i++){
 
-    // read reference field value
-    if(_myDField)
-      val = _myDField->getValueIJ(i,1);
-    else
-      val = (double)_myIField->getValueIJ(i,1);
+    // if reference field is on elements get reference field on nodes
+    switch(typ){
+    case MED_CELL:
+      // listElements contains elements which contains a node of element i
+      listElements.clear();
+
+      for(int j=indC[i-1];j<indC[i];j++){
+	// c element contains node i
+	int c=revC[j-1];
+	listElements.insert(c);
+      }
+
+      // calculate field value on node i 
+      sigmaV = 0.;
+      val = 0.;
+      for(elemIt=listElements.begin();elemIt!=listElements.end();elemIt++){
+	int elem = *elemIt;
+	double vol = volume->getValueIJ(elem,1);
+	if( vol != 0. ){
+	  sigmaV += 1./vol;
+	  if(_myDField)
+	    val += _myDField->getValueIJ(elem,1)/vol;
+	  else
+	    val += ((double)_myIField->getValueIJ(elem,1))/vol;
+	}
+      }
+      val /= sigmaV;
+      break;
+    case MED_FACE:
+      throw SALOME_Exception("Filter doesn't run on reference field on faces");
+      break;
+    case MED_EDGE:
+      throw SALOME_Exception("Filter doesn't run on reference field on edges");
+      break;
+    case MED_NODE:
+      // read reference field value
+      if(_myDField)
+	val = _myDField->getValueIJ(i,1);
+      else
+	val = (double)_myIField->getValueIJ(i,1);
+      break;
+    case MED_ALL_ENTITIES:
+      throw SALOME_Exception("Filter doesn't run on reference field on all entities");
+      break;
+    }
 
     // set criteria field value
     if( _myOneThresh->isChecked() ){
-- 
2.39.2