Index: /issm/trunk/src/c/Bamgx/Bamgx.cpp =================================================================== --- /issm/trunk/src/c/Bamgx/Bamgx.cpp (revision 2793) +++ /issm/trunk/src/c/Bamgx/Bamgx.cpp (revision 2794) @@ -109,27 +109,5 @@ //Build output - if (verbosity>1) printf(" Processing Output...\n"); - NumVerticesOut=Th.nbv; - VerticesOut=(double*)xmalloc(3*NumVerticesOut*sizeof(double)); - for (i=0;iNumTriangles=NumTrianglesOut; - xfree((void**)&bamgmesh->Triangles); - bamgmesh->Triangles=TrianglesOut; - - bamgmesh->NumVertices=NumVerticesOut; - xfree((void**)&bamgmesh->Vertices); - bamgmesh->Vertices=VerticesOut; - /*No error return*/ if (verbosity>1) printf(" Exiting Bamg.\n"); Index: /issm/trunk/src/c/Bamgx/Mesh2.h =================================================================== --- /issm/trunk/src/c/Bamgx/Mesh2.h (revision 2793) +++ /issm/trunk/src/c/Bamgx/Mesh2.h (revision 2794) @@ -899,4 +899,5 @@ void Read(MeshIstream &,int version,Real8 cutoffradian); void ReadFromMatlabMesh(BamgMesh* bamgmesh, BamgOpts* bamgopts); + void WriteBamgMesh(BamgMesh* bamgmesh,BamgOpts* bamgopts); void Read_am_fmt(MeshIstream &); void Read_amdba(MeshIstream &); Index: /issm/trunk/src/c/Bamgx/MeshRead.cpp =================================================================== --- /issm/trunk/src/c/Bamgx/MeshRead.cpp (revision 2793) +++ /issm/trunk/src/c/Bamgx/MeshRead.cpp (revision 2794) @@ -132,16 +132,15 @@ } - //VertexOnGeometricEdge - if(bamgmesh->VertexOnGeometricEdge){ - if(verbose>3) printf(" processing VertexOnGeometricEdge\n"); - NbVerticesOnGeomEdge=bamgmesh->NumVertexOnGeometricEdge; + //VerVerticesOnGeometricEdge + if(bamgmesh->VerticesOnGeometricEdge){ + if(verbose>3) printf(" processing VerticesOnGeometricEdge\n"); + NbVerticesOnGeomEdge=bamgmesh->NumVerticesOnGeometricEdge; VerticesOnGeomEdge= new VertexOnGeom[NbVerticesOnGeomEdge] ; for (i=0;iVertexOnGeometricEdge[i*3+0]-1; //for C indexing - i2=(Int4)bamgmesh->VertexOnGeometricEdge[i*3+1]-1; //for C indexing - s =(Int4)bamgmesh->VertexOnGeometricEdge[i*3+2]; + i1=(Int4)bamgmesh->VerticesOnGeometricEdge[i*3+0]-1; //for C indexing + i2=(Int4)bamgmesh->VerticesOnGeometricEdge[i*3+1]-1; //for C indexing + s =(Int4)bamgmesh->VerticesOnGeometricEdge[i*3+2]; VerticesOnGeomEdge[i]=VertexOnGeom(vertices[i1],Gh.edges[i2],s); } @@ -224,13 +223,13 @@ //EdgeOnGeometricEdge - if(bamgmesh->EdgeOnGeometricEdge){ - if(verbose>3) printf(" processing EdgeOnGeometricEdge\n"); + if(bamgmesh->EdgesOnGeometricEdge){ + if(verbose>3) printf(" processing EdgesOnGeometricEdge\n"); int i1,i2,i,j; - i2=bamgmesh->NumEdgeOnGeometricEdge; + i2=bamgmesh->NumEdgesOnGeometricEdge; for (i1=0;i1EdgeOnGeometricEdge[i*2+0]; - j=(int)bamgmesh->EdgeOnGeometricEdge[i*2+1]; + i=(int)bamgmesh->EdgesOnGeometricEdge[i*2+0]; + j=(int)bamgmesh->EdgesOnGeometricEdge[i*2+1]; if(!(i>0 && j >0 && i <= nbe && j <= Gh.nbe)) { - throw ErrorException(__FUNCT__,exprintf("EdgeOnGeometricEdge error: We must have : (i>0 && j >0 && i <= nbe && j <= Gh.nbe)")); + throw ErrorException(__FUNCT__,exprintf("EdgesOnGeometricEdge error: We must have : (i>0 && j >0 && i <= nbe && j <= Gh.nbe)")); } edges[i-1].on = Gh.edges + j-1; @@ -238,17 +237,17 @@ } else{ - if(verbose>3) printf(" no EdgeOnGeometricEdge found\n"); + if(verbose>3) printf(" no EdgesOnGeometricEdge found\n"); } //SubDomain - if(bamgmesh->SubDomain){ + if(bamgmesh->SubDomains){ Int4 i3,head,sens; - if(verbose>3) printf(" processing SubDomain\n"); - NbSubDomains=bamgmesh->NumSubDomain; + if(verbose>3) printf(" processing SubDomains\n"); + NbSubDomains=bamgmesh->NumSubDomains; subdomains = new SubDomain [ NbSubDomains ]; for (i=0;iSubDomain[i*3+0]; - head=(int)bamgmesh->SubDomain[i*3+1]-1;//C indexing - sens=(int)bamgmesh->SubDomain[i*3+2]; + i3 =(int)bamgmesh->SubDomains[i*3+0]; + head=(int)bamgmesh->SubDomains[i*3+1]-1;//C indexing + sens=(int)bamgmesh->SubDomains[i*3+2]; if (i3!=23) throw ErrorException(__FUNCT__,exprintf("Bad Subdomain definition: first number should be 3")); if (head<0 || head>=nbt) throw ErrorException(__FUNCT__,exprintf("Bad Subdomain definition: head should in [1 %i] (triangle number)",nbt)); @@ -257,5 +256,5 @@ } else{ - if(verbose>3) printf(" no SubDomain found\n"); + if(verbose>3) printf(" no SubDomains found\n"); } @@ -1391,14 +1390,14 @@ //SubDomain - if(bamggeom->SubDomain){ + if(bamggeom->SubDomains){ Int4 i0,i1,i2,i3; - if(verbose>3) printf(" processing SubDomain\n"); - NbSubDomains=bamggeom->NumSubDomain; + if(verbose>3) printf(" processing SubDomains\n"); + NbSubDomains=bamggeom->NumSubDomains; subdomains = new GeometricalSubDomain[NbSubDomains]; for (i=0;iSubDomain[i*4+0]; - i1=(int)bamggeom->SubDomain[i*4+1]; - i2=(int)bamggeom->SubDomain[i*4+2]; - i3=(int)bamggeom->SubDomain[i*4+3]; + i0=(int)bamggeom->SubDomains[i*4+0]; + i1=(int)bamggeom->SubDomains[i*4+1]; + i2=(int)bamggeom->SubDomains[i*4+2]; + i3=(int)bamggeom->SubDomains[i*4+3]; if (i0!=2) throw ErrorException(__FUNCT__,exprintf("Bad Subdomain definition: first number should be 2 (for Edges)")); if (i1>nbe || i1<=0) throw ErrorException(__FUNCT__,exprintf("Bad Subdomain definition: second number should in [1 %i] (edge number)",nbe)); @@ -1409,5 +1408,5 @@ } else{ - if(verbose>3) printf(" no SubDomain found\n"); + if(verbose>3) printf(" no SubDomains found\n"); } } Index: /issm/trunk/src/c/Bamgx/MeshWrite.cpp =================================================================== --- /issm/trunk/src/c/Bamgx/MeshWrite.cpp (revision 2793) +++ /issm/trunk/src/c/Bamgx/MeshWrite.cpp (revision 2794) @@ -1,29 +1,6 @@ -// -*- Mode : c++ -*- -// -// SUMMARY : -// USAGE : -// ORG : -// AUTHOR : Frederic Hecht -// E-MAIL : hecht@ann.jussieu.fr -// - -/* - - This file is part of Freefem++ - - Freefem++ is free software; you can redistribute it and/or modify - it under the terms of the GNU Lesser General Public License as published by - the Free Software Foundation; either version 2.1 of the License, or - (at your option) any later version. - - Freefem++ is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU Lesser General Public License for more details. - - You should have received a copy of the GNU Lesser General Public License - along with Freefem++; if not, write to the Free Software - Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA - */ +#include "../shared/shared.h" +#include "../include/macros.h" +#include "../toolkits/toolkits.h" + #include #include @@ -553,4 +530,204 @@ } +void Triangles::WriteBamgMesh(BamgMesh* bamgmesh,BamgOpts* bamgopts){ + + int i,j; + int verbose; + + verbose=bamgopts->verbose; + Int4 *reft = new Int4[nbt]; + Int4 nbInT = ConsRefTriangle(reft); + + //Vertices + if(verbose>3) printf(" writing Vertices\n"); + bamgmesh->NumVertices=nbv; + xfree((void**)&bamgmesh->Vertices); + if (nbv){ + bamgmesh->Vertices=(double*)xmalloc(3*nbv*sizeof(double)); + for (i=0;iVertices[i*3+0]=vertices[i].r.x; + bamgmesh->Vertices[i*3+1]=vertices[i].r.y; + bamgmesh->Vertices[i*3+2]=vertices[i].ref(); + } + } + else{ + bamgmesh->Vertices=NULL; + } + + //Edges + if(verbose>3) printf(" writing Edges\n"); + bamgmesh->NumEdges=nbe; + xfree((void**)&bamgmesh->Edges); + if (nbe){ + bamgmesh->Edges=(double*)xmalloc(3*nbe*sizeof(double)); + for (i=0;iEdges[i*3+0]=Number(edges[i][0])+1; //back to M indexing + bamgmesh->Edges[i*3+1]=Number(edges[i][1])+1; //back to M indexing + bamgmesh->Edges[i*3+2]=edges[i].ref; + } + } + else{ + bamgmesh->Edges=NULL; + } + + //CrackedEdges + if(verbose>3) printf(" writing CrackedEdges\n"); + bamgmesh->NumCrackedEdges=NbCrackedEdges; + xfree((void**)&bamgmesh->CrackedEdges); + if (NbCrackedEdges){ + bamgmesh->CrackedEdges=(double*)xmalloc(2*NbCrackedEdges*sizeof(double)); + for (i=0;iCrackedEdges[i*2+0]=Number(CrackedEdges[i].a.edge)+1; //back to M indexing + bamgmesh->CrackedEdges[i*2+1]=Number(CrackedEdges[i].b.edge)+1; //back to M indexing + } + } + else{ + bamgmesh->CrackedEdges=NULL; + } + + //Triangles + if(verbose>3) printf(" writing Triangles\n"); + Int4 k=nbInT-NbOfQuad*2; + Int4 num=0; + bamgmesh->NumTriangles=k; + xfree((void**)&bamgmesh->Triangles); + if (k){ + bamgmesh->Triangles=(double*)xmalloc(4*k*sizeof(double)); + for (i=0;i=0 && !( t.Hidden(0) || t.Hidden(1) || t.Hidden(2) )){ + bamgmesh->Triangles[num*4+0]=Number(t[0])+1; //back to M indexing + bamgmesh->Triangles[num*4+1]=Number(t[1])+1; //back to M indexing + bamgmesh->Triangles[num*4+2]=Number(t[2])+1; //back to M indexing + bamgmesh->Triangles[num*4+3]=subdomains[reft[i]].ref; + num=num+1; + } + } + } + else{ + bamgmesh->Triangles=NULL; + } + + //Quadrilaterals + if(verbose>3) printf(" writing Quadrilaterals\n"); + bamgmesh->NumQuadrilaterals=NbOfQuad; + xfree((void**)&bamgmesh->Quadrilaterals); + if (NbOfQuad){ + bamgmesh->Quadrilaterals=(double*)xmalloc(5*NbOfQuad*sizeof(double)); + for (i=0;iQuadrilaterals[i*5+0]=Number(v0)+1; //back to M indexing + bamgmesh->Quadrilaterals[i*5+1]=Number(v1)+1; //back to M indexing + bamgmesh->Quadrilaterals[i*5+2]=Number(v2)+1; //back to M indexing + bamgmesh->Quadrilaterals[i*5+3]=Number(v3)+1; //back to M indexing + bamgmesh->Quadrilaterals[i*5+4]=subdomains[reft[i]].ref; + } + } + } + else{ + bamgmesh->Quadrilaterals=NULL; + } + + //SubDomains + if(verbose>3) printf(" writing SubDomains\n"); + bamgmesh->NumSubDomains=NbSubDomains; + xfree((void**)&bamgmesh->SubDomains); + if (NbSubDomains){ + bamgmesh->SubDomains=(double*)xmalloc(4*NbSubDomains*sizeof(double)); + for (i=0;iSubDomains[i*4+0]=3; + bamgmesh->SubDomains[i*4+1]=reft[Number(subdomains[i].head)]; + bamgmesh->SubDomains[i*4+2]=1; + bamgmesh->SubDomains[i*4+3]=subdomains[i].ref; + } + } + else{ + bamgmesh->SubDomains=NULL; + } + + //SubDomainsFromGeom + if(verbose>3) printf(" writing SubDomainsFromGeom\n"); + bamgmesh->NumSubDomainsFromGeom=Gh.NbSubDomains; + xfree((void**)&bamgmesh->SubDomainsFromGeom); + if (Gh.NbSubDomains){ + bamgmesh->SubDomainsFromGeom=(double*)xmalloc(4*Gh.NbSubDomains*sizeof(double)); + for (i=0;iSubDomainsFromGeom[i*4+0]=2; + bamgmesh->SubDomainsFromGeom[i*4+1]=Number(subdomains[i].edge)+1; //back to Matlab indexing + bamgmesh->SubDomainsFromGeom[i*4+2]=subdomains[i].sens; + bamgmesh->SubDomainsFromGeom[i*4+3]=Gh.subdomains[i].ref; + } + } + else{ + bamgmesh->SubDomainsFromGeom=NULL; + } + + //VerticesOnGeomVertex + if(verbose>3) printf(" writing VerticesOnGeometricVertex\n"); + bamgmesh->NumVerticesOnGeometricVertex=NbVerticesOnGeomVertex; + xfree((void**)&bamgmesh->VerticesOnGeometricVertex); + if (NbVerticesOnGeomVertex){ + bamgmesh->VerticesOnGeometricVertex=(double*)xmalloc(2*NbVerticesOnGeomVertex*sizeof(double)); + for (i=0;iVerticesOnGeometricVertex[i*2+0]=Number((Vertex*)v)+1; //back to Matlab indexing + bamgmesh->VerticesOnGeometricVertex[i*2+1]=Gh.Number(( GeometricalVertex*)v)+1; //back to Matlab indexing + } + } + else{ + bamgmesh->VerticesOnGeometricVertex=NULL; + } + + //VertexOnGeometricEdge + if(verbose>3) printf(" writing VerticesOnGeometricEdge\n"); + bamgmesh->NumVerticesOnGeometricEdge=NbVerticesOnGeomEdge; + xfree((void**)&bamgmesh->VerticesOnGeometricEdge); + if (NbVerticesOnGeomEdge){ + bamgmesh->VerticesOnGeometricEdge=(double*)xmalloc(3*NbVerticesOnGeomEdge*sizeof(double)); + for (i=0;iVerticesOnGeometricEdge[i*3+0]=Number((Vertex*)v)+1; //back to Matlab indexing + bamgmesh->VerticesOnGeometricEdge[i*3+1]=Gh.Number((const GeometricalEdge*)v)+1; //back to Matlab indexing + bamgmesh->VerticesOnGeometricEdge[i*3+2]=(double)v; + } + } + else{ + bamgmesh->VerticesOnGeometricEdge=NULL; + } + + //EdgesOnGeometricEdge + if(verbose>3) printf(" writing EdgesOnGeometricEdge\n"); + k=0; + for (i=0;iNumEdgesOnGeometricEdge=k; + xfree((void**)&bamgmesh->EdgesOnGeometricEdge); + if (k){ + bamgmesh->EdgesOnGeometricEdge=(double*)xmalloc(2*(int)k*sizeof(double)); + int count=0; + for (i=0;iEdgesOnGeometricEdge[count*2+0]=(double)i+1; //back to Matlab indexing + bamgmesh->EdgesOnGeometricEdge[count*2+1]=(double)Gh.Number(edges[i].on)+1; //back to Matlab indexing + count=count+1; + } + } + } + else{ + bamgmesh->EdgesOnGeometricEdge=NULL; + } +} + void Triangles::Write(const char * filename) { @@ -719,7 +896,7 @@ f << "\nEdgeOnGeometricEdge\n"<< k << endl; for (i0=0;i01, clockwise=-1; - bamg_geometry.SubDomain=[2 count+1 clockwise 1]; + bamg_geometry.SubDomains=[2 count+1 clockwise 1]; end count=count+nods; @@ -44,5 +44,5 @@ bamg_geometry.NumVertices=size(bamg_geometry.Vertices,1); bamg_geometry.NumEdges=size(bamg_geometry.Edges,1); - bamg_geometry.NumSubDomain=size(bamg_geometry.SubDomain,1); + bamg_geometry.NumSubDomains=size(bamg_geometry.SubDomains,1); end end Index: /issm/trunk/src/mex/Bamg/Bamg.cpp =================================================================== --- /issm/trunk/src/mex/Bamg/Bamg.cpp (revision 2793) +++ /issm/trunk/src/mex/Bamg/Bamg.cpp (revision 2794) @@ -26,6 +26,6 @@ double* hVerticesGeom=NULL; double MaximalAngleOfCorner; - int NumSubDomainGeom; - double* SubDomainGeom=NULL; + int NumSubDomainsGeom; + double* SubDomainsGeom=NULL; /*Options inputs*/ @@ -63,8 +63,8 @@ bamggeom.NumRequiredEdges=0; bamggeom.RequiredEdges=NULL; - FetchData(&NumSubDomainGeom,mxGetField(BAMGGEOMETRY,0,"NumSubDomain")); - bamggeom.NumSubDomain=NumSubDomainGeom; - FetchData(&SubDomainGeom,NULL,NULL,mxGetField(BAMGGEOMETRY,0,"SubDomain")); - bamggeom.SubDomain=SubDomainGeom; + FetchData(&NumSubDomainsGeom,mxGetField(BAMGGEOMETRY,0,"NumSubDomains")); + bamggeom.NumSubDomains=NumSubDomainsGeom; + FetchData(&SubDomainsGeom,NULL,NULL,mxGetField(BAMGGEOMETRY,0,"SubDomains")); + bamggeom.SubDomains=SubDomainsGeom; /*create bamg mesh input*/ @@ -80,12 +80,18 @@ bamgmesh.NumQuadrilaterals=0; bamgmesh.Quadrilaterals=NULL; - bamgmesh.NumVertexOnGeometricEdge=0; - bamgmesh.VertexOnGeometricEdge=NULL; - bamgmesh.NumEdgeOnGeometricEdge=0; - bamgmesh.EdgeOnGeometricEdge=NULL; + bamgmesh.NumVerticesOnGeometricVertex=0; + bamgmesh.VerticesOnGeometricVertex=NULL; + bamgmesh.NumVerticesOnGeometricEdge=0; + bamgmesh.VerticesOnGeometricEdge=NULL; + bamgmesh.NumEdgesOnGeometricEdge=0; + bamgmesh.EdgesOnGeometricEdge=NULL; bamgmesh.NumEdges=0; bamgmesh.Edges=NULL; - bamgmesh.NumSubDomain=0; - bamgmesh.SubDomain=NULL; + bamgmesh.NumCrackedEdges=0; + bamgmesh.CrackedEdges=NULL; + bamgmesh.NumSubDomains=0; + bamgmesh.SubDomains=NULL; + bamgmesh.NumSubDomainsFromGeom=0; + bamgmesh.SubDomainsFromGeom=NULL; /*create bamg options input*/ @@ -122,5 +128,5 @@ xfree((void**)&EdgesGeom); xfree((void**)&hVerticesGeom); - xfree((void**)&SubDomainGeom); + xfree((void**)&SubDomainsGeom); //xfree((void**)&TrianglesMesh); //xfree((void**)&VerticesMesh);