![(please configure the [header_logo] section in trac.ini)](/trac/issm/chrome/common/trac_banner.png)
Changeset 5531
- Timestamp:
- 08/24/10 08:40:17 (15 years ago)
- File:
-
- 1 edited
-
issm/trunk/src/c/objects/Bamg/Mesh.cpp (modified) (16 diffs)
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/c/objects/Bamg/Mesh.cpp
r5489 r5531 9 9 10 10 static const Direction NoDirOfSearch=Direction(); 11 long NbUnSwap =0;12 11 13 12 /*Constructors/Destructors*/ … … 4895 4894 int verbose=bamgopts->verbose; 4896 4895 4897 //initialize Mesh4898 nbv=0;4899 NbVerticesOnGeomVertex=0;4900 NbVerticesOnGeomEdge=0;4901 4902 4896 //build background mesh flag (1 if background, else 0) 4903 4897 bool background=(&BTh != this); … … 4912 4906 VerticesOnGeomVertex = new VertexOnGeom[NbVerticesOnGeomVertex]; 4913 4907 if(NbVerticesOnGeomVertex >= maxnbv) ISSMERROR("too many vertices on geometry: %i >= %i",NbVerticesOnGeomVertex,maxnbv); 4908 ISSMASSERT(nbv==0); 4914 4909 //Build VerticesOnGeomVertex 4915 nbv=0;4916 4910 for (i=0;i<Gh.nbv;i++){ 4917 4911 /* Add vertex only if required*/ … … 4941 4935 * to the one present in the geometry provided. 4942 4936 * We proceed in 2 steps 4943 * -step 1: we count all the edges4944 * we allocate the number of edges at the end of step 14945 * -step 2: the edges are created */4937 * -step 0: we count all the edges 4938 * we allocate the number of edges at the end of step 0 4939 * -step 1: the edges are created */ 4946 4940 for (int step=0;step<2;step++){ 4947 4941 … … 4962 4956 for(int j=0;j<2;j++) { 4963 4957 4964 /*Take only required vertices (corner )*/4958 /*Take only required vertices (corner->beginning of a new curve)*/ 4965 4959 if (!ei.Mark() && ei[j].Required()){ 4966 4960 … … 4970 4964 4971 4965 /*If Edge is required (do that only once for the 2 vertices)*/ 4972 if(ei.Required() && j==0){ 4973 //do not create internal points if required (take it as is) 4974 if(step==0) nbe++; 4975 else{ 4976 e=&ei; 4977 a=ei(0); 4978 b=ei(1); 4979 4980 //check that edges has been allocated 4981 ISSMASSERT(edges); 4982 edges[nbe].v[0]=a->MeshVertexHook; 4983 edges[nbe].v[1]=b->MeshVertexHook;; 4984 edges[nbe].ReferenceNumber = e->ReferenceNumber; 4985 edges[nbe].GeometricalEdgeHook = e; 4986 edges[nbe].adj[0] = 0; 4987 edges[nbe].adj[1] = 0; 4988 nbe++; 4966 if(ei.Required()){ 4967 if (j==0){ 4968 //do not create internal points if required (take it as is) 4969 if(step==0) nbe++; 4970 else{ 4971 e=&ei; 4972 a=ei(0); 4973 b=ei(1); 4974 4975 //check that edges has been allocated 4976 ISSMASSERT(edges); 4977 edges[nbe].v[0]=a->MeshVertexHook; 4978 edges[nbe].v[1]=b->MeshVertexHook;; 4979 edges[nbe].ReferenceNumber = e->ReferenceNumber; 4980 edges[nbe].GeometricalEdgeHook = e; 4981 edges[nbe].adj[0] = 0; 4982 edges[nbe].adj[1] = 0; 4983 nbe++; 4984 } 4989 4985 } 4990 4986 } … … 4993 4989 else { 4994 4990 for (int kstep=0;kstep<=step;kstep++){ 4995 //kstep=0, do nothing 4996 //kstep=1, compute the length of the curve 4997 //kstep=2 create the points and edge 4991 //kstep=0, compute number of edges (discretize curve) 4992 //kstep=1 create the points and edge 4998 4993 PreviousNewEdge=0; 4999 4994 NbNewPoints=0; … … 5001 4996 if (nbvend>=maxnbv) ISSMERROR("maximum number of vertices too low! Check the domain outline or increase maxnbv"); 5002 4997 lcurve =0; 5003 s = lstep; 4998 s = lstep; //-1 initially, then length of each sub edge 5004 4999 5005 5000 /*reminder: i = edge number, j=[0;1] vertex index in edge*/ … … 5010 5005 e->SetMark(); // Mark edge 5011 5006 5012 /* If we have a Background mesh, we can use it to discretizethe curve*/5007 /*Loop until we reach the end of the curve*/ 5013 5008 for(;;){ 5014 5009 k = 1-k; // other vertx index of the curve … … 5016 5011 AB= b->r - a->r; // AB = vector of the current edge 5017 5012 Metric MA = background ? BTh.MetricAt(a->r) :a->m ; //Get metric associated to A 5018 Metric MB = background ? BTh.MetricAt(b->r) :b->m ; //Get metric associated to B5019 double ledge = (MA(AB) + MB(AB))/2; //Get edge length in metric5013 Metric MB = background ? BTh.MetricAt(b->r) :b->m ; //Get metric associated to B 5014 double ledge = (MA(AB) + MB(AB))/2; //Get edge length in metric 5020 5015 5021 5016 /* We are now creating the mesh edges from the geometrical edge selected above. … … 5036 5031 //else, divide the edge 5037 5032 else { 5038 //compute number of subedges (division of the edge) 5033 //compute number of subedges (division of the edge), Maximum is 10 5039 5034 NbSubEdge = Min( MaxSubEdge, (int) (ledge +0.5)); 5040 //A and B are the position of points on the edge 5035 /*Now, we are going to divide the edge according to the metric. 5036 * Get segment by sement along the edge. 5037 * Build lSubEdge, which holds the distance between the first vertex 5038 * of the edge and the next point on the edge according to the 5039 * discretization (each SubEdge is AB)*/ 5041 5040 R2 A,B; 5042 5041 A=a->r; … … 5049 5048 MBs= background ? BTh.MetricAt(B) : Metric(1-x,MA,x,MB); 5050 5049 AB = A-B; 5051 lSubEdge[kk]= 5050 lSubEdge[kk]=(ledge+=(MAs(AB)+MBs(AB))/2); 5052 5051 } 5053 5052 } … … 5056 5055 5057 5056 /*Now, create corresponding points*/ 5058 while (lcurve<=s&& s<=lcurveb && nbv<nbvend){5057 while(s>=lcurve && s<=lcurveb && nbv<nbvend){ 5059 5058 5060 5059 double ss = s-lcurve; … … 5095 5094 va = vb; 5096 5095 } 5096 5097 /*We just added one edge to the curve: Go to the next one*/ 5097 5098 lcurve = lcurveb; 5098 5099 e->SetMark(); 5099 5100 a=b; 5100 if (b->Required() ) break; 5101 5102 /*If b is required, we are on a new curve->break*/ 5103 if (b->Required()) break; 5101 5104 int kprev=k; 5102 5105 k = e->AdjVertexIndex[kprev];// next vertices … … 5105 5108 }// for(;;) 5106 5109 vb = b->MeshVertexHook; 5110 5111 /*Number of edges in the last disretized curve*/ 5107 5112 NbEdgeCurve = Max((long) (lcurve +0.5), (long) 1); 5113 /*Number of internal vertices in the last disretized curve*/ 5108 5114 NbNewPoints = NbEdgeCurve-1; 5109 5115 if(!kstep){ … … 5112 5118 } 5113 5119 nbvend=nbv+NbNewPoints; 5114 lstep = lcurve / NbEdgeCurve; 5120 lstep = lcurve / NbEdgeCurve; //approximately one 5115 5121 }// end of curve -- 5116 5122 if (edges) { // last edges of the curves
Note:
See TracChangeset
for help on using the changeset viewer.
