From a9b14a1d20ce46a4afadbe3598e6545eee2bc01f Mon Sep 17 00:00:00 2001
From: Jean-Francois Remacle <jean-francois.remacle@uclouvain.be>
Date: Mon, 22 Nov 2004 11:35:26 +0000
Subject: [PATCH] *** empty log message ***
---
Common/AdaptiveViews.cpp | 30 ++++++++++++++++++++----------
Plugin/Integrate.cpp | 29 +++++++++++++++++++++++++----
2 files changed, 45 insertions(+), 14 deletions(-)
diff --git a/Common/AdaptiveViews.cpp b/Common/AdaptiveViews.cpp
index 9401c08435..15a3022434 100644
--- a/Common/AdaptiveViews.cpp
+++ b/Common/AdaptiveViews.cpp
@@ -80,6 +80,19 @@ void _triangle::Recur_Create (_triangle *t, int maxlevel, int level , Double_Mat
all_triangles.push_back(t);
if (level++ >= maxlevel)
return;
+
+ /*
+
+ p3
+
+
+ p13 p23
+
+
+ p1 p12 p2
+
+
+ */
_point *p1 = t->p[0];
_point *p2 = t->p[1];
@@ -87,13 +100,13 @@ void _triangle::Recur_Create (_triangle *t, int maxlevel, int level , Double_Mat
_point *p12 = _point::New ( (p1->x+p2->x)*0.5,(p1->y+p2->y)*0.5,0, coeffs, eexps);
_point *p13 = _point::New ( (p1->x+p3->x)*0.5,(p1->y+p3->y)*0.5,0, coeffs, eexps);
_point *p23 = _point::New ( (p3->x+p2->x)*0.5,(p3->y+p2->y)*0.5,0, coeffs, eexps);
- _triangle *t1 = new _triangle (p1,p12,p13);
+ _triangle *t1 = new _triangle (p1,p13,p12);
Recur_Create (t1, maxlevel,level,coeffs,eexps);
_triangle *t2 = new _triangle (p12,p23,p2);
Recur_Create (t2, maxlevel,level,coeffs,eexps);
_triangle *t3 = new _triangle (p23,p13,p3);
Recur_Create (t3, maxlevel,level,coeffs,eexps);
- _triangle *t4 = new _triangle (p12,p23,p13);
+ _triangle *t4 = new _triangle (p12,p13,p23);
Recur_Create (t4, maxlevel,level,coeffs,eexps);
t->t[0]=t1;t->t[1]=t2;t->t[2]=t3;t->t[3]=t4;
@@ -215,18 +228,15 @@ void Adaptive_Post_View:: zoomElement (Post_View * view ,
{
_point *p = (_point*) &(*it);
p->val = res(kk);
- // for ( int k=0;k<N;++k)
- // {
- // p->val += p->shape_functions[k] * (*_STval )( ielem , k );
- // }
- // p->X = (*_STposX) ( ielem , 0 ) * ( 1.-p->x-p->y) + (*_STposX) ( ielem , 1 ) * p->x + (*_STposX) ( ielem , 2 ) * p->y;
- // p->Y = (*_STposY) ( ielem , 0 ) * ( 1.-p->x-p->y) + (*_STposY) ( ielem , 1 ) * p->x + (*_STposY) ( ielem , 2 ) * p->y;
- // p->Z = (*_STposZ) ( ielem , 0 ) * ( 1.-p->x-p->y) + (*_STposZ) ( ielem , 1 ) * p->x + (*_STposZ) ( ielem , 2 ) * p->y;
-
p->X = XYZ(kk,0);
p->Y = XYZ(kk,1);
p->Z = XYZ(kk,2);
+ /// -----------------------------------------------
+ /// TEST TEST ZZUT CHIOTTE VIREZ MOI CA VITE FAIT
+ /// -----------------------------------------------
+ // p->val = p->X * p->X + p->Y*p->Y - 1;
+
if (min > p->val) min = p->val;
if (max < p->val) max = p->val;
kk++;
diff --git a/Plugin/Integrate.cpp b/Plugin/Integrate.cpp
index f9f8112dd4..52485c211c 100644
--- a/Plugin/Integrate.cpp
+++ b/Plugin/Integrate.cpp
@@ -1,4 +1,4 @@
-// $Id: Integrate.cpp,v 1.1 2004-11-13 23:57:28 geuzaine Exp $
+// $Id: Integrate.cpp,v 1.2 2004-11-22 11:35:26 remacle Exp $
//
// Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
//
@@ -30,7 +30,8 @@
extern Context_T CTX;
StringXNumber IntegrateOptions_Number[] = {
- {GMSH_FULLRC, "iView", NULL, -1.}
+ {GMSH_FULLRC, "iView", NULL, -1.},
+ {GMSH_FULLRC, "computeLevelsetPositive", NULL, 0.}
};
extern "C"
@@ -63,7 +64,11 @@ void GMSH_IntegratePlugin::getInfos(char *author, char *copyright,
"line/surface elements. If `iView' < 0, the plugin\n"
"is run on the current view.\n"
"\n"
- "Plugin(Integrate) creates one new view.\n");
+ "If computeLevelsetPositive is not 0, then the plugin \n"
+ "computes the positive area (volume)of the map\n"
+ "\n"
+ "Plugin(Integrate) creates one new view.\n"
+);
}
int GMSH_IntegratePlugin::getNbOptions() const
@@ -86,6 +91,8 @@ static double integrate(int nbList, List_T *list, int dim,
{
if(!nbList) return 0.;
+ const int levelsetPositive = (int)IntegrateOptions_Number[1].def;
+
double res = 0.;
int nb = List_Nbr(list) / nbList;
for(int i = 0; i < List_Nbr(list); i += nb) {
@@ -110,7 +117,20 @@ static double integrate(int nbList, List_T *list, int dim,
else if(dim == 2){
if(nbNod == 3){
triangle t(x, y, z);
- if(nbComp == 1) res += t.integrate(v);
+ if(nbComp == 1)
+ {
+ if ( ! levelsetPositive )
+ res += t.integrate(v);
+ else
+ {
+ double ONES[] = { 1.0 , 1.0 , 1.0 };
+ const double area = t.integrate (ONES);
+ const double SUM = v[0] + v[1] + v[2];
+ const double SUMABS = fabs(v[0]) + fabs(v[1]) + fabs (v[2]);
+ const double XI = SUM / SUMABS;
+ res += area * (1 - XI) * 0.5 ;
+ }
+ }
else if(nbComp == 3) res += t.integrateFlux(v);
}
else if(nbNod == 4){
@@ -138,6 +158,7 @@ static double integrate(int nbList, List_T *list, int dim,
}
}
}
+ printf("integration res = %22.15E\n",res);
return res;
}
--
GitLab