source: issm/trunk/src/c/shared/TriMesh/SplitMeshForRifts.cpp@ 15396

Last change on this file since 15396 was 15396, checked in by Mathieu Morlighem, 12 years ago

merged trunk-jpl and trunk for revision 15394
File size: 2.9 KB
Line 
1/*
2 * SplitMeshForRifts.c:
3 */
4#include "./trimesh.h"
5#include "../MemOps/MemOps.h"
6
7int SplitMeshForRifts(int* pnel,int** pindex,int* pnods,double** px,double** py,int* pnsegs,int** psegments,int** psegmentmarkerlist){
8
9 /*Some notes on dimensions:
10 index of size nelx3
11 x and y of size nodsx1
12 segments of size nsegsx3*/
13
14 int i,j,k,l;
15 int node;
16 int el;
17 int nriftsegs;
18 int* riftsegments=NULL;
19 int* flags=NULL;
20 int NumGridElementListOnOneSideOfRift;
21 int* GridElementListOnOneSideOfRift=NULL;
22
23 /*Recover input: */
24 int nel = *pnel;
25 int *index = *pindex;
26 int nods = *pnods;
27 double *x = *px;
28 double *y = *py;
29 int nsegs = *pnsegs;
30 int *segments = *psegments;
31 int *segmentmarkerlist = *psegmentmarkerlist;
32
33 /*Establish list of segments that belong to a rift: */
34 /*riftsegments of size nriftsegsx4 (4 for first element on segment,second element,first node and second snode)*/
35 RiftSegmentsFromSegments(&nriftsegs,&riftsegments,nel,index,nsegs,segments);
36
37 /*Go through all nodes of the rift segments, and start splitting the mesh: */
38 flags=xNewZeroInit<int>(nods); //to make sure we don't split the same nodes twice!
39 for (i=0;i<nriftsegs;i++){
40 for (j=0;j<2;j++){
41
42 node=riftsegments[4*i+j+2];
43 if(flags[node-1]){
44 /*This node was already split, skip:*/
45 continue;
46 }
47 else{
48 flags[node-1]=1;
49 }
50
51 if(IsGridOnRift(riftsegments,nriftsegs,node)){
52
53 DetermineGridElementListOnOneSideOfRift(&NumGridElementListOnOneSideOfRift,&GridElementListOnOneSideOfRift,i,nriftsegs,riftsegments,node,index,nel);
54
55 /*Summary: we have for node, a list of elements
56 * (GridElementListOnOneSideOfRift, of size
57 * NumGridElementListOnOneSideOfRift) that all contain node
58 *and that are on the same side of the rift. For all these
59 elements, we clone node into another node, and we swap all
60 instances of node in the triangulation *for those elements, to the
61 new node.*/
62
63 //create new node
64 x=xReNew<double>(x,nods,nods+1);
65 y=xReNew<double>(y,nods,nods+1);
66 x[nods]=x[node-1]; //matlab indexing
67 y[nods]=y[node-1]; //matlab indexing
68
69 //augment number of nodes
70 nods++;
71
72 //change elements owning this node
73 for (k=0;k<NumGridElementListOnOneSideOfRift;k++){
74 el=GridElementListOnOneSideOfRift[k];
75 for (l=0;l<3;l++){
76 if (index[3*el+l]==node) index[3*el+l]=nods; //again, matlab indexing.
77 }
78 }
79 }// if(IsGridOnRift(riftsegments,nriftsegs,node))
80 } //for(j=0;j<2;j++)
81 } //for (i=0;i<nriftsegs;i++)
82
83 /*update segments: they got modified completely by adding new nodes.*/
84 UpdateSegments(&segments,&segmentmarkerlist, &nsegs,index,x,y,riftsegments,nriftsegs,nods,nel);
85
86 /*Assign output pointers: */
87 *pnel=nel;
88 *pindex=index;
89 *pnods=nods;
90 *px=x;
91 *py=y;
92 *pnsegs=nsegs;
93 *psegments=segments;
94 *psegmentmarkerlist=segmentmarkerlist;
95 return 1;
96}
Note: See TracBrowser for help on using the repository browser.