Skip to content
Snippets Groups Projects
Commit f55145d2 authored by Christophe Geuzaine's avatar Christophe Geuzaine
Browse files

forgot to indent these two
parent d3aba987
No related branches found
No related tags found
No related merge requests found
......@@ -18,9 +18,6 @@
//
// Please report all bugs and problems to <gmsh@geuz.org>.
// Don't indent this file
// *INDENT-OFF*
#include <stdio.h>
#include <math.h>
#include <list>
......@@ -30,7 +27,8 @@
// A recursive effective implementation
void computeShapeFunctions ( Double_Matrix *coeffs, Double_Matrix *eexps , double u, double v, double w,double *sf);
void computeShapeFunctions(Double_Matrix * coeffs, Double_Matrix * eexps,
double u, double v, double w, double *sf);
std::set < adapt_point > adapt_point::all_points;
std::list < adapt_edge * >adapt_edge::all_elems;
......@@ -45,13 +43,15 @@ int adapt_tet::nbNod = 4;
int adapt_quad::nbNod = 4;
int adapt_hex::nbNod = 8;
adapt_point * adapt_point::New ( double x, double y, double z, Double_Matrix *coeffs, Double_Matrix *eexps)
adapt_point *adapt_point::New(double x, double y, double z,
Double_Matrix * coeffs, Double_Matrix * eexps)
{
adapt_point p;
p.x=x; p.y=y; p.z=z;
p.x = x;
p.y = y;
p.z = z;
std::set < adapt_point >::iterator it = all_points.find(p);
if ( it == all_points.end() )
{
if(it == all_points.end()) {
all_points.insert(p);
it = all_points.find(p);
double *kkk = (double *)(it->shape_functions);
......@@ -62,7 +62,8 @@ adapt_point * adapt_point::New ( double x, double y, double z, Double_Matrix *co
return (adapt_point *) & (*it);
}
void adapt_edge::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_edge::Create(int maxlevel, Double_Matrix * coeffs,
Double_Matrix * eexps)
{
int level = 0;
cleanElement < adapt_edge > ();
......@@ -72,8 +73,8 @@ void adapt_edge::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eex
Recur_Create(t, maxlevel, level, coeffs, eexps);
}
void adapt_triangle::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_triangle::Create(int maxlevel, Double_Matrix * coeffs,
Double_Matrix * eexps)
{
int level = 0;
cleanElement < adapt_triangle > ();
......@@ -84,7 +85,8 @@ void adapt_triangle::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix
Recur_Create(t, maxlevel, level, coeffs, eexps);
}
void adapt_quad::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_quad::Create(int maxlevel, Double_Matrix * coeffs,
Double_Matrix * eexps)
{
int level = 0;
cleanElement < adapt_quad > ();
......@@ -96,7 +98,8 @@ void adapt_quad::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eex
Recur_Create(q, maxlevel, level, coeffs, eexps);
}
void adapt_tet::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_tet::Create(int maxlevel, Double_Matrix * coeffs,
Double_Matrix * eexps)
{
int level = 0;
cleanElement < adapt_tet > ();
......@@ -108,7 +111,8 @@ void adapt_tet::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexp
Recur_Create(t, maxlevel, level, coeffs, eexps);
}
void adapt_hex::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_hex::Create(int maxlevel, Double_Matrix * coeffs,
Double_Matrix * eexps)
{
int level = 0;
cleanElement < adapt_hex > ();
......@@ -124,8 +128,8 @@ void adapt_hex::Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexp
Recur_Create(h, maxlevel, level, coeffs, eexps);
}
void adapt_edge::Recur_Create (adapt_edge *e, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_edge::Recur_Create(adapt_edge * e, int maxlevel, int level,
Double_Matrix * coeffs, Double_Matrix * eexps)
{
all_elems.push_back(e);
if(level++ >= maxlevel)
......@@ -137,40 +141,45 @@ void adapt_edge::Recur_Create (adapt_edge *e, int maxlevel, int level , Double_M
adapt_point *p1 = e->p[0];
adapt_point *p2 = e->p[1];
adapt_point *p12 = adapt_point::New ( (p1->x+p2->x)*0.5,(p1->y+p2->y)*0.5, (p1->z+p2->z)*0.5, coeffs, eexps);
adapt_point *p12 =
adapt_point::New((p1->x + p2->x) * 0.5, (p1->y + p2->y) * 0.5,
(p1->z + p2->z) * 0.5, coeffs, eexps);
adapt_edge *e1 = new adapt_edge(p1, p12);
Recur_Create(e1, maxlevel, level, coeffs, eexps);
adapt_edge *e2 = new adapt_edge(p12, p2);
Recur_Create(e2, maxlevel, level, coeffs, eexps);
e->e[0]=e1;e->e[1]=e2;
e->e[0] = e1;
e->e[1] = e2;
}
void adapt_triangle::Recur_Create (adapt_triangle *t, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_triangle::Recur_Create(adapt_triangle * t, int maxlevel, int level,
Double_Matrix * coeffs,
Double_Matrix * eexps)
{
all_elems.push_back(t);
if(level++ >= maxlevel)
return;
/*
p3
p13 p23
p1 p12 p2
*/
adapt_point *p1 = t->p[0];
adapt_point *p2 = t->p[1];
adapt_point *p3 = t->p[2];
adapt_point *p12 = adapt_point::New ( (p1->x+p2->x)*0.5,(p1->y+p2->y)*0.5,0, coeffs, eexps);
adapt_point *p13 = adapt_point::New ( (p1->x+p3->x)*0.5,(p1->y+p3->y)*0.5,0, coeffs, eexps);
adapt_point *p23 = adapt_point::New ( (p3->x+p2->x)*0.5,(p3->y+p2->y)*0.5,0, coeffs, eexps);
adapt_point *p12 =
adapt_point::New((p1->x + p2->x) * 0.5, (p1->y + p2->y) * 0.5, 0, coeffs,
eexps);
adapt_point *p13 =
adapt_point::New((p1->x + p3->x) * 0.5, (p1->y + p3->y) * 0.5, 0, coeffs,
eexps);
adapt_point *p23 =
adapt_point::New((p3->x + p2->x) * 0.5, (p3->y + p2->y) * 0.5, 0, coeffs,
eexps);
adapt_triangle *t1 = new adapt_triangle(p1, p13, p12);
Recur_Create(t1, maxlevel, level, coeffs, eexps);
adapt_triangle *t2 = new adapt_triangle(p12, p23, p2);
......@@ -179,38 +188,47 @@ void adapt_triangle::Recur_Create (adapt_triangle *t, int maxlevel, int level ,
Recur_Create(t3, maxlevel, level, coeffs, eexps);
adapt_triangle *t4 = new adapt_triangle(p12, p13, p23);
Recur_Create(t4, maxlevel, level, coeffs, eexps);
t->e[0]=t1;t->e[1]=t2;t->e[2]=t3;t->e[3]=t4;
t->e[0] = t1;
t->e[1] = t2;
t->e[2] = t3;
t->e[3] = t4;
}
void adapt_quad::Recur_Create (adapt_quad *q, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_quad::Recur_Create(adapt_quad * q, int maxlevel, int level,
Double_Matrix * coeffs, Double_Matrix * eexps)
{
all_elems.push_back(q);
if(level++ >= maxlevel)
return;
/*
p2 p23 p3
p12 pc p34
p1 p14 p4
*/
adapt_point *p1 = q->p[0];
adapt_point *p2 = q->p[1];
adapt_point *p3 = q->p[2];
adapt_point *p4 = q->p[3];
adapt_point *p12 = adapt_point::New ( (p1->x+p2->x)*0.5,(p1->y+p2->y)*0.5,0, coeffs, eexps);
adapt_point *p23 = adapt_point::New ( (p2->x+p3->x)*0.5,(p2->y+p3->y)*0.5,0, coeffs, eexps);
adapt_point *p34 = adapt_point::New ( (p4->x+p3->x)*0.5,(p4->y+p3->y)*0.5,0, coeffs, eexps);
adapt_point *p14 = adapt_point::New ( (p1->x+p4->x)*0.5,(p1->y+p4->y)*0.5,0, coeffs, eexps);
adapt_point *pc = adapt_point::New ( (p1->x+p2->x+p3->x+p4->x)*0.25,(p1->y+p2->y+p3->y+p4->y)*0.25,0, coeffs, eexps);
adapt_point *p12 =
adapt_point::New((p1->x + p2->x) * 0.5, (p1->y + p2->y) * 0.5, 0, coeffs,
eexps);
adapt_point *p23 =
adapt_point::New((p2->x + p3->x) * 0.5, (p2->y + p3->y) * 0.5, 0, coeffs,
eexps);
adapt_point *p34 =
adapt_point::New((p4->x + p3->x) * 0.5, (p4->y + p3->y) * 0.5, 0, coeffs,
eexps);
adapt_point *p14 =
adapt_point::New((p1->x + p4->x) * 0.5, (p1->y + p4->y) * 0.5, 0, coeffs,
eexps);
adapt_point *pc =
adapt_point::New((p1->x + p2->x + p3->x + p4->x) * 0.25,
(p1->y + p2->y + p3->y + p4->y) * 0.25, 0, coeffs,
eexps);
adapt_quad *q1 = new adapt_quad(p1, p12, pc, p14);
Recur_Create(q1, maxlevel, level, coeffs, eexps);
adapt_quad *q2 = new adapt_quad(p12, p2, p23, pc);
......@@ -219,11 +237,14 @@ void adapt_quad::Recur_Create (adapt_quad *q, int maxlevel, int level , Double_M
Recur_Create(q3, maxlevel, level, coeffs, eexps);
adapt_quad *q4 = new adapt_quad(p14, pc, p34, p4);
Recur_Create(q4, maxlevel, level, coeffs, eexps);
q->e[0]=q1;q->e[1]=q2;q->e[2]=q3;q->e[3]=q4;
q->e[0] = q1;
q->e[1] = q2;
q->e[2] = q3;
q->e[3] = q4;
}
void adapt_tet::Recur_Create (adapt_tet *t, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_tet::Recur_Create(adapt_tet * t, int maxlevel, int level,
Double_Matrix * coeffs, Double_Matrix * eexps)
{
all_elems.push_back(t);
if(level++ >= maxlevel)
......@@ -233,12 +254,24 @@ void adapt_tet::Recur_Create (adapt_tet *t, int maxlevel, int level , Double_Mat
adapt_point *p1 = t->p[1];
adapt_point *p2 = t->p[2];
adapt_point *p3 = t->p[3];
adapt_point *pe0 = adapt_point::New ( (p0->x+p1->x)*0.5,(p0->y+p1->y)*0.5,(p0->z+p1->z)*0.5,coeffs, eexps);
adapt_point *pe1 = adapt_point::New ( (p0->x+p2->x)*0.5,(p0->y+p2->y)*0.5,(p0->z+p2->z)*0.5,coeffs, eexps);
adapt_point *pe2 = adapt_point::New ( (p0->x+p3->x)*0.5,(p0->y+p3->y)*0.5,(p0->z+p3->z)*0.5,coeffs, eexps);
adapt_point *pe3 = adapt_point::New ( (p1->x+p2->x)*0.5,(p1->y+p2->y)*0.5,(p1->z+p2->z)*0.5,coeffs, eexps);
adapt_point *pe4 = adapt_point::New ( (p1->x+p3->x)*0.5,(p1->y+p3->y)*0.5,(p1->z+p3->z)*0.5,coeffs, eexps);
adapt_point *pe5 = adapt_point::New ( (p2->x+p3->x)*0.5,(p2->y+p3->y)*0.5,(p2->z+p3->z)*0.5,coeffs, eexps);
adapt_point *pe0 =
adapt_point::New((p0->x + p1->x) * 0.5, (p0->y + p1->y) * 0.5,
(p0->z + p1->z) * 0.5, coeffs, eexps);
adapt_point *pe1 =
adapt_point::New((p0->x + p2->x) * 0.5, (p0->y + p2->y) * 0.5,
(p0->z + p2->z) * 0.5, coeffs, eexps);
adapt_point *pe2 =
adapt_point::New((p0->x + p3->x) * 0.5, (p0->y + p3->y) * 0.5,
(p0->z + p3->z) * 0.5, coeffs, eexps);
adapt_point *pe3 =
adapt_point::New((p1->x + p2->x) * 0.5, (p1->y + p2->y) * 0.5,
(p1->z + p2->z) * 0.5, coeffs, eexps);
adapt_point *pe4 =
adapt_point::New((p1->x + p3->x) * 0.5, (p1->y + p3->y) * 0.5,
(p1->z + p3->z) * 0.5, coeffs, eexps);
adapt_point *pe5 =
adapt_point::New((p2->x + p3->x) * 0.5, (p2->y + p3->y) * 0.5,
(p2->z + p3->z) * 0.5, coeffs, eexps);
adapt_tet *t1 = new adapt_tet(p0, pe0, pe2, pe1);
Recur_Create(t1, maxlevel, level, coeffs, eexps);
......@@ -258,12 +291,18 @@ void adapt_tet::Recur_Create (adapt_tet *t, int maxlevel, int level , Double_Mat
adapt_tet *t8 = new adapt_tet(pe0, pe2, pe3, pe1);
Recur_Create(t8, maxlevel, level, coeffs, eexps);
t->e[0]=t1;t->e[1]=t2;t->e[2]=t3;t->e[3]=t4;
t->e[4]=t5;t->e[5]=t6;t->e[6]=t7;t->e[7]=t8;
t->e[0] = t1;
t->e[1] = t2;
t->e[2] = t3;
t->e[3] = t4;
t->e[4] = t5;
t->e[5] = t6;
t->e[6] = t7;
t->e[7] = t8;
}
void adapt_hex::Recur_Create (adapt_hex *h, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps)
void adapt_hex::Recur_Create(adapt_hex * h, int maxlevel, int level,
Double_Matrix * coeffs, Double_Matrix * eexps)
{
all_elems.push_back(h);
if(level++ >= maxlevel)
......@@ -277,29 +316,70 @@ void adapt_hex::Recur_Create (adapt_hex *h, int maxlevel, int level , Double_Mat
adapt_point *p5 = h->p[5];
adapt_point *p6 = h->p[6];
adapt_point *p7 = h->p[7];
adapt_point *p01 = adapt_point::New ( (p0->x+p1->x)*0.5,(p0->y+p1->y)*0.5,(p0->z+p1->z)*0.5,coeffs, eexps);
adapt_point *p12 = adapt_point::New ( (p1->x+p2->x)*0.5,(p1->y+p2->y)*0.5,(p1->z+p2->z)*0.5,coeffs, eexps);
adapt_point *p23 = adapt_point::New ( (p2->x+p3->x)*0.5,(p2->y+p3->y)*0.5,(p2->z+p3->z)*0.5,coeffs, eexps);
adapt_point *p03 = adapt_point::New ( (p3->x+p0->x)*0.5,(p3->y+p0->y)*0.5,(p3->z+p0->z)*0.5,coeffs, eexps);
adapt_point *p45 = adapt_point::New ( (p4->x+p5->x)*0.5,(p4->y+p5->y)*0.5,(p4->z+p5->z)*0.5,coeffs, eexps);
adapt_point *p56 = adapt_point::New ( (p5->x+p6->x)*0.5,(p5->y+p6->y)*0.5,(p5->z+p6->z)*0.5,coeffs, eexps);
adapt_point *p67 = adapt_point::New ( (p6->x+p7->x)*0.5,(p6->y+p7->y)*0.5,(p6->z+p7->z)*0.5,coeffs, eexps);
adapt_point *p47 = adapt_point::New ( (p7->x+p4->x)*0.5,(p7->y+p4->y)*0.5,(p7->z+p4->z)*0.5,coeffs, eexps);
adapt_point *p04 = adapt_point::New ( (p4->x+p0->x)*0.5,(p4->y+p0->y)*0.5,(p4->z+p0->z)*0.5,coeffs, eexps);
adapt_point *p15 = adapt_point::New ( (p5->x+p1->x)*0.5,(p5->y+p1->y)*0.5,(p5->z+p1->z)*0.5,coeffs, eexps);
adapt_point *p26 = adapt_point::New ( (p6->x+p2->x)*0.5,(p6->y+p2->y)*0.5,(p6->z+p2->z)*0.5,coeffs, eexps);
adapt_point *p37 = adapt_point::New ( (p7->x+p3->x)*0.5,(p7->y+p3->y)*0.5,(p7->z+p3->z)*0.5,coeffs, eexps);
adapt_point *p0145 = adapt_point::New ( (p45->x+p01->x)*0.5,(p45->y+p01->y)*0.5,(p45->z+p01->z)*0.5,coeffs, eexps);
adapt_point *p1256 = adapt_point::New ( (p12->x+p56->x)*0.5,(p12->y+p56->y)*0.5,(p12->z+p56->z)*0.5,coeffs, eexps);
adapt_point *p2367 = adapt_point::New ( (p23->x+p67->x)*0.5,(p23->y+p67->y)*0.5,(p23->z+p67->z)*0.5,coeffs, eexps);
adapt_point *p0347 = adapt_point::New ( (p03->x+p47->x)*0.5,(p03->y+p47->y)*0.5,(p03->z+p47->z)*0.5,coeffs, eexps);
adapt_point *p4756 = adapt_point::New ( (p47->x+p56->x)*0.5,(p47->y+p56->y)*0.5,(p47->z+p56->z)*0.5,coeffs, eexps);
adapt_point *p0312 = adapt_point::New ( (p03->x+p12->x)*0.5,(p03->y+p12->y)*0.5,(p03->z+p12->z)*0.5,coeffs, eexps);
adapt_point *pc = adapt_point::New ( (p0->x+p1->x+p2->x+p3->x+p4->x+p5->x+p6->x+p7->x)*0.125,
(p0->y+p1->y+p2->y+p3->y+p4->y+p5->y+p6->y+p7->y)*0.125,
(p0->z+p1->z+p2->z+p3->z+p4->z+p5->z+p6->z+p7->z)*0.125,
adapt_point *p01 =
adapt_point::New((p0->x + p1->x) * 0.5, (p0->y + p1->y) * 0.5,
(p0->z + p1->z) * 0.5, coeffs, eexps);
adapt_point *p12 =
adapt_point::New((p1->x + p2->x) * 0.5, (p1->y + p2->y) * 0.5,
(p1->z + p2->z) * 0.5, coeffs, eexps);
adapt_point *p23 =
adapt_point::New((p2->x + p3->x) * 0.5, (p2->y + p3->y) * 0.5,
(p2->z + p3->z) * 0.5, coeffs, eexps);
adapt_point *p03 =
adapt_point::New((p3->x + p0->x) * 0.5, (p3->y + p0->y) * 0.5,
(p3->z + p0->z) * 0.5, coeffs, eexps);
adapt_point *p45 =
adapt_point::New((p4->x + p5->x) * 0.5, (p4->y + p5->y) * 0.5,
(p4->z + p5->z) * 0.5, coeffs, eexps);
adapt_point *p56 =
adapt_point::New((p5->x + p6->x) * 0.5, (p5->y + p6->y) * 0.5,
(p5->z + p6->z) * 0.5, coeffs, eexps);
adapt_point *p67 =
adapt_point::New((p6->x + p7->x) * 0.5, (p6->y + p7->y) * 0.5,
(p6->z + p7->z) * 0.5, coeffs, eexps);
adapt_point *p47 =
adapt_point::New((p7->x + p4->x) * 0.5, (p7->y + p4->y) * 0.5,
(p7->z + p4->z) * 0.5, coeffs, eexps);
adapt_point *p04 =
adapt_point::New((p4->x + p0->x) * 0.5, (p4->y + p0->y) * 0.5,
(p4->z + p0->z) * 0.5, coeffs, eexps);
adapt_point *p15 =
adapt_point::New((p5->x + p1->x) * 0.5, (p5->y + p1->y) * 0.5,
(p5->z + p1->z) * 0.5, coeffs, eexps);
adapt_point *p26 =
adapt_point::New((p6->x + p2->x) * 0.5, (p6->y + p2->y) * 0.5,
(p6->z + p2->z) * 0.5, coeffs, eexps);
adapt_point *p37 =
adapt_point::New((p7->x + p3->x) * 0.5, (p7->y + p3->y) * 0.5,
(p7->z + p3->z) * 0.5, coeffs, eexps);
adapt_point *p0145 =
adapt_point::New((p45->x + p01->x) * 0.5, (p45->y + p01->y) * 0.5,
(p45->z + p01->z) * 0.5, coeffs, eexps);
adapt_point *p1256 =
adapt_point::New((p12->x + p56->x) * 0.5, (p12->y + p56->y) * 0.5,
(p12->z + p56->z) * 0.5, coeffs, eexps);
adapt_point *p2367 =
adapt_point::New((p23->x + p67->x) * 0.5, (p23->y + p67->y) * 0.5,
(p23->z + p67->z) * 0.5, coeffs, eexps);
adapt_point *p0347 =
adapt_point::New((p03->x + p47->x) * 0.5, (p03->y + p47->y) * 0.5,
(p03->z + p47->z) * 0.5, coeffs, eexps);
adapt_point *p4756 =
adapt_point::New((p47->x + p56->x) * 0.5, (p47->y + p56->y) * 0.5,
(p47->z + p56->z) * 0.5, coeffs, eexps);
adapt_point *p0312 =
adapt_point::New((p03->x + p12->x) * 0.5, (p03->y + p12->y) * 0.5,
(p03->z + p12->z) * 0.5, coeffs, eexps);
adapt_point *pc =
adapt_point::
New((p0->x + p1->x + p2->x + p3->x + p4->x + p5->x + p6->x +
p7->x) * 0.125,
(p0->y + p1->y + p2->y + p3->y + p4->y + p5->y + p6->y +
p7->y) * 0.125,
(p0->z + p1->z + p2->z + p3->z + p4->z + p5->z + p6->z +
p7->z) * 0.125,
coeffs, eexps);
adapt_hex *h1 = new adapt_hex(p0, p01, p0312, p03, p04, p0145, pc, p0347); //p0
......@@ -319,8 +399,14 @@ void adapt_hex::Recur_Create (adapt_hex *h, int maxlevel, int level , Double_Mat
adapt_hex *h8 = new adapt_hex(p0312, pc, p1256, p12, p23, p2367, p26, p2); //p2
Recur_Create(h8, maxlevel, level, coeffs, eexps);
h->e[0]=h1;h->e[1]=h2;h->e[2]=h3;h->e[3]=h4;
h->e[4]=h5;h->e[5]=h6;h->e[6]=h7;h->e[7]=h8;
h->e[0] = h1;
h->e[1] = h2;
h->e[2] = h3;
h->e[3] = h4;
h->e[4] = h5;
h->e[5] = h6;
h->e[6] = h7;
h->e[7] = h8;
}
void adapt_edge::Error(double AVG, double tol)
......@@ -354,12 +440,11 @@ void adapt_hex::Error ( double AVG , double tol )
void adapt_edge::Recur_Error(adapt_edge * e, double AVG, double tol)
{
if(!e->e[0])e->visible = true;
else
{
if(!e->e[0])
e->visible = true;
else {
double vr;
if (!e->e[0]->e[0])
{
if(!e->e[0]->e[0]) {
double v1 = e->e[0]->V();
double v2 = e->e[1]->V();
vr = (v1 + v2) / 2.;
......@@ -374,8 +459,7 @@ void adapt_edge::Recur_Error ( adapt_edge *e, double AVG, double tol )
else
e->visible = true;
}
else
{
else {
double v11 = e->e[0]->e[0]->V();
double v12 = e->e[0]->e[1]->V();
......@@ -387,8 +471,7 @@ void adapt_edge::Recur_Error ( adapt_edge *e, double AVG, double tol )
vr = (vr1 + vr2) / 2.;
if(fabs(e->e[0]->V() - vr1) > AVG * tol ||
fabs(e->e[1]->V() - vr2) > AVG * tol ||
fabs(e->V() - vr) > AVG * tol )
{
fabs(e->V() - vr) > AVG * tol) {
e->visible = false;
Recur_Error(e->e[0], AVG, tol);
Recur_Error(e->e[1], AVG, tol);
......@@ -401,12 +484,11 @@ void adapt_edge::Recur_Error ( adapt_edge *e, double AVG, double tol )
void adapt_triangle::Recur_Error(adapt_triangle * t, double AVG, double tol)
{
if(!t->e[0])t->visible = true;
else
{
if(!t->e[0])
t->visible = true;
else {
double vr;
if (!t->e[0]->e[0])
{
if(!t->e[0]->e[0]) {
double v1 = t->e[0]->V();
double v2 = t->e[1]->V();
double v3 = t->e[2]->V();
......@@ -425,8 +507,7 @@ void adapt_triangle::Recur_Error ( adapt_triangle *t, double AVG, double tol )
else
t->visible = true;
}
else
{
else {
double v11 = t->e[0]->e[0]->V();
double v12 = t->e[0]->e[1]->V();
double v13 = t->e[0]->e[2]->V();
......@@ -473,12 +554,11 @@ void adapt_triangle::Recur_Error ( adapt_triangle *t, double AVG, double tol )
void adapt_quad::Recur_Error(adapt_quad * q, double AVG, double tol)
{
if(!q->e[0])q->visible = true;
else
{
if(!q->e[0])
q->visible = true;
else {
double vr;
if (!q->e[0]->e[0])
{
if(!q->e[0]->e[0]) {
double v1 = q->e[0]->V();
double v2 = q->e[1]->V();
double v3 = q->e[2]->V();
......@@ -497,8 +577,7 @@ void adapt_quad::Recur_Error ( adapt_quad *q, double AVG, double tol )
else
q->visible = true;
}
else
{
else {
double v11 = q->e[0]->e[0]->V();
double v12 = q->e[0]->e[1]->V();
double v13 = q->e[0]->e[2]->V();
......@@ -543,12 +622,11 @@ void adapt_quad::Recur_Error ( adapt_quad *q, double AVG, double tol )
}
}
void adapt_tet::Recur_Error(adapt_tet * t, double AVG, double tol)
{
if(!t->e[0])t->visible = true;
else
{
if(!t->e[0])
t->visible = true;
else {
const double v1 = t->e[0]->V();
const double v2 = t->e[1]->V();
const double v3 = t->e[2]->V();
......@@ -559,10 +637,8 @@ void adapt_tet::Recur_Error ( adapt_tet *t, double AVG, double tol )
const double v8 = t->e[7]->V();
const double vr = (v1 + v2 + v3 + v4 + v5 + v6 + v7 + v8) * .125;
const double v = t->V();
if (!t->e[0]->e[0])
{
if ( fabs(v - vr) > AVG * tol )
{
if(!t->e[0]->e[0]) {
if(fabs(v - vr) > AVG * tol) {
t->visible = false;
Recur_Error(t->e[0], AVG, tol);
Recur_Error(t->e[1], AVG, tol);
......@@ -576,18 +652,15 @@ void adapt_tet::Recur_Error ( adapt_tet *t, double AVG, double tol )
else
t->visible = true;
}
else
{
else {
double vi[8][8];
for(int k = 0; k < 8; k++)
for(int l = 0; l < 8; l++)
vi[k][l] = t->e[k]->e[l]->V();
double vri[8];
for (int k=0;k<8;k++)
{
for(int k = 0; k < 8; k++) {
vri[k] = 0.0;
for (int l=0;l<8;l++)
{
for(int l = 0; l < 8; l++) {
vri[k] += vi[k][l];
}
vri[k] /= 8.0;
......@@ -600,8 +673,7 @@ void adapt_tet::Recur_Error ( adapt_tet *t, double AVG, double tol )
fabs(t->e[5]->V() - vri[5]) > AVG * tol ||
fabs(t->e[6]->V() - vri[6]) > AVG * tol ||
fabs(t->e[7]->V() - vri[7]) > AVG * tol ||
fabs(v - vr) > AVG * tol )
{
fabs(v - vr) > AVG * tol) {
t->visible = false;
Recur_Error(t->e[0], AVG, tol);
Recur_Error(t->e[1], AVG, tol);
......@@ -620,9 +692,9 @@ void adapt_tet::Recur_Error ( adapt_tet *t, double AVG, double tol )
void adapt_hex::Recur_Error(adapt_hex * h, double AVG, double tol)
{
if(!h->e[0])h->visible = true;
else
{
if(!h->e[0])
h->visible = true;
else {
double vr;
double v1 = h->e[0]->V();
double v2 = h->e[1]->V();
......@@ -634,8 +706,7 @@ void adapt_hex::Recur_Error ( adapt_hex *h, double AVG, double tol )
double v8 = h->e[7]->V();
vr = (v1 + v2 + v3 + v4 + v5 + v6 + v7 + v8) * .125;
double v = h->V();
if ( fabs(v - vr) > AVG * tol )
{
if(fabs(v - vr) > AVG * tol) {
h->visible = false;
Recur_Error(h->e[0], AVG, tol);
Recur_Error(h->e[1], AVG, tol);
......@@ -659,14 +730,14 @@ int Adaptive_Post_View:: zoomElement (Post_View * view ,
int level,
int levelmax,
GMSH_Post_Plugin * plug,
List_T *theList,
int *counter)
List_T * theList, int *counter)
{
const int nbNod = ELEM::nbNod;
typename std::set<adapt_point>::iterator it = adapt_point::all_points.begin();
typename std::set<adapt_point>::iterator ite = adapt_point::all_points.end();
typename std::set < adapt_point >::iterator it =
adapt_point::all_points.begin();
typename std::set < adapt_point >::iterator ite =
adapt_point::all_points.end();
double c0 = Cpu();
......@@ -676,15 +747,13 @@ int Adaptive_Post_View:: zoomElement (Post_View * view ,
Double_Matrix xyz(nbNod, 3);
Double_Matrix XYZ(adapt_point::all_points.size(), 3);
for ( int k=0;k<nbNod;++k)
{
for(int k = 0; k < nbNod; ++k) {
xyz(k, 0) = (*_STposX) (ielem, k);
xyz(k, 1) = (*_STposY) (ielem, k);
xyz(k, 2) = (*_STposZ) (ielem, k);
}
for ( int k=0;k<N;++k)
{
for(int k = 0; k < N; ++k) {
val(k) = (*_STval) (ielem, k);
}
_Interpolate->mult(val, res);
......@@ -692,15 +761,16 @@ int Adaptive_Post_View:: zoomElement (Post_View * view ,
double c1 = Cpu();
int kk = 0;
for ( ; it !=ite ; ++it)
{
for(; it != ite; ++it) {
adapt_point *p = (adapt_point *) & (*it);
p->val = res(kk);
p->X = XYZ(kk, 0);
p->Y = XYZ(kk, 1);
p->Z = XYZ(kk, 2);
if (min > p->val) min = p->val;
if (max < p->val) max = p->val;
if(min > p->val)
min = p->val;
if(max < p->val)
max = p->val;
kk++;
}
double c2 = Cpu();
......@@ -708,14 +778,12 @@ int Adaptive_Post_View:: zoomElement (Post_View * view ,
typename std::list < ELEM * >::iterator itt = ELEM::all_elems.begin();
typename std::list < ELEM * >::iterator itte = ELEM::all_elems.end();
for ( ;itt != itte ; itt++)
{
for(; itt != itte; itt++) {
(*itt)->visible = false;
}
if (!plug || tol != 0.0)
{
if(!plug || tol != 0.0) {
ELEM::Error(max - min, tol);
}
if(plug)
......@@ -726,24 +794,26 @@ int Adaptive_Post_View:: zoomElement (Post_View * view ,
itt = ELEM::all_elems.begin();
for ( ;itt != itte ; itt++)
{
if ((*itt)->visible && !(*itt)->e[0] && level != levelmax )return 0;
for(; itt != itte; itt++) {
if((*itt)->visible && !(*itt)->e[0] && level != levelmax)
return 0;
}
itt = ELEM::all_elems.begin();
adapt_point **p;
for ( ;itt != itte ; itt++)
{
if ((*itt)->visible)
{
for(; itt != itte; itt++) {
if((*itt)->visible) {
p = (*itt)->p;
for (int k=0;k<nbNod;++k)List_Add ( theList , &p[k]->X );
for (int k=0;k<nbNod;++k)List_Add ( theList , &p[k]->Y );
for (int k=0;k<nbNod;++k)List_Add ( theList , &p[k]->Z );
for (int k=0;k<nbNod;++k)List_Add ( theList , &p[k]->val );
for(int k = 0; k < nbNod; ++k)
List_Add(theList, &p[k]->X);
for(int k = 0; k < nbNod; ++k)
List_Add(theList, &p[k]->Y);
for(int k = 0; k < nbNod; ++k)
List_Add(theList, &p[k]->Z);
for(int k = 0; k < nbNod; ++k)
List_Add(theList, &p[k]->val);
(*counter)++;
}
}
......@@ -765,56 +835,52 @@ int Adaptive_Post_View:: zoomElement (Post_View * view ,
We compute first the matrix at max recursion level (those should be stored
on disk once in the GMSHPLUGINSHOME directory, i'll do that at some point).
*/
void Adaptive_Post_View:: setAdaptiveResolutionLevel (Post_View * view , int level, GMSH_Post_Plugin *plug)
void Adaptive_Post_View::setAdaptiveResolutionLevel(Post_View * view,
int level,
GMSH_Post_Plugin * plug)
{
if (presentTol==tol && presentZoomLevel == level && !plug)return;
if(presentTol == tol && presentZoomLevel == level && !plug)
return;
printf("calling setAdaptive with level = %d and plug = %p tol %g presentTol %g\n",
printf
("calling setAdaptive with level = %d and plug = %p tol %g presentTol %g\n",
level, (void *)plug, tol, presentTol);
int *done = new int[_STposX->size1()];
for (int i=0;i<_STposX->size1();++i)done[i] = 0;
for(int i = 0; i < _STposX->size1(); ++i)
done[i] = 0;
int level_act = (level > 2) ? 2 : level;
int nbelm = _STposX->size1();
int TYP = 0;
if (view->NbSL)
{
if(view->NbSL) {
TYP = 5;
List_Delete(view->SL);
view->NbSL = 0;
view->SL = List_Create(nbelm * 8, nbelm, sizeof(double));
}
if (view->NbST)
{
if(view->NbST) {
TYP = 1;
List_Delete(view->ST);
view->NbST = 0;
view->ST = List_Create(nbelm * 4, nbelm, sizeof(double));
}
if (view->NbSS)
{
if(view->NbSS) {
TYP = 3;
List_Delete(view->SS);
view->NbSS = 0;
view->SS = List_Create(nbelm * 4, nbelm, sizeof(double));
}
if (view->NbSQ)
{
if(view->NbSQ) {
TYP = 2;
List_Delete(view->SQ);
view->NbSQ = 0;
view->SQ = List_Create(nbelm * 20, nbelm * 20, sizeof(double));
}
if (view->NbSH)
{
if(view->NbSH) {
TYP = 4;
List_Delete(view->SH);
view->NbSH = 0;
......@@ -823,21 +889,40 @@ void Adaptive_Post_View:: setAdaptiveResolutionLevel (Post_View * view , int lev
view->NbTimeStep = 1;
t0 = t1 = t2 = t3 = 0;
while (1)
{
if (TYP == 5)setAdaptiveResolutionLevel_TEMPL <adapt_edge> ( view,level_act,level, plug,&(view->SL),&(view->NbSL),done) ;
if (TYP == 1)setAdaptiveResolutionLevel_TEMPL <adapt_triangle> ( view,level_act,level, plug,&(view->ST),&(view->NbST),done) ;
if (TYP == 2)setAdaptiveResolutionLevel_TEMPL <adapt_quad> ( view,level_act,level, plug,&(view->SQ),&(view->NbSQ),done) ;
if (TYP == 4)setAdaptiveResolutionLevel_TEMPL <adapt_hex> ( view,level_act,level, plug,&(view->SH),&(view->NbSH),done) ;
if (TYP == 3)setAdaptiveResolutionLevel_TEMPL <adapt_tet> ( view,level_act,level, plug,&(view->SS),&(view->NbSS),done) ;
while(1) {
if(TYP == 5)
setAdaptiveResolutionLevel_TEMPL < adapt_edge > (view, level_act, level,
plug, &(view->SL),
&(view->NbSL), done);
if(TYP == 1)
setAdaptiveResolutionLevel_TEMPL < adapt_triangle > (view, level_act,
level, plug,
&(view->ST),
&(view->NbST),
done);
if(TYP == 2)
setAdaptiveResolutionLevel_TEMPL < adapt_quad > (view, level_act, level,
plug, &(view->SQ),
&(view->NbSQ), done);
if(TYP == 4)
setAdaptiveResolutionLevel_TEMPL < adapt_hex > (view, level_act, level,
plug, &(view->SH),
&(view->NbSH), done);
if(TYP == 3)
setAdaptiveResolutionLevel_TEMPL < adapt_tet > (view, level_act, level,
plug, &(view->SS),
&(view->NbSS), done);
int nbDone = 0;
for (int i=0;i<_STposX->size1();++i)nbDone += done[i];
printf("adaptive %d %d %d %d\n",level, level_act, nbDone, _STposX->size1());
if (nbDone ==_STposX->size1()) break;
if (level_act >= level) break;
for(int i = 0; i < _STposX->size1(); ++i)
nbDone += done[i];
printf("adaptive %d %d %d %d\n", level, level_act, nbDone,
_STposX->size1());
if(nbDone == _STposX->size1())
break;
if(level_act >= level)
break;
level_act++;
}
......@@ -853,7 +938,8 @@ template <class ELEM>
void Adaptive_Post_View::setAdaptiveResolutionLevel_TEMPL(Post_View * view,
int level,
int levelmax,
GMSH_Post_Plugin *plug,
GMSH_Post_Plugin *
plug,
List_T ** myList,
int *counter,
int *done)
......@@ -875,42 +961,39 @@ void Adaptive_Post_View:: setAdaptiveResolutionLevel_TEMPL (Post_View * view ,
_Geometry = new Double_Matrix(adapt_point::all_points.size(), ELEM::nbNod);
int kk = 0;
for ( ; it !=ite ; ++it)
{
for(; it != ite; ++it) {
adapt_point *p = (adapt_point *) & (*it);
for ( int k=0;k<N;++k)
{
for(int k = 0; k < N; ++k) {
(*_Interpolate) (kk, k) = p->shape_functions[k];
}
ELEM::GSF(p->x, p->y, p->z, sf);
for (int k=0;k<ELEM::nbNod;k++)(*_Geometry)(kk,k) = sf[k];
for(int k = 0; k < ELEM::nbNod; k++)
(*_Geometry) (kk, k) = sf[k];
kk++;
}
for ( int i=0;i<nbelm;++i)
{
for(int i = 0; i < nbelm; ++i) {
if(!done[i])
done[i] = zoomElement<ELEM> ( view , i , level, levelmax, plug,*myList,counter);
done[i] =
zoomElement < ELEM > (view, i, level, levelmax, plug, *myList,
counter);
}
}
void computeShapeFunctions ( Double_Matrix *coeffs, Double_Matrix *eexps , double u, double v, double w,double *sf)
void computeShapeFunctions(Double_Matrix * coeffs, Double_Matrix * eexps,
double u, double v, double w, double *sf)
{
static double powsuvw[256];
for (int j=0;j<coeffs->size2();++j)
{
for(int j = 0; j < coeffs->size2(); ++j) {
double powu = (*eexps) (j, 0);
double powv = (*eexps) (j, 1);
double poww = (*eexps) (j, 2);
powsuvw[j] = pow(u, powu) * pow(v, powv) * pow(w, poww);
}
for (int i=0;i<coeffs->size1();++i)
{
for(int i = 0; i < coeffs->size1(); ++i) {
sf[i] = 0.0;
for (int j=0;j<coeffs->size2();++j)
{
for(int j = 0; j < coeffs->size2(); ++j) {
sf[i] += (*coeffs) (i, j) * powsuvw[j];
}
}
......@@ -918,41 +1001,36 @@ void computeShapeFunctions ( Double_Matrix *coeffs, Double_Matrix *eexps , doubl
void Adaptive_Post_View::initWithLowResolution(Post_View * view)
{
List_T *myList;
int nbelm;
int nbnod;
if (view->NbSL)
{
if(view->NbSL) {
myList = view->SL;
nbelm = view->NbSL;
nbnod = 2;
}
else if (view->NbST)
{
else if(view->NbST) {
myList = view->ST;
nbelm = view->NbST;
nbnod = 3;
}
else if (view->NbSS)
{
else if(view->NbSS) {
myList = view->SS;
nbelm = view->NbSS;
nbnod = 4;
}
else if (view->NbSQ)
{
else if(view->NbSQ) {
myList = view->SQ;
nbelm = view->NbSQ;
nbnod = 4;
}
else if (view->NbSH)
{
else if(view->NbSH) {
myList = view->SH;
nbelm = view->NbSH;
nbnod = 8;
}
else return;
else
return;
min = VAL_INF;
max = -VAL_INF;
......@@ -966,13 +1044,11 @@ void Adaptive_Post_View:: initWithLowResolution (Post_View *view)
/// Store non interpolated data
int k = 0;
for (int i=0;i<List_Nbr(myList);i+=nb)
{
for(int i = 0; i < List_Nbr(myList); i += nb) {
double *x = (double *)List_Pointer_Fast(myList, i);
double *y = (double *)List_Pointer_Fast(myList, i + nbnod);
double *z = (double *)List_Pointer_Fast(myList, i + 2 * nbnod);
for (int NN=0;NN<nbnod;NN++)
{
for(int NN = 0; NN < nbnod; NN++) {
(*_STposX) (k, NN) = x[NN];
(*_STposY) (k, NN) = y[NN];
(*_STposZ) (k, NN) = z[NN];
......@@ -986,7 +1062,8 @@ void Adaptive_Post_View:: initWithLowResolution (Post_View *view)
setAdaptiveResolutionLevel(view, 0);
}
Adaptive_Post_View:: Adaptive_Post_View (Post_View *view, List_T *_c , List_T *_pol)
Adaptive_Post_View::Adaptive_Post_View(Post_View * view, List_T * _c,
List_T * _pol)
:tol(1.e-3)
{
......@@ -994,8 +1071,7 @@ Adaptive_Post_View:: Adaptive_Post_View (Post_View *view, List_T *_c , List_T *_
_coefs = new Double_Matrix(List_Nbr(_c), List_Nbr(_c));
_eexps = new Double_Matrix(List_Nbr(_c), 3);
for (int i=0; i< List_Nbr ( _c ); ++i)
{
for(int i = 0; i < List_Nbr(_c); ++i) {
List_T **line = (List_T **) List_Pointer_Fast(_c, i);
List_T **eexp = (List_T **) List_Pointer_Fast(_pol, i);
......@@ -1009,8 +1085,7 @@ Adaptive_Post_View:: Adaptive_Post_View (Post_View *view, List_T *_c , List_T *_
(*_eexps) (i, 1) = dpowv;
(*_eexps) (i, 2) = dpoww;
for (int j=0;j < List_Nbr ( *line ); ++j)
{
for(int j = 0; j < List_Nbr(*line); ++j) {
double val;
List_Read(*line, j, &val);
(*_coefs) (i, j) = val;
......@@ -1027,8 +1102,10 @@ Adaptive_Post_View::~Adaptive_Post_View()
delete _STposY;
delete _STposZ;
delete _STval;
if(_Interpolate)delete _Interpolate;
if(_Geometry)delete _Geometry;
if(_Interpolate)
delete _Interpolate;
if(_Geometry)
delete _Geometry;
cleanElement < adapt_edge > ();
cleanElement < adapt_triangle > ();
cleanElement < adapt_tet > ();
......
......@@ -29,18 +29,19 @@
class Post_View;
class GMSH_Post_Plugin;
// On a triangle, we suppose that there exists an
// interpolation scheme such that u = \sum_i u_i \phi_i
// phi_i being polynomials of order p, i goes from 1...(p+1)(p+2)/2
// and phi_i = \sum_j coeffs_{ij} monomials_j and
// monomials are 1,x,y,x^2,xy,y^2,x^3,x^2y,xy^2,y^3...
// On a triangle, we suppose that there exists an interpolation scheme
// such that u = \sum_i u_i \phi_i, phi_i being polynomials of order
// p, i goes from 1...(p+1)(p+2)/2 and phi_i = \sum_j coeffs_{ij}
// monomials_j and monomials are 1,x,y,x^2,xy,y^2,x^3,x^2y,xy^2,y^3...
class adapt_point
{
public :
double x, y, z;
double X, Y, Z, val;
double shape_functions[128];
static adapt_point * New ( double x, double y, double z, Double_Matrix *coeffs, Double_Matrix *eexps);
static adapt_point * New (double x, double y, double z,
Double_Matrix *coeffs, Double_Matrix *eexps);
void print () const
{
printf("p %g %g\n" ,x,y);
......@@ -81,7 +82,8 @@ public:
printf("p1 %g %g p2 %g %g\n", p[0]->x, p[0]->y, p[1]->x, p[1]->y);
}
static void Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps) ;
static void Recur_Create (adapt_edge *e, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_edge *e, int maxlevel, int level,
Double_Matrix *coeffs, Double_Matrix *eexps);
static void Error (double AVG , double tol );
static void Recur_Error (adapt_edge *e, double AVG, double tol);
bool visible;
......@@ -91,8 +93,6 @@ public:
static int nbNod;
};
class adapt_triangle
{
public:
......@@ -104,7 +104,6 @@ public:
p[2] = p3;
e[0] = e[1] = e[2] = e[3] = 0;
}
inline double V () const
{
return (p[0]->val + p[1]->val + p[2]->val)/3.;
......@@ -117,10 +116,12 @@ public:
}
void print ()
{
printf ("p1 %g %g p2 %g %g p3 %g %g \n",p[0]->x,p[0]->y,p[1]->x,p[1]->y,p[2]->x,p[2]->y);
printf("p1 %g %g p2 %g %g p3 %g %g\n",
p[0]->x, p[0]->y, p[1]->x, p[1]->y, p[2]->x, p[2]->y);
}
static void Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_triangle *t, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_triangle *t, int maxlevel, int level,
Double_Matrix *coeffs, Double_Matrix *eexps);
static void Error (double AVG, double tol);
static void Recur_Error (adapt_triangle *t, double AVG, double tol);
bool visible;
......@@ -142,7 +143,6 @@ public:
p[3] = p4;
e[0] = e[1] = e[2] = e[3] = 0;
}
inline double V () const
{
return (p[0]->val + p[1]->val + p[2]->val + p[3]->val)/4.;
......@@ -156,10 +156,12 @@ public:
}
void print ()
{
printf ("p1 %g %g p2 %g %g p3 %g %g \n",p[0]->x,p[0]->y,p[1]->x,p[1]->y,p[2]->x,p[2]->y);
printf("p1 %g %g p2 %g %g p3 %g %g\n",
p[0]->x, p[0]->y, p[1]->x, p[1]->y, p[2]->x, p[2]->y);
}
static void Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_quad *q, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_quad *q, int maxlevel, int level,
Double_Matrix *coeffs, Double_Matrix *eexps);
static void Error (double AVG, double tol);
static void Recur_Error (adapt_quad *q, double AVG, double tol);
bool visible;
......@@ -182,7 +184,6 @@ public:
e[0] = e[1] = e[2] = e[3] = 0;
e[4] = e[5] = e[6] = e[7] = 0;
}
inline static void GSF (const double u, const double v, double w, double sf[])
{
sf[0] = 1-u-v-w;
......@@ -196,10 +197,12 @@ public:
}
void print ()
{
printf ("p1 %g %g p2 %g %g p3 %g %g \n",p[0]->x,p[0]->y,p[1]->x,p[1]->y,p[2]->x,p[2]->y);
printf("p1 %g %g p2 %g %g p3 %g %g\n",
p[0]->x, p[0]->y, p[1]->x, p[1]->y, p[2]->x, p[2]->y);
}
static void Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_tet *t, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_tet *t, int maxlevel, int level,
Double_Matrix *coeffs, Double_Matrix *eexps);
static void Error (double AVG, double tol);
static void Recur_Error (adapt_tet *t, double AVG, double tol);
bool visible;
......@@ -212,7 +215,8 @@ public:
class adapt_hex
{
public:
adapt_hex (adapt_point *p1,adapt_point *p2,adapt_point *p3,adapt_point *p4,adapt_point *p5,adapt_point *p6,adapt_point *p7,adapt_point *p8)
adapt_hex (adapt_point *p1, adapt_point *p2, adapt_point *p3, adapt_point *p4,
adapt_point *p5, adapt_point *p6, adapt_point *p7, adapt_point *p8)
: visible (false)
{
p[0] = p1;
......@@ -226,7 +230,6 @@ public:
e[0] = e[1] = e[2] = e[3] = 0;
e[4] = e[5] = e[6] = e[7] = 0;
}
inline static void GSF (const double u, const double v, double w, double sf[])
{
sf[0] = 0.125 * (1 - u) * (1 - v) * (1 - w);
......@@ -240,14 +243,17 @@ public:
}
inline double V () const
{
return (p[0]->val + p[1]->val + p[2]->val+ p[3]->val+p[4]->val + p[5]->val + p[6]->val+ p[7]->val)/8.;
return (p[0]->val + p[1]->val + p[2]->val+ p[3]->val +
p[4]->val + p[5]->val + p[6]->val+ p[7]->val)/8.;
}
void print ()
{
printf ("p1 %g %g p2 %g %g p3 %g %g \n",p[0]->x,p[0]->y,p[1]->x,p[1]->y,p[2]->x,p[2]->y);
printf("p1 %g %g p2 %g %g p3 %g %g\n",
p[0]->x, p[0]->y, p[1]->x, p[1]->y, p[2]->x, p[2]->y);
}
static void Create (int maxlevel, Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_hex *h, int maxlevel, int level , Double_Matrix *coeffs, Double_Matrix *eexps);
static void Recur_Create (adapt_hex *h, int maxlevel, int level,
Double_Matrix *coeffs, Double_Matrix *eexps);
static void Error (double AVG, double tol);
static void Recur_Error (adapt_hex *h, double AVG, double tol);
bool visible;
......@@ -280,14 +286,17 @@ public:
setAdaptiveResolutionLevel(view, level);
}
template <class ELEM>
void setAdaptiveResolutionLevel_TEMPL (Post_View * view , int level, int lemvelmax, GMSH_Post_Plugin *plug, List_T **myList, int *counter, int *done);
void setAdaptiveResolutionLevel ( Post_View * view , int levelmax, GMSH_Post_Plugin *plug = 0);
void setAdaptiveResolutionLevel_TEMPL(Post_View * view, int level, int lemvelmax,
GMSH_Post_Plugin *plug, List_T **myList,
int *counter, int *done);
void setAdaptiveResolutionLevel (Post_View * view , int levelmax,
GMSH_Post_Plugin *plug = 0);
void initWithLowResolution (Post_View *view);
void setTolerance (const double eps) { tol=eps; }
double getTolerance () const { return tol; }
template <class ELEM>
int zoomElement (Post_View * view ,
int ielem, int level, int levelmax, GMSH_Post_Plugin *plug, List_T *theList, int *counter);
int zoomElement (Post_View * view, int ielem, int level, int levelmax,
GMSH_Post_Plugin *plug, List_T *theList, int *counter);
};
template <class ELEM>
......@@ -295,14 +304,11 @@ void cleanElement ()
{
typename std::list<ELEM*>::iterator it = ELEM::all_elems.begin();
typename std::list<ELEM*>::iterator ite = ELEM::all_elems.end();
for (;it!=ite;++it)
{
for (; it != ite; ++it){
delete *it;
}
ELEM::all_elems.clear();
adapt_point::all_points.clear();
}
#endif
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment