Changeset 13634

Timestamp:

10/12/12 10:45:13 (12 years ago)

Author:

Mathieu Morlighem

Message:

NEW: added TriMeshProcessRiftsx, x layer of TriMeshProcessRifts for python port

Location:

issm/trunk-jpl/src

Files:

• c/Makefile.am (modified) (2 diffs)
• c/classes/DofIndexing.cpp (modified) (1 diff)
• c/classes/FemModel.h (modified) (1 diff)
• c/classes/RiftStruct.cpp (added)
• c/classes/RiftStruct.h (added)
• c/classes/Update.h (modified) (1 diff)
• c/classes/classes.h (modified) (1 diff)
• c/classes/objects/Vertex.h (modified) (1 diff)
• c/modules/TriMeshProcessRiftsx (added)
• c/modules/TriMeshProcessRiftsx/TriMeshProcessRiftsx.cpp (added)
• c/modules/TriMeshProcessRiftsx/TriMeshProcessRiftsx.h (added)
• c/modules/modules.h (modified) (1 diff)
• c/shared/TriMesh/TriMeshUtils.cpp (modified) (1 diff)
• c/shared/shared.h (modified) (1 diff)
• modules/TriMeshProcessRifts/TriMeshProcessRifts.cpp (modified) (8 diffs)

issm/trunk-jpl/src/c/Makefile.am ¶

@@ -360 +360 @@
                                         ./classes/objects/Options/GenericOption.h\
                                         ./classes/objects/Options/OptionUtilities.cpp\
-                                        ./classes/objects/Options/OptionUtilities.h
-
+                                        ./classes/objects/Options/OptionUtilities.h\
+                                        ./classes/RiftStruct.cpp\
+                                        ./classes/RiftStruct.h
 #}}}
 #DAKOTA sources  {{{
@@ -848 +849 @@
                         ./modules/TriMeshx/TriMeshx.h\
                         ./modules/TriMeshx/TriMeshx.cpp\
+                        ./modules/TriMeshProcessRiftsx/TriMeshProcessRiftsx.h\
+                        ./modules/TriMeshProcessRiftsx/TriMeshProcessRiftsx.cpp\
                         ./modules/PointCloudFindNeighborsx/PointCloudFindNeighborsx.h\
                         ./modules/PointCloudFindNeighborsx/PointCloudFindNeighborsx.cpp\

issm/trunk-jpl/src/c/classes/DofIndexing.cpp ¶

@@ -10 +10 @@
 
 #include <stdio.h>
+#include <string.h>
 #include "./classes.h"
-#include <string.h>
 #include "../EnumDefinitions/EnumDefinitions.h"
 #include "../shared/shared.h"

issm/trunk-jpl/src/c/classes/FemModel.h ¶

@@ -20 +20 @@
 class Materials;
 class Profiler;
+class Elements;
 /*}}}*/
 

issm/trunk-jpl/src/c/classes/Update.h ¶

@@ -7 +7 @@
 
 /*Headers:*/
-/*{{{*/
 #include "../shared/shared.h"
-/*}}}*/
+class IoModel;
 
 class Update{

issm/trunk-jpl/src/c/classes/classes.h ¶

@@ -32 +32 @@
 #include "./AdolcEdf.h"
 #include "./IssmComm.h"
+#include "./RiftStruct.h"
 
 #endif

issm/trunk-jpl/src/c/classes/objects/Vertex.h ¶

@@ -9 +9 @@
 /*{{{*/
 #include "../classes.h"
-class IoModel;
-template <class doubletype> class Vector;
-class Parameters;
 #include "../../shared/Exceptions/exceptions.h"
 #include "../../toolkits/toolkits.h"
 #include "../../include/include.h"
-
+template <class doubletype> class Vector;
+class Parameters;
+class IoModel;
 /*}}}*/
 

issm/trunk-jpl/src/c/modules/modules.h ¶

@@ -117 +117 @@
 #include "./TriaSearchx/TriaSearchx.h"
 #include "./TriMeshx/TriMeshx.h"
+#include "./TriMeshProcessRiftsx/TriMeshProcessRiftsx.h"
 #include "./ThicknessAbsMisfitx/ThicknessAbsMisfitx.h"
 #include "./ThicknessAbsGradientx/ThicknessAbsGradientx.h"

issm/trunk-jpl/src/c/shared/TriMesh/TriMeshUtils.cpp ¶

@@ -266 +266 @@
         /*Reallocate segments: */
         segments         =xReNew<double>(segments,         nsegs*3,(nsegs+nriftsegs)*3);
-        segmentmarkerlist=xReNew<double>(segmentmarkerlist,nsegs*3,(nsegs+nriftsegs)*3);
+        segmentmarkerlist=xReNew<double>(segmentmarkerlist,nsegs,(nsegs+nriftsegs));
 
         /*First, update the existing segments to the new nodes :*/

issm/trunk-jpl/src/c/shared/shared.h ¶

@@ -6 +6 @@
 #define _SHARED_H_
 
-#include "Alloc/alloc.h"
-#include "Alloc/xNewDelete.h"
-#include "Bamg/shared.h"
-#include "Elements/elements.h"
-#include "Exceptions/exceptions.h"
-#include "Exp/exp.h"
-#include "Matrix/matrix.h"
-#include "MemOps/xMemCpy.h"
-#include "Numerics/numerics.h"
-#include "Sorting/sorting.h"
-#include "String/sharedstring.h"
-#include "Threads/issm_threads.h"
-#include "TriMesh/trimesh.h"
-#include "Wrapper/wrappershared.h"
+#include "./Alloc/alloc.h"
+#include "./Alloc/xNewDelete.h"
+#include "./Bamg/shared.h"
+#include "./Elements/elements.h"
+#include "./Exceptions/exceptions.h"
+#include "./Exp/exp.h"
+#include "./Matrix/matrix.h"
+#include "./MemOps/xMemCpy.h"
+#include "./Numerics/numerics.h"
+#include "./Sorting/sorting.h"
+#include "./String/sharedstring.h"
+#include "./Threads/issm_threads.h"
+#include "./TriMesh/trimesh.h"
+#include "./Wrapper/wrappershared.h"
 
 #endif

issm/trunk-jpl/src/modules/TriMeshProcessRifts/TriMeshProcessRifts.cpp ¶

@@ -11 +11 @@
         _printLine_("      index2,x2,y2,segments2,segmentmarkers2,rifts2 is the resulting triangulation where rifts have been processed.");
 }/*}}}*/
-void TriMeshProcessRiftsx(double** pindex, int* pnel,double** px,double** py,int* pnods,double** psegments,double** psegmentmarkers,int *pnum_seg,
-                        int  *pnumrifts, int  **priftsnumsegments, double ***priftssegments,
-                        int  **priftsnumpairs, double ***priftspairs, double **priftstips,
-                        double ***priftspenaltypairs, int  **priftsnumpenaltypairs);
 WRAPPER(TriMeshProcessRifts){
 
         /* returned quantities: */
-        int      numrifts;
-        int     *riftsnumsegments     = NULL;
-        double **riftssegments        = NULL;
-        int     *riftsnumpairs        = NULL;
-        double **riftspairs           = NULL;
-        double  *riftstips            = NULL;
-        double **riftspenaltypairs    = NULL;
-        int     *riftsnumpenaltypairs = NULL;
+        RiftStruct *riftstruct = NULL;
 
         /* input: */
@@ -51 +40 @@
         /*call x layer*/
         TriMeshProcessRiftsx(&index,&nel,&x,&y,&nods,&segments,&segmentmarkers,&num_seg,
-                                &numrifts,&riftsnumsegments,&riftssegments,&riftsnumpairs,&riftspairs,&riftstips,&riftspenaltypairs,&riftsnumpenaltypairs);
+                                &riftstruct);
 
         /*Output : */
@@ -60 +49 @@
         WriteData(SEGMENTMARKERSOUT,segmentmarkers,num_seg,1);
 
-        mxArray *pmxa_array  = NULL;
-        mxArray *pmxa_array2 = NULL;
-        mxArray *pmxa_array3 = NULL;
+        mxArray    *pmxa_array    = NULL;
+        mxArray    *pmxa_array2   = NULL;
+        mxArray    *pmxa_array3   = NULL;
         const char *fnames[10];
         mwSize      ndim          = 2;
         mwSize      dimensions[2] = {1,1};
         double     *pair          = NULL;
-
 
         /*Create a structure rifts where if i is a rift number, we have the following fields rifts(i).segments and rifts(i).numsegs: */
@@ -81 +69 @@
         fnames[9] = "state";
 
-        dimensions[0]=numrifts;
-        pmxa_array=mxCreateStructArray( ndim,dimensions,10,fnames);
+        dimensions[0]=riftstruct->numrifts;
+        pmxa_array=mxCreateStructArray(ndim,dimensions,10,fnames);
 
-        for(int i=0;i<numrifts;i++){
+        for(int i=0;i<riftstruct->numrifts;i++){
 
                 /*Segments: */
-                WriteData(&pmxa_array3,riftssegments[i],riftsnumsegments[i],3);
+                WriteData(&pmxa_array3,riftstruct->riftssegments[i],riftstruct->riftsnumsegments[i],3);
                 mxSetField(pmxa_array,i,"segments",pmxa_array3);
-                mxSetField(pmxa_array,i,"numsegs",mxCreateDoubleScalar((double)riftsnumsegments[i]));
+                mxSetField(pmxa_array,i,"numsegs",mxCreateDoubleScalar((double)riftstruct->riftsnumsegments[i]));
 
                 /*Element pairs: */
-                WriteData(&pmxa_array3,riftspairs[i],riftsnumpairs[i],2);
+                WriteData(&pmxa_array3,riftstruct->riftspairs[i],riftstruct->riftsnumpairs[i],2);
                 mxSetField(pmxa_array,i,"pairs",pmxa_array3);
 
@@ -99 +87 @@
                 mxSetM(pmxa_array2,1);
                 pair=(double*)mxMalloc(2*sizeof(double));
-                pair[0]=*(riftstips+2*i+0);
-                pair[1]=*(riftstips+2*i+1);
+                pair[0]=*(riftstruct->riftstips+2*i+0);
+                pair[1]=*(riftstruct->riftstips+2*i+1);
                 mxSetN(pmxa_array2,2);
                 mxSetPr(pmxa_array2,pair);
@@ -106 +94 @@
 
                 /*Penalty pairs: */
-                WriteData(&pmxa_array3,riftspenaltypairs[i],riftsnumpenaltypairs[i],7);
+                WriteData(&pmxa_array3,riftstruct->riftspenaltypairs[i],riftstruct->riftsnumpenaltypairs[i],7);
                 mxSetField(pmxa_array,i,"penaltypairs",pmxa_array3);
 
@@ -116 +104 @@
 
                 /*State: */
-                double *state=(double*)mxMalloc(riftsnumpenaltypairs[i]*sizeof(double));
-                for(int j=0;j<riftsnumpenaltypairs[i];j++) state[j]=FreeEnum;
                 pmxa_array2= mxCreateDoubleMatrix(0,0,mxREAL);
                 mxSetM(pmxa_array2,1);
-                mxSetN(pmxa_array2,riftsnumpenaltypairs[i]);
-                mxSetPr(pmxa_array2,state);
+                mxSetN(pmxa_array2,riftstruct->riftsnumpenaltypairs[i]);
+                mxSetPr(pmxa_array2,riftstruct->state[i]);
                 mexCallMATLAB( 1, &pmxa_array3, 1, &pmxa_array2, "transpose");
 
@@ -132 +118 @@
         MODULEEND();
 }
-void TriMeshProcessRiftsx(double** pindex, int* pnel,double** px,double** py,int* pnods,double** psegments,double** psegmentmarkers,int *pnum_seg,
-                        int  *pnumrifts, int  **priftsnumsegments, double ***priftssegments,
-                        int  **priftsnumpairs, double ***priftspairs, double **priftstips,
-                        double ***priftspenaltypairs, int  **priftsnumpenaltypairs){
-
-        /*Output*/
-        int      numrifts,numrifts0;
-        int     *riftsnumsegments     = NULL;
-        double **riftssegments        = NULL;
-        int     *riftsnumpairs        = NULL;
-        double **riftspairs           = NULL;
-        double  *riftstips            = NULL;
-        double **riftspenaltypairs    = NULL;
-        int     *riftsnumpenaltypairs = NULL;
-
-        /*Recover initial mesh*/
-        int     nel            = *pnel;
-        double *index          = *pindex;
-        double *x              = *px;
-        double *y              = *py;
-        int     nods           = *pnods;
-        double *segments       = *psegments;
-        double *segmentmarkers = *psegmentmarkers;
-        int     num_seg        = *pnum_seg;
-
-        /*Intermediary*/
-        int     riftflag;
-
-        /*First, do some fixing on the existing mesh: we do not want any element belonging entirely to the segment list (ie: 
-         *all the nodes of this element belong to the segments (tends to happen when there are corners: */
-        RemoveCornersFromRifts(&index,&nel,&x,&y,&nods,segments,segmentmarkers,num_seg);
-
-        /*Figure out if we have rifts, and how many: */
-        IsRiftPresent(&riftflag,&numrifts0,segmentmarkers,num_seg);
-
-        if(!riftflag) _error_("No rift present in mesh");
-
-        /*Split mesh*/
-        SplitMeshForRifts(&nel,&index,&nods,&x,&y,&num_seg,&segments,&segmentmarkers);
-
-        /*Order segments so that their normals point outside the domain: */
-        OrderSegments(&segments,num_seg, index,nel);
-
-        /*We do not want to output segments mixed with rift segments: wring out the rifts from the segments, using the 
-         *segmentmarkerlist:*/
-        SplitRiftSegments(&segments,&segmentmarkers,&num_seg,&numrifts,&riftsnumsegments,&riftssegments,numrifts0,nods,nel);
-
-        /*Using rift segments, associate rift faces in pairs, each pair face representing opposite flanks of the rifts facing one another directly: */
-        PairRiftElements(&riftsnumpairs,&riftspairs,numrifts,riftsnumsegments,riftssegments,x,y);
-
-        /*Order rifts so that they start from one tip, go to the other tip, and back: */
-        OrderRifts(&riftstips,riftssegments,riftspairs,numrifts,riftsnumsegments,x,y,nods,nel);
-
-        /*Create penalty pairs, used by Imp: */
-        PenaltyPairs(&riftspenaltypairs,&riftsnumpenaltypairs,numrifts,riftssegments,riftsnumsegments,riftspairs,riftstips,x,y);
-
-        /*Assign output pointers for mesh*/
-        *pnel            = nel;
-        *pindex          = index;
-        *px              = x;
-        *py              = y;
-        *pnods           = nods;
-        *psegments       = segments;
-        *psegmentmarkers = segmentmarkers;
-        *pnum_seg        = num_seg;
-
-        /*Assign output pointers for rifts*/
-        *pnumrifts             = numrifts;
-        *priftsnumsegments     = riftsnumsegments;
-        *priftssegments        = riftssegments;
-        *priftsnumpairs        = riftsnumpairs;
-        *priftspairs           = riftspairs;
-        *priftstips            = riftstips;
-        *priftspenaltypairs    = riftspenaltypairs;
-        *priftsnumpenaltypairs = riftsnumpenaltypairs;
-}