From 018b76dabfa8c314d12428f6b2c9a11b77be2724 Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Thu, 26 Jul 2001 21:26:35 +0000
Subject: [PATCH] Single MeanPlane Better tolerance in sys3x3_with_tol 1.22
---
Graphics/Geom.cpp | 218 ++++--------------------------------
Makefile | 4 +-
Mesh/2D_Mesh.cpp | 206 ++--------------------------------
Mesh/Numeric.cpp | 257 ++++++++++++++++++++++++++++++++++++++++++-
Mesh/Numeric.h | 2 +
doc/VERSIONS | 4 +-
doc/gmsh.1 | 4 +-
utils/gmsh_fltk.spec | 4 +-
www/gmsh.html | 20 ++--
9 files changed, 303 insertions(+), 416 deletions(-)
diff --git a/Graphics/Geom.cpp b/Graphics/Geom.cpp
index 08d189a2e8..a7dd9b9fdb 100644
--- a/Graphics/Geom.cpp
+++ b/Graphics/Geom.cpp
@@ -1,4 +1,4 @@
-// $Id: Geom.cpp,v 1.23 2001-07-18 07:36:36 geuzaine Exp $
+// $Id: Geom.cpp,v 1.24 2001-07-26 21:26:34 geuzaine Exp $
#include "Gmsh.h"
#include "GmshUI.h"
@@ -216,209 +216,36 @@ int isPointOnPlanarSurface (Surface *S, double X, double Y, double Z, double n[3
}
-void Plan_SurfPlane (void *data,void *dum){
- Surface **pS , *s;
- Curve *pC;
- Vertex *v;
- int ix,iy,iz,i,j,N;
- double det,sys[3][3],b[3],res[3],mod,t1[3],t2[3],ex[3],s2s[2][2],r2[2],X,Y,Z;
-
- static List_T *points;
- static int deb=1;
-
- pS = (Surface**)data;
- s = *pS;
-
- if(s->Typ != MSH_SURF_PLAN) return ;
-
- if(deb){
- points = List_Create(10,10,sizeof(Vertex*));
- deb = 0;
- }
- else
- List_Reset(points);
-
- for(i=0;i<List_Nbr(s->Generatrices);i++){
- List_Read(s->Generatrices,i,&pC);
- for(j=0;j<List_Nbr(pC->Control_Points);j++){
- List_Add(points,List_Pointer(pC->Control_Points,j));
- }
- }
-
- N = List_Nbr(points);
-
- for(i=0;i<3;i++){
- b[i] = 0.0;
- for(j=0;j<3;j++){
- sys[i][j] = 0.0;
- }
- }
+void Draw_Plane_Surface (Surface *s){
+ int i, j, k;
+ Curve *c;
+ double minx, miny, maxx, maxy, t, n[3], nn;
+ Vertex P1, P2, P3, V[4], vv, vv1, vv2;
+ char Num[100];
- /* ax + by + cz = 1 */
+ static List_T *points;
+ static int deb=1;
- ix=iy=iz=0;
+ if(!s->Orientations){
- // TOLERANCE ! WARNING WARNING
- double eps = 1.e-6 * CTX.lc;
+ s->Orientations = List_Create(20,2,sizeof(Vertex));
- for(i=0;i<N;i++){
- List_Read(points,i,&v);
- if(!i){
- X = v->Pos.X;
- Y = v->Pos.Y;
- Z = v->Pos.Z;
- }
- else{
- if(fabs(X-v->Pos.X) > eps) ix = 1;
- if(fabs(Y-v->Pos.Y) > eps) iy = 1;
- if(fabs(Z-v->Pos.Z) > eps) iz = 1;
+ if(deb){
+ points = List_Create(10,10,sizeof(Vertex*));
+ deb = 0;
}
+ else
+ List_Reset(points);
- sys[0][0] += v->Pos.X * v->Pos.X;
- sys[1][1] += v->Pos.Y * v->Pos.Y;
- sys[2][2] += v->Pos.Z * v->Pos.Z;
- sys[0][1] += v->Pos.X * v->Pos.Y;
- sys[0][2] += v->Pos.X * v->Pos.Z;
- sys[1][2] += v->Pos.Y * v->Pos.Z;
- sys[2][1] = sys[1][2];
- sys[1][0] = sys[0][1];
- sys[2][0] = sys[0][2];
- b[0] += v->Pos.X;
- b[1] += v->Pos.Y;
- b[2] += v->Pos.Z;
- }
-
- s->d = 1.0;
-
- /* x = X */
-
- if(!ix){
- s->d = X;
- res[0] = 1.;
- res[1] = res[2] = 0.0;
- }
-
- /* y = Y */
-
- else if(!iy){
- s->d = Y;
- res[1] = 1.;
- res[0] = res[2] = 0.0;
- }
-
- /* z = Z */
-
- else if(!iz){
- s->d = Z;
- res[2] = 1.;
- res[1] = res[0] = 0.0;
- }
-
- /* by + cz = -x */
-
- else if (!sys3x3_with_tol(sys,b,res,&det) ){
-
- s->d = 0.0;
- s2s[0][0] = sys[1][1];
- s2s[0][1] = sys[1][2];
- s2s[1][0] = sys[1][2];
- s2s[1][1] = sys[2][2];
- b[0] = -sys[0][1];
- b[1] = -sys[0][2];
- if(sys2x2(s2s,b,r2)){
- res[0] = 1.;
- res[1] = r2[0];
- res[2] = r2[1];
- }
- /* ax + cz = -y */
- else {
- s->d = 0.0;
- s2s[0][0] = sys[0][0];
- s2s[0][1] = sys[0][2];
- s2s[1][0] = sys[0][2];
- s2s[1][1] = sys[2][2];
- b[0] = -sys[0][1];
- b[1] = -sys[1][2];
- if(sys2x2(s2s,b,r2)){
- res[0] = r2[0];
- res[1] = 1.;
- res[2] = r2[1];
- }
- /* ax + by = -z */
- else {
- s->d = 1.0;
- s2s[0][0] = sys[0][0];
- s2s[0][1] = sys[0][1];
- s2s[1][0] = sys[0][1];
- s2s[1][1] = sys[1][1];
- b[0] = -sys[0][2];
- b[1] = -sys[1][2];
- if(sys2x2(s2s,b,r2)){
- res[0] = r2[0];
- res[1] = r2[1];
- res[2] = 1.;
- }
- else {
- Msg(GERROR, "Draw Geometry (Plan_SurfPlane)");
- }
+ for(i=0;i<List_Nbr(s->Generatrices);i++){
+ List_Read(s->Generatrices,i,&c);
+ for(j=0;j<List_Nbr(c->Control_Points);j++){
+ List_Add(points,List_Pointer(c->Control_Points,j));
}
}
- }
-
- s->a = res[0];
- s->b = res[1];
- s->c = res[2];
- mod = sqrt (res[0] * res[0] + res[1] * res[1] + res[2] * res[2]);
- for(i=0;i<3;i++) res[i]/=mod;
-
- /* L'axe n'est pas l'axe des x */
- if(res[0] > res[1]){
- ex[0] = 0.;
- ex[1] = 1.;
- ex[2] = 0.;
- }
- else{
- ex[0] = 1.;
- ex[1] = 0.;
- ex[2] = 0.;
- }
-
- prodve(res,ex,t1);
-
- mod = sqrt (t1[0] * t1[0] + t1[1] * t1[1] + t1[2] * t1[2] ) ;
- for(i=0;i<3;i++) t1[i]/=mod;
-
- prodve(t1,res,t2);
-
- mod = sqrt (t2[0] * t2[0] + t2[1] * t2[1] + t2[2] * t2[2] ) ;
- for(i=0;i<3;i++) t2[i]/=mod;
-
- for(i=0;i<3;i++)s->plan[0][i] = t1[i];
- for(i=0;i<3;i++)s->plan[1][i] = t2[i];
- for(i=0;i<3;i++)s->plan[2][i] = res[i];
-
- /* Matrice orthogonale */
-
- for(i=0;i<3;i++){
- for(j=0;j<3;j++){
- s->invplan[i][j] = s->plan[j][i];
- }
- }
-
-}
+ MeanPlane(points, s);
-void Draw_Plane_Surface (Surface *s){
- int i,k;
- Curve *c;
- double minx,miny,maxx,maxy,t,n[3],nn;
- Vertex P1,P2,P3,V[4],vv,vv1,vv2;
- char Num[100];
-
- if(!s->Orientations){
-
- s->Orientations = List_Create(20,2,sizeof(Vertex));
- Plan_SurfPlane(&s,NULL);
k = 0;
for(i=0;i<List_Nbr(s->Generatrices);i++){
@@ -682,8 +509,7 @@ void Draw_Surface (void *a, void *b){
glDisable(GL_LINE_STIPPLE);
Tree_Action(s->STL->Simplexes,Draw_Simplex_Surfaces);
}
- else if(s->Typ == MSH_SURF_DISCRETE)
- {
+ else if(s->Typ == MSH_SURF_DISCRETE){
glDisable(GL_LINE_STIPPLE);
Tree_Action(s->Simplexes,Draw_Simplex_Surfaces);
}
diff --git a/Makefile b/Makefile
index 77bb7d1cd0..e37a468774 100644
--- a/Makefile
+++ b/Makefile
@@ -1,9 +1,9 @@
-# $Id: Makefile,v 1.106 2001-07-26 18:47:59 remacle Exp $
+# $Id: Makefile,v 1.107 2001-07-26 21:26:34 geuzaine Exp $
# ----------------------------------------------------------------------
# Makefile for Gmsh
# ----------------------------------------------------------------------
- GMSH_RELEASE = 1.21
+ GMSH_RELEASE = 1.22
MAKE = make
CC = c++
diff --git a/Mesh/2D_Mesh.cpp b/Mesh/2D_Mesh.cpp
index dbe5d1042b..3b9d9655a9 100644
--- a/Mesh/2D_Mesh.cpp
+++ b/Mesh/2D_Mesh.cpp
@@ -1,4 +1,4 @@
-// $Id: 2D_Mesh.cpp,v 1.29 2001-06-25 18:34:59 remacle Exp $
+// $Id: 2D_Mesh.cpp,v 1.30 2001-07-26 21:26:34 geuzaine Exp $
/*
Maillage Delaunay d'une surface (Point insertion Technique)
@@ -55,14 +55,12 @@ void Projette_Inverse (void *a, void *b){
}
void Plan_Moyen (void *data, void *dum){
- int ix, iy, iz, i, j, N;
- static List_T *points;
- Curve *pC;
- static int deb = 1;
- double det, sys[3][3], b[3], res[3], mod, t1[3], t2[3], ex[3];
- double s2s[2][2], r2[2], X, Y, Z;
- Vertex *v;
- Surface **pS, *s;
+ int i, j;
+ Vertex *v;
+ Curve *pC;
+ Surface **pS, *s;
+ static List_T *points;
+ static int deb = 1;
pS = (Surface **) data;
s = *pS;
@@ -91,195 +89,7 @@ void Plan_Moyen (void *data, void *dum){
break;
}
- N = List_Nbr (points);
-
- for (i = 0; i < 3; i++){
- b[i] = 0.0;
- for (j = 0; j < 3; j++){
- sys[i][j] = 0.0;
- }
- }
-
- /* ax + by + cz = 1 */
-
- ix = iy = iz = 0;
-
- // TOLERANCE ! WARNING WARNING
- double eps = 1.e-6 * CTX.lc;
-
- for (i = 0; i < N; i++){
- List_Read (points, i, &v);
-
- if (!i){
- X = v->Pos.X;
- Y = v->Pos.Y;
- Z = v->Pos.Z;
- }
- else{
- if(fabs(X-v->Pos.X) > eps) ix = 1;
- if(fabs(Y-v->Pos.Y) > eps) iy = 1;
- if(fabs(Z-v->Pos.Z) > eps) iz = 1;
- }
-
- sys[0][0] += v->Pos.X * v->Pos.X;
- sys[1][1] += v->Pos.Y * v->Pos.Y;
- sys[2][2] += v->Pos.Z * v->Pos.Z;
- sys[0][1] += v->Pos.X * v->Pos.Y;
- sys[0][2] += v->Pos.X * v->Pos.Z;
- sys[1][2] += v->Pos.Y * v->Pos.Z;
- sys[2][1] = sys[1][2];
- sys[1][0] = sys[0][1];
- sys[2][0] = sys[0][2];
- b[0] += v->Pos.X;
- b[1] += v->Pos.Y;
- b[2] += v->Pos.Z;
- }
-
- s->d = 1.0;
-
- /* x = X */
-
- if (!ix){
- s->d = X;
- res[0] = 1.;
- res[1] = res[2] = 0.0;
- Msg(DEBUG, "Mean plane of type 'x = c'");
- }
-
- /* y = Y */
-
- else if (!iy){
- s->d = Y;
- res[1] = 1.;
- res[0] = res[2] = 0.0;
- Msg(DEBUG, "Mean plane of type 'y = c'");
- }
-
- /* z = Z */
-
- else if (!iz){
- s->d = Z;
- res[2] = 1.;
- res[1] = res[0] = 0.0;
- Msg(DEBUG, "Mean plane of type 'y = c'");
- }
-
- /* by + cz = -x */
-
- else if (!sys3x3_with_tol (sys, b, res, &det)){
- s->d = 0.0;
- s2s[0][0] = sys[1][1];
- s2s[0][1] = sys[1][2];
- s2s[1][0] = sys[1][2];
- s2s[1][1] = sys[2][2];
- b[0] = -sys[0][1];
- b[1] = -sys[0][2];
- if (sys2x2 (s2s, b, r2)){
- res[0] = 1.;
- res[1] = r2[0];
- res[2] = r2[1];
- Msg(DEBUG, "Mean plane of type 'by + cz = -x'");
- }
-
- /* ax + cz = -y */
-
- else{
- s->d = 0.0;
- s2s[0][0] = sys[0][0];
- s2s[0][1] = sys[0][2];
- s2s[1][0] = sys[0][2];
- s2s[1][1] = sys[2][2];
- b[0] = -sys[0][1];
- b[1] = -sys[1][2];
- if (sys2x2 (s2s, b, r2)){
- res[0] = r2[0];
- res[1] = 1.;
- res[2] = r2[1];
- Msg(DEBUG, "Mean plane of type 'ax + cz = -y'");
- }
-
- /* ax + by = -z */
-
- else{
- s->d = 1.0;
- s2s[0][0] = sys[0][0];
- s2s[0][1] = sys[0][1];
- s2s[1][0] = sys[0][1];
- s2s[1][1] = sys[1][1];
- b[0] = -sys[0][2];
- b[1] = -sys[1][2];
- if (sys2x2 (s2s, b, r2)){
- res[0] = r2[0];
- res[1] = r2[1];
- res[2] = 1.;
- Msg(DEBUG, "Mean plane of type 'ax + by = -z'");
- }
- else{
- Msg(GERROR, "Problem in mean plane computation");
- }
- }
- }
- }
-
- s->a = res[0];
- s->b = res[1];
- s->c = res[2];
- mod = sqrt (res[0] * res[0] + res[1] * res[1] + res[2] * res[2]);
- for (i = 0; i < 3; i++)
- res[i] /= mod;
-
- /* L'axe n'est pas l'axe des x */
-
- ex[0] = ex[1] = ex[2] = 0.0;
- if(res[0] == 0.0)
- ex[0] = 1.0;
- else if(res[1] == 0.0)
- ex[1] = 1.0;
- else
- ex[2] = 1.0;
-
- prodve (res, ex, t1);
-
- mod = sqrt (t1[0] * t1[0] + t1[1] * t1[1] + t1[2] * t1[2]);
- for (i = 0; i < 3; i++)
- t1[i] /= mod;
-
- prodve (t1, res, t2);
-
- mod = sqrt (t2[0] * t2[0] + t2[1] * t2[1] + t2[2] * t2[2]);
- for (i = 0; i < 3; i++)
- t2[i] /= mod;
-
- for (i = 0; i < 3; i++)
- s->plan[0][i] = t1[i];
- for (i = 0; i < 3; i++)
- s->plan[1][i] = t2[i];
- for (i = 0; i < 3; i++)
- s->plan[2][i] = res[i];
-
- Msg(DEBUG1, "Plane : (%g x + %g y + %g z = %g)", s->a, s->b, s->c, s->d);
- Msg(DEBUG2, "Normal : (%g , %g , %g )", res[0], res[1], res[2]);
- Msg(DEBUG2, "t1 : (%g , %g , %g )", t1[0], t1[1], t1[2]);
- Msg(DEBUG3, "t2 : (%g , %g , %g )", t2[0], t2[1], t2[2]);
-
- /* Matrice orthogonale */
-
- if (!iz){
- for (i = 0; i < 3; i++){
- for (j = 0; j < 3; j++){
- s->invplan[i][j] = (i == j) ? 1. : 0.;
- s->plan[i][j] = (i == j) ? 1. : 0.;
- }
- }
- }
- else{
- for (i = 0; i < 3; i++){
- for (j = 0; j < 3; j++){
- s->invplan[i][j] = s->plan[j][i];
- }
- }
- }
-
+ MeanPlane(points, s);
}
diff --git a/Mesh/Numeric.cpp b/Mesh/Numeric.cpp
index 7772de5ffe..7e89139c7f 100644
--- a/Mesh/Numeric.cpp
+++ b/Mesh/Numeric.cpp
@@ -1,4 +1,4 @@
-// $Id: Numeric.cpp,v 1.15 2001-06-02 13:09:14 geuzaine Exp $
+// $Id: Numeric.cpp,v 1.16 2001-07-26 21:26:34 geuzaine Exp $
#include "Gmsh.h"
#include "Const.h"
@@ -8,6 +8,9 @@
#include "Numeric.h"
#include "Interpolation.h"
#include "nrutil.h"
+#include "Context.h"
+
+extern Context_T CTX;
double myatan2 (double a, double b){
if (a == 0.0 && b == 0)
@@ -67,7 +70,9 @@ int sys2x2 (double mat[2][2], double b[2], double res[2]){
norm = DSQR (mat[0][0]) + DSQR (mat[1][1]) + DSQR (mat[0][1]) + DSQR (mat[1][0]);
det = mat[0][0] * mat[1][1] - mat[1][0] * mat[0][1];
+ // TOLERANCE ! WARNING WARNING
if (norm == 0.0 || fabs (det) / norm < 1.e-07){
+ Msg(DEBUG, "Assuming 2x2 matrix is singular (det/norm == %g)", fabs(det)/norm);
res[0] = res[1] = 0.0 ;
return 0;
}
@@ -116,13 +121,19 @@ int sys3x3 (double mat[3][3], double b[3], double res[3], double *det){
int sys3x3_with_tol (double mat[3][3], double b[3], double res[3], double *det){
int out;
+ double norm;
out = sys3x3(mat,b,res,det);
-
- if (fabs(*det) < 1.e-14){
- Msg(DEBUG, "Assuming 3x3 matrix is singular (det == %g)", fabs(*det));
+ norm =
+ DSQR (mat[0][0]) + DSQR (mat[0][1]) + DSQR (mat[0][2]) +
+ DSQR (mat[1][0]) + DSQR (mat[1][1]) + DSQR (mat[1][2]) +
+ DSQR (mat[2][0]) + DSQR (mat[2][1]) + DSQR (mat[2][2]) ;
+
+ // TOLERANCE ! WARNING WARNING
+ if (norm == 0.0 || fabs (*det) / norm < 1.e-07){
+ Msg(DEBUG, "Assuming 3x3 matrix is singular (det/norm == %g)", fabs(*det)/norm);
res[0] = res[1] = res[2] = 0.0 ;
- return (0);
+ return 0;
}
return out ;
@@ -161,6 +172,242 @@ int inv3x3 (double mat[3][3], double inv[3][3], double *det){
}
+void MeanPlane(List_T *points, Surface *s){
+ int i, j, ix, iy, iz, N;
+ double det,sys[3][3],b[3],res[3],mod,t1[3],t2[3],ex[3],s2s[2][2],r2[2],X,Y,Z;
+ Vertex *v;
+
+ N = List_Nbr (points);
+
+ for (i = 0; i < 3; i++){
+ b[i] = 0.0;
+ for (j = 0; j < 3; j++){
+ sys[i][j] = 0.0;
+ }
+ }
+
+ /* ax + by + cz = 1 */
+
+ ix = iy = iz = 0;
+
+ // TOLERANCE ! WARNING WARNING
+ double eps = 1.e-6 * CTX.lc;
+
+ for (i = 0; i < N; i++){
+ List_Read (points, i, &v);
+
+ if (!i){
+ X = v->Pos.X;
+ Y = v->Pos.Y;
+ Z = v->Pos.Z;
+ }
+ else{
+ if(fabs(X-v->Pos.X) > eps) ix = 1;
+ if(fabs(Y-v->Pos.Y) > eps) iy = 1;
+ if(fabs(Z-v->Pos.Z) > eps) iz = 1;
+ }
+
+ sys[0][0] += v->Pos.X * v->Pos.X;
+ sys[1][1] += v->Pos.Y * v->Pos.Y;
+ sys[2][2] += v->Pos.Z * v->Pos.Z;
+ sys[0][1] += v->Pos.X * v->Pos.Y;
+ sys[0][2] += v->Pos.X * v->Pos.Z;
+ sys[1][2] += v->Pos.Y * v->Pos.Z;
+ sys[2][1] = sys[1][2];
+ sys[1][0] = sys[0][1];
+ sys[2][0] = sys[0][2];
+ b[0] += v->Pos.X;
+ b[1] += v->Pos.Y;
+ b[2] += v->Pos.Z;
+ }
+
+ s->d = 1.0;
+
+ /* x = X */
+
+ if (!ix){
+ s->d = X;
+ res[0] = 1.;
+ res[1] = res[2] = 0.0;
+ Msg(DEBUG, "Mean plane of type 'x = c'");
+ }
+
+ /* y = Y */
+
+ else if (!iy){
+ s->d = Y;
+ res[1] = 1.;
+ res[0] = res[2] = 0.0;
+ Msg(DEBUG, "Mean plane of type 'y = c'");
+ }
+
+ /* z = Z */
+
+ else if (!iz){
+ s->d = Z;
+ res[2] = 1.;
+ res[1] = res[0] = 0.0;
+ Msg(DEBUG, "Mean plane of type 'z = c'");
+ }
+
+ /* by + cz = -x */
+
+ else if (!sys3x3_with_tol (sys, b, res, &det)){
+ s->d = 0.0;
+ s2s[0][0] = sys[1][1];
+ s2s[0][1] = sys[1][2];
+ s2s[1][0] = sys[1][2];
+ s2s[1][1] = sys[2][2];
+ b[0] = -sys[0][1];
+ b[1] = -sys[0][2];
+ if (sys2x2 (s2s, b, r2)){
+ res[0] = 1.;
+ res[1] = r2[0];
+ res[2] = r2[1];
+ Msg(DEBUG, "Mean plane of type 'by + cz = -x'");
+ }
+
+ /* ax + cz = -y */
+
+ else{
+ s->d = 0.0;
+ s2s[0][0] = sys[0][0];
+ s2s[0][1] = sys[0][2];
+ s2s[1][0] = sys[0][2];
+ s2s[1][1] = sys[2][2];
+ b[0] = -sys[0][1];
+ b[1] = -sys[1][2];
+ if (sys2x2 (s2s, b, r2)){
+ res[0] = r2[0];
+ res[1] = 1.;
+ res[2] = r2[1];
+ Msg(DEBUG, "Mean plane of type 'ax + cz = -y'");
+ }
+
+ /* ax + by = -z */
+
+ else{
+ s->d = 1.0;
+ s2s[0][0] = sys[0][0];
+ s2s[0][1] = sys[0][1];
+ s2s[1][0] = sys[0][1];
+ s2s[1][1] = sys[1][1];
+ b[0] = -sys[0][2];
+ b[1] = -sys[1][2];
+ if (sys2x2 (s2s, b, r2)){
+ res[0] = r2[0];
+ res[1] = r2[1];
+ res[2] = 1.;
+ Msg(DEBUG, "Mean plane of type 'ax + by = -z'");
+ }
+ else{
+ Msg(GERROR, "Problem in mean plane computation");
+ }
+ }
+ }
+ }
+
+ s->a = res[0];
+ s->b = res[1];
+ s->c = res[2];
+ mod = sqrt (res[0] * res[0] + res[1] * res[1] + res[2] * res[2]);
+ for (i = 0; i < 3; i++)
+ res[i] /= mod;
+
+ /* L'axe n'est pas l'axe des x */
+
+ ex[0] = ex[1] = ex[2] = 0.0;
+ if(res[0] == 0.0)
+ ex[0] = 1.0;
+ else if(res[1] == 0.0)
+ ex[1] = 1.0;
+ else
+ ex[2] = 1.0;
+
+ prodve (res, ex, t1);
+
+ mod = sqrt (t1[0] * t1[0] + t1[1] * t1[1] + t1[2] * t1[2]);
+ for (i = 0; i < 3; i++)
+ t1[i] /= mod;
+
+ prodve (t1, res, t2);
+
+ mod = sqrt (t2[0] * t2[0] + t2[1] * t2[1] + t2[2] * t2[2]);
+ for (i = 0; i < 3; i++)
+ t2[i] /= mod;
+
+ for (i = 0; i < 3; i++)
+ s->plan[0][i] = t1[i];
+ for (i = 0; i < 3; i++)
+ s->plan[1][i] = t2[i];
+ for (i = 0; i < 3; i++)
+ s->plan[2][i] = res[i];
+
+ Msg(DEBUG1, "Plane : (%g x + %g y + %g z = %g)", s->a, s->b, s->c, s->d);
+ Msg(DEBUG2, "Normal : (%g , %g , %g )", res[0], res[1], res[2]);
+ Msg(DEBUG2, "t1 : (%g , %g , %g )", t1[0], t1[1], t1[2]);
+ Msg(DEBUG3, "t2 : (%g , %g , %g )", t2[0], t2[1], t2[2]);
+
+ /* Matrice orthogonale */
+
+ if (!iz){
+ for (i = 0; i < 3; i++){
+ for (j = 0; j < 3; j++){
+ s->invplan[i][j] = (i == j) ? 1. : 0.;
+ s->plan[i][j] = (i == j) ? 1. : 0.;
+ }
+ }
+ }
+ else{
+ for (i = 0; i < 3; i++){
+ for (j = 0; j < 3; j++){
+ s->invplan[i][j] = s->plan[j][i];
+ }
+ }
+ }
+
+
+// this is the end of the algo as it was used for surface drawing:
+
+#if 0
+ /* L'axe n'est pas l'axe des x */
+ if(res[0] > res[1]){
+ ex[0] = 0.;
+ ex[1] = 1.;
+ ex[2] = 0.;
+ }
+ else{
+ ex[0] = 1.;
+ ex[1] = 0.;
+ ex[2] = 0.;
+ }
+
+ prodve(res,ex,t1);
+
+ mod = sqrt (t1[0] * t1[0] + t1[1] * t1[1] + t1[2] * t1[2] ) ;
+ for(i=0;i<3;i++) t1[i]/=mod;
+
+ prodve(t1,res,t2);
+
+ mod = sqrt (t2[0] * t2[0] + t2[1] * t2[1] + t2[2] * t2[2] ) ;
+ for(i=0;i<3;i++) t2[i]/=mod;
+
+ for(i=0;i<3;i++)s->plan[0][i] = t1[i];
+ for(i=0;i<3;i++)s->plan[1][i] = t2[i];
+ for(i=0;i<3;i++)s->plan[2][i] = res[i];
+
+ /* Matrice orthogonale */
+
+ for(i=0;i<3;i++){
+ for(j=0;j<3;j++){
+ s->invplan[i][j] = s->plan[j][i];
+ }
+ }
+#endif
+}
+
+
+
#define Precision 1.e-10
#define MaxIter 20
diff --git a/Mesh/Numeric.h b/Mesh/Numeric.h
index 7d7a05aad4..4208df365c 100644
--- a/Mesh/Numeric.h
+++ b/Mesh/Numeric.h
@@ -13,6 +13,8 @@ int sys3x3 (double mat[3][3], double b[3], double res[3], double *det);
int sys3x3_with_tol (double mat[3][3], double b[3], double res[3], double *det);
int det3x3 (double mat[3][3], double *det);
int inv3x3 (double mat[3][3], double inv[3][3], double *det);
+void MeanPlane(List_T *point, Surface *s);
+
void find_bestuv (Surface * s, double X, double Y,
double *U, double *V, double *Z, int N);
void XYtoUV (Surface * s, double *X, double *Y,
diff --git a/doc/VERSIONS b/doc/VERSIONS
index 1ffa8141ef..bdfbd58aa1 100644
--- a/doc/VERSIONS
+++ b/doc/VERSIONS
@@ -1,4 +1,6 @@
-$Id: VERSIONS,v 1.32 2001-07-26 15:25:36 geuzaine Exp $
+$Id: VERSIONS,v 1.33 2001-07-26 21:26:35 geuzaine Exp $
+
+New in 1.22: Fixed (yet another) bug for 2D mesh in the mean plane;
New in 1.21: Fixed more memory leaks; added -opt command line option
to parse definitions directly from the command line; fixed missing
diff --git a/doc/gmsh.1 b/doc/gmsh.1
index b3c2d96275..de93e21487 100644
--- a/doc/gmsh.1
+++ b/doc/gmsh.1
@@ -5,7 +5,7 @@
.\" Copyright (c) 2000-2001 J.-F. Remacle, C. Geuzaine
.\"
.\" ======================================================================
-.TH Gmsh 1 "15 February 2001" "Version 1.21" "Gmsh Manual Pages"
+.TH Gmsh 1 "15 February 2001" "Version 1.22" "Gmsh Manual Pages"
.UC 4
.\" ======================================================================
.SH NAME
@@ -220,6 +220,6 @@ Remacle (Remacle@scorec.rpi.edu).
.SH SEE ALSO
.BR getdp (1),
.br
-Gmsh examples (\fI/usr/doc/gmsh-1.21/\fR),
+Gmsh examples (\fI/usr/doc/gmsh-1.22/\fR),
.br
Gmsh homepage (\fIhttp://www.geuz.org/gmsh/\fR).
diff --git a/utils/gmsh_fltk.spec b/utils/gmsh_fltk.spec
index 745eb3b537..4eb52ba57c 100644
--- a/utils/gmsh_fltk.spec
+++ b/utils/gmsh_fltk.spec
@@ -1,7 +1,7 @@
Summary: A 3D mesh generator with pre- and post-processing facilities
Name: gmsh
-Version: 1.21
-Source: gmsh-1.21.tar.gz
+Version: 1.22
+Source: gmsh-1.22.tar.gz
Release: 1
Copyright: distributable
Group: Applications/Engineering
diff --git a/www/gmsh.html b/www/gmsh.html
index a3eb14088b..ab54aa60cf 100644
--- a/www/gmsh.html
+++ b/www/gmsh.html
@@ -46,7 +46,7 @@ cat << EOM
page requests since<br>1998/05/24<p>
ENDSCRIPT--->
-<!---BEGINDATE$Date: 2001-07-25 17:53:55 $ENDDATE--->
+<!---BEGINDATE$Date: 2001-07-26 21:26:35 $ENDDATE--->
Copyright © 1998-2001<br>
Jean-François Remacle and
@@ -324,7 +324,7 @@ ENDSCRIPT--->
<td bgcolor="#ededed"><font face="Helvetica, Arial" size=-1>
-<b>Latest Release: 1.21 (July 26, 2001)</b>
+<b>Latest Release: 1.22 (August 10, 2001)</b>
<p>
Executable versions of Gmsh are available for Windows and for most of
the classical UNIX platforms. These versions are free, and are all
@@ -333,14 +333,14 @@ name="opengl-footmark"><sup>1</sup></a>. The only thing required if
you use Gmsh is to mention it in your work. The tutorial and demo
files are included in the archives.
<ul>
-<li><A href="/gmsh/bin/gmsh-1.21-Windows.zip">Windows zip archive (95/98/NT)</A>
-<li><A href="/gmsh/bin/gmsh-1.21-1.i386.rpm">Linux RPM (Red Hat 6.2 and compatible, i386, glibc 2.1)</A>
-<li><A href="/gmsh/bin/gmsh-1.21-Linux.tgz">Linux tarball (i386, glibc 2.1)</A>
-<li><A href="/gmsh/bin/gmsh-1.21-OSF1.tgz">Compaq Tru64 tarball (OSF 4.0)</A>
-<li><A href="/gmsh/bin/gmsh-1.21-SunOS.tgz">Sun tarball (SunOS 5.5)</A>
-<li><A href="/gmsh/bin/gmsh-1.21-AIX.tgz">IBM tarball (AIX)</A>
-<li><A href="/gmsh/bin/gmsh-1.21-IRIX.tgz">SGI IRIX tarball (IRIX 6.5)</A>
-<li><A href="/gmsh/bin/gmsh-1.21-HP-UX.tgz">HP tarball (HPUX 10.20)</A>
+<li><A href="/gmsh/bin/gmsh-1.22-Windows.zip">Windows zip archive (95/98/NT)</A>
+<li><A href="/gmsh/bin/gmsh-1.22-1.i386.rpm">Linux RPM (Red Hat 6.2 and compatible, i386, glibc 2.1)</A>
+<li><A href="/gmsh/bin/gmsh-1.22-Linux.tgz">Linux tarball (i386, glibc 2.1)</A>
+<li><A href="/gmsh/bin/gmsh-1.22-OSF1.tgz">Compaq Tru64 tarball (OSF 4.0)</A>
+<li><A href="/gmsh/bin/gmsh-1.22-SunOS.tgz">Sun tarball (SunOS 5.5)</A>
+<li><A href="/gmsh/bin/gmsh-1.22-AIX.tgz">IBM tarball (AIX)</A>
+<li><A href="/gmsh/bin/gmsh-1.22-IRIX.tgz">SGI IRIX tarball (IRIX 6.5)</A>
+<li><A href="/gmsh/bin/gmsh-1.22-HP-UX.tgz">HP tarball (HPUX 10.20)</A>
</ul>
--
GitLab