Index: /issm/trunk/configs/linux64/linux64.sh.petsc2 =================================================================== --- /issm/trunk/configs/linux64/linux64.sh.petsc2 (revision 4772) +++ /issm/trunk/configs/linux64/linux64.sh.petsc2 (revision 4773) @@ -1,3 +1,3 @@ #!/bin/sh -./configure --prefix=$ISSM_DIR --with-matlab-dir=$MATLAB_DIR --with-triangle-dir=$ISSM_DIR/externalpackages/triangle/install --with-metis-dir=$ISSM_DIR/externalpackages/metis/install --with-petsc-dir=$ISSM_DIR/externalpackages/petsc/install --with-mpi-include=$ISSM_DIR/externalpackages/mpich2/install/include --with-mpi-lib="-L$ISSM_DIR/externalpackages/mpich2/install/lib/ -lmpich" --with-petsc-arch=$ISSM_ARCH --with-dakota-dir=$ISSM_DIR/externalpackages/dakota/install --with-blas-lapack-dir=$ISSM_DIR/externalpackages/petsc/install/externalpackages/fblaslapack/$ISSM_ARCH --with-plapack-lib="-L$ISSM_DIR/externalpackages/petsc/install/externalpackages/PLAPACKR32-hg/$ISSM_ARCH -lPLAPACK" --with-plapack-include="-I$ISSM_DIR/externalpackages/petsc/install/externalpackages/PLAPACKR32-hg/$ISSM_ARCH/INCLUDE" --with-blacs-dir=$ISSM_DIR/externalpackages/petsc/install/externalpackages/blacs-dev/$ISSM_ARCH --with-scalapack-dir=$ISSM_DIR/externalpackages/petsc/install/externalpackages/SCALAPACK/$ISSM_ARCH --with-mumps-dir=$ISSM_DIR/externalpackages/petsc/install/externalpackages/MUMPS_4.6.3/$ISSM_ARCH --with-fortran-lib="-L/usr/lib/gcc/x86_64-redhat-linux/4.1.1/ -lgfortran" --with-graphics-lib=/usr/lib64/libX11.so --with-cxxoptflags="-march=opteron -O2" --with-numthreads=32 --with-petsc-version=2 --with-chaco-dir=$ISSM_DIR/externalpackages/chaco/install +./configure --prefix=$ISSM_DIR --with-matlab-dir=$MATLAB_DIR --with-triangle-dir=$ISSM_DIR/externalpackages/triangle/install --with-metis-dir=$ISSM_DIR/externalpackages/metis/install --with-petsc-dir=$ISSM_DIR/externalpackages/petsc/install --with-mpi-include=$ISSM_DIR/externalpackages/mpich2/install/include --with-mpi-lib="-L$ISSM_DIR/externalpackages/mpich2/install/lib/ -lmpich" --with-petsc-arch=$ISSM_ARCH --with-dakota-dir=$ISSM_DIR/externalpackages/dakota/install --with-blas-lapack-dir=$ISSM_DIR/externalpackages/petsc/install/externalpackages/fblaslapack/$ISSM_ARCH --with-plapack-lib="-L$ISSM_DIR/externalpackages/petsc/install/externalpackages/PLAPACKR32-hg/$ISSM_ARCH -lPLAPACK" --with-plapack-include="-I$ISSM_DIR/externalpackages/petsc/install/externalpackages/PLAPACKR32-hg/$ISSM_ARCH/INCLUDE" --with-blacs-dir=$ISSM_DIR/externalpackages/petsc/install/externalpackages/blacs-dev/$ISSM_ARCH --with-scalapack-dir=$ISSM_DIR/externalpackages/petsc/install/externalpackages/SCALAPACK/$ISSM_ARCH --with-mumps-dir=$ISSM_DIR/externalpackages/petsc/install/externalpackages/MUMPS_4.6.3/$ISSM_ARCH --with-fortran-lib="-L/usr/lib/gcc/x86_64-redhat-linux/4.1.1/ -lgfortran" --with-graphics-lib=/usr/lib64/libX11.so --with-cxxoptflags="-march=opteron -O2" --with-numthreads=32 --with-petsc-version=2 --with-scotch-dir=$ISSM_DIR/externalpackages/scotch/install --with-chaco-dir=$ISSM_DIR/externalpackages/chaco/install Index: /issm/trunk/etc/cluster.rc =================================================================== --- /issm/trunk/etc/cluster.rc (revision 4772) +++ /issm/trunk/etc/cluster.rc (revision 4773) @@ -54,5 +54,5 @@ cluster_name=larsen cluster_codepath=/u/astrid-r1b/larour/issm/trunk/bin -cluster_executionpath=/u/astrid-r1b/larour/Testing/Execution +cluster_executionpath=/u/astrid-r1b/larour/issm/execution cluster_login=larour cluster_port=0 Index: /issm/trunk/src/c/Makefile.am =================================================================== --- /issm/trunk/src/c/Makefile.am (revision 4772) +++ /issm/trunk/src/c/Makefile.am (revision 4773) @@ -94,4 +94,6 @@ ./objects/Node.h\ ./objects/Node.cpp\ + ./objects/Segment.h\ + ./objects/Segment.cpp\ ./objects/Vertex.h\ ./objects/Vertex.cpp\ @@ -499,4 +501,9 @@ ./modules/MeshPartitionx/MeshPartitionx.cpp\ ./modules/MeshPartitionx/MeshPartitionx.h\ + ./modules/MeshProfileIntersectionx/MeshProfileIntersectionx.cpp\ + ./modules/MeshProfileIntersectionx/MeshProfileIntersectionx.h\ + ./modules/MeshProfileIntersectionx/MeshSegmentsIntersection.cpp\ + ./modules/MeshProfileIntersectionx/ElementSegmentsIntersection.cpp\ + ./modules/MeshProfileIntersectionx/ElementSegment.cpp\ ./modules/ContourToMeshx/ContourToMeshx.cpp\ ./modules/ContourToMeshx/ContourToMeshxt.cpp\ @@ -653,4 +660,6 @@ ./objects/Node.h\ ./objects/Node.cpp\ + ./objects/Segment.h\ + ./objects/Segment.cpp\ ./objects/Vertex.h\ ./objects/Vertex.cpp\ @@ -1051,4 +1060,9 @@ ./modules/MeshPartitionx/MeshPartitionx.cpp\ ./modules/MeshPartitionx/MeshPartitionx.h\ + ./modules/MeshProfileIntersectionx/MeshProfileIntersectionx.cpp\ + ./modules/MeshProfileIntersectionx/MeshProfileIntersectionx.h\ + ./modules/MeshProfileIntersectionx/MeshSegmentsIntersection.cpp\ + ./modules/MeshProfileIntersectionx/ElementSegmentsIntersection.cpp\ + ./modules/MeshProfileIntersectionx/ElementSegment.cpp\ ./modules/ContourToMeshx/ContourToMeshx.cpp\ ./modules/ContourToMeshx/ContourToMeshxt.cpp\ @@ -1146,5 +1160,5 @@ endif -LDADD = ./libpISSM.a $(PETSCLIB) $(FLIBS) $(PLAPACKLIB) $(MUMPSLIB) $(SCALAPACKLIB) $(BLACSLIB) $(DAKOTALIB) $(METISLIB) $(CHACOLIB) $(SCOTCHLIB) $(BLASLAPACKLIB) $(MKLLIB) $(MPILIB) $(MATHLIB) $(FORTRANLIB) $(GRAPHICSLIB) libOverload.a $(MULTITHREADINGLIB) +LDADD = ./libpISSM.a $(PETSCLIB) $(FLIBS) $(PLAPACKLIB) $(MUMPSLIB) $(SCALAPACKLIB) $(BLACSLIB) $(DAKOTALIB) $(METISLIB) $(CHACOLIB) $(SCOTCHLIB) $(BLASLAPACKLIB) $(MKLLIB) $(MPILIB) $(MATHLIB) $(FORTRANLIB) $(GRAPHICSLIB) libOverload.a $(MULTITHREADINGLIB) issm_exe_SOURCES = solutions/issm.cpp Index: /issm/trunk/src/c/modules/DakotaResponsesx/DakotaResponsesx.cpp =================================================================== --- /issm/trunk/src/c/modules/DakotaResponsesx/DakotaResponsesx.cpp (revision 4772) +++ /issm/trunk/src/c/modules/DakotaResponsesx/DakotaResponsesx.cpp (revision 4773) @@ -34,53 +34,56 @@ /*Compute response for this response_descriptor:*/ - if(strcmp(response_descriptor,"min_vel")==0){ + if(strcmp(response_descriptor,"MinVel")==0){ MinVelx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"max_vel")==0){ + else if(strcmp(response_descriptor,"MaxVel")==0){ MaxVelx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"min_vx")==0){ + else if(strcmp(response_descriptor,"MinVx")==0){ MinVxx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"max_vx")==0){ + else if(strcmp(response_descriptor,"MaxVx")==0){ MaxVxx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"max_abs_vx")==0){ + else if(strcmp(response_descriptor,"MaxAbsVx")==0){ MaxAbsVxx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"min_vy")==0){ + else if(strcmp(response_descriptor,"MinVy")==0){ MinVyx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"max_vy")==0){ + else if(strcmp(response_descriptor,"MaxVy")==0){ MaxVyx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"max_abs_vy")==0){ + else if(strcmp(response_descriptor,"MaxAbsVy")==0){ MaxAbsVyx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"min_vz")==0){ + else if(strcmp(response_descriptor,"MinVz")==0){ MinVzx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"max_vz")==0){ + else if(strcmp(response_descriptor,"MaxVz")==0){ MaxVzx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"max_abs_vz")==0){ + else if(strcmp(response_descriptor,"MaxAbsVz")==0){ MaxAbsVzx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"misfit")==0){ + else if(strcmp(response_descriptor,"Misfit")==0){ Misfitx( &femmodel_response, elements,nodes, vertices, loads, materials, parameters); } - else if(strcmp(response_descriptor,"mass_flux")==0){ + else if(strlen(response_descriptor)>=8){ - double* segments=NULL; - int num_segments; - - /*retrieve qmu_mass_flux_segments: */ - parameters->FindParam(&segments,&num_segments,&dummy,QmuMassFluxSegmentsEnum); + if(strncmp(response_descriptor,"MassFlux",8)==0){ - /*call mass flux module: */ - MassFluxx(&femmodel_response,elements,nodes,vertices,loads,materials,parameters,segments,num_segments); + double* segments=NULL; + int num_segments; - /*Free ressources:*/ - xfree((void**)&segments); + /*retrieve qmu_mass_flux_segments: */ + parameters->FindParam(&segments,&num_segments,&dummy,QmuMassFluxSegmentsEnum); + + /*call mass flux module: */ + MassFluxx(&femmodel_response,elements,nodes,vertices,loads,materials,parameters,segments,num_segments); + + /*Free ressources:*/ + xfree((void**)&segments); + } } Index: /issm/trunk/src/c/modules/MeshProfileIntersectionx/ElementSegment.cpp =================================================================== --- /issm/trunk/src/c/modules/MeshProfileIntersectionx/ElementSegment.cpp (revision 4773) +++ /issm/trunk/src/c/modules/MeshProfileIntersectionx/ElementSegment.cpp (revision 4773) @@ -0,0 +1,14 @@ +/*! \file ElementSegment.cpp + */ + +#include "./MeshProfileIntersectionx.h" + +void ElementSegment(DataSet* segments_dataset,double* xgrids,double* ygrids,double* xsegment,double* ysegment){ + + + /*We have a tria element (xgrids,ygrids) and a segment (xsegment,ysegment). Find whether they intersect. + * If they do, create a Segment object with the intersection, and add to segments_dataset dataset: */ + + ISSMERROR("not supported yet!"); + +} Index: /issm/trunk/src/c/modules/MeshProfileIntersectionx/ElementSegmentsIntersection.cpp =================================================================== --- /issm/trunk/src/c/modules/MeshProfileIntersectionx/ElementSegmentsIntersection.cpp (revision 4773) +++ /issm/trunk/src/c/modules/MeshProfileIntersectionx/ElementSegmentsIntersection.cpp (revision 4773) @@ -0,0 +1,23 @@ +/*! \file ElementSegmentsIntersection.cpp + */ + +#include "./MeshProfileIntersectionx.h" + +void ElementSegmentsIntersection(DataSet* segments_dataset,double* xgrids,double* ygrids,double* xc,double* yc,int numgrids){ + + int i; + double xsegment[2]; + double ysegment[2]; + + /*Loop through contour: */ + for(i=0;inods; + xc=contouri->x; + yc=contouri->y; + + /*determine segmentsi and numsegsi for this contour and the mesh intersection: */ + MeshSegmentsIntersection(&segmentsi,&numsegsi,index,x,y,nel,nods,xc,yc,numgrids); + + /*save segmentsi: */ + allsegments[i]=segmentsi; + allnumsegs[i]=numsegsi; + } + + /*total number of segments? */ + numsegs=0; + for(i=0;iSize(); + segments=(double*)xmalloc(5*numsegs*sizeof(double)); + for(i=0;iGetObjectByOffset(i); + + /*x1,y1,x2,y2 then element_id: */ + *(segments+5*i+0)=segment->x1; + *(segments+5*i+1)=segment->y1; + *(segments+5*i+2)=segment->x2; + *(segments+5*i+3)=segment->y2; + *(segments+5*i+4)=(double)segment->eid; + } + + /*Free ressources:*/ + delete segments_dataset; + + /*Assign output pointers:*/ + *psegments=segments; + *pnumsegs=numsegs; +} Index: /issm/trunk/src/c/modules/modules.h =================================================================== --- /issm/trunk/src/c/modules/modules.h (revision 4772) +++ /issm/trunk/src/c/modules/modules.h (revision 4773) @@ -56,4 +56,5 @@ #include "./Mergesolutionfromftogx/Mergesolutionfromftogx.h" #include "./MeshPartitionx/MeshPartitionx.h" +#include "./MeshProfileIntersectionx/MeshProfileIntersectionx.h" #include "./MinVelx/MinVelx.h" #include "./MinVxx/MinVxx.h" Index: /issm/trunk/src/c/objects/Segment.cpp =================================================================== --- /issm/trunk/src/c/objects/Segment.cpp (revision 4773) +++ /issm/trunk/src/c/objects/Segment.cpp (revision 4773) @@ -0,0 +1,102 @@ +/*!\file Segment.c + * \brief: implementation of the Segment object + */ + +/*Include files: {{{1*/ +#ifdef HAVE_CONFIG_H + #include "config.h" +#else +#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!" +#endif + +#include "stdio.h" +#include +#include "./objects.h" +#include "../Container/Container.h" +#include "../EnumDefinitions/EnumDefinitions.h" +#include "../shared/shared.h" +#include "../include/include.h" +/*}}}*/ + +/*Segment constructors and destructors:*/ +/*FUNCTION Segment::Segment() default constructor {{{1*/ +Segment::Segment(){ + this->eid=UNDEF; + this->x1=UNDEF; + this->y1=UNDEF; + this->x2=UNDEF; + this->y2=UNDEF; +} +/*}}}*/ +/*FUNCTION Segment::Segment(int eid, double x1,double y1,double x2, double y2){{{1*/ +Segment::Segment(int segment_eid, double segment_x1,double segment_y1,double segment_x2, double segment_y2){ + + this->eid=segment_eid; + this->x1=segment_x1; + this->y1=segment_y1; + this->x2=segment_x2; + this->y2=segment_y2; + +} +/*}}}*/ +/*FUNCTION Segment::~Segment(){{{1*/ +Segment::~Segment(){ +} +/*}}}*/ + +/*Object virtual functions definitions:*/ +/*FUNCTION Segment::Echo{{{1*/ +void Segment::Echo(void){ + + printf("Segment:\n"); + printf(" eid: %i\n",eid); + printf(" grid 1: %g|%g\n",this->x1,this->y1); + printf(" grid 2: %g|%g\n",this->x2,this->y2); + +} +/*}}}*/ +/*FUNCTION Segment::DeepEcho{{{1*/ +void Segment::DeepEcho(void){ + this->Echo(); +} +/*}}}*/ +/*FUNCTION Segment::Id{{{1*/ +int Segment::Id(void){ return eid; } +/*}}}*/ +/*FUNCTION Segment::MyRank{{{1*/ +int Segment::MyRank(void){ + extern int my_rank; + + return my_rank; +} +/*}}}*/ +/*FUNCTION Segment::Marshall{{{1*/ +void Segment::Marshall(char** pmarshalled_dataset){ + + ISSMERROR(" not supported yet!"); +} +/*}}}*/ +/*FUNCTION Segment::MarshallSize{{{1*/ +int Segment::MarshallSize(){ + ISSMERROR(" not supported yet!"); +} +/*}}}*/ +/*FUNCTION Segment::Demarshall{{{1*/ +void Segment::Demarshall(char** pmarshalled_dataset){ + ISSMERROR(" not supported yet!"); +} +/*}}}*/ +/*FUNCTION Segment::Enum{{{1*/ +int Segment::Enum(void){ + + return SegmentEnum; + +} +/*}}}*/ +/*FUNCTION Segment::copy {{{1*/ +Object* Segment::copy() { + return new Segment(this->eid,this->x1,this->y1,this->x2,this->y2); + +} + +/*}}}*/ Index: /issm/trunk/src/c/objects/Segment.h =================================================================== --- /issm/trunk/src/c/objects/Segment.h (revision 4773) +++ /issm/trunk/src/c/objects/Segment.h (revision 4773) @@ -0,0 +1,40 @@ +/*!\file Segment.h + * \brief: header file for node object + */ + +#ifndef _SEGMENT_H_ +#define _SEGMENT_H_ + +/*Headers:*/ +/*{{{1*/ +#include "./Object.h" +/*}}}*/ + +class Segment: public Object{ + + public: + int eid; + double x1; + double y1; + double x2; + double y2; + + /*Segment constructors, destructors {{{1*/ + Segment(); + Segment(int eid,double x1,double y1, double x2, double y2); + ~Segment(); + /*}}}*/ + /*Object virtual functions definitions:{{{1 */ + void Echo(); + void DeepEcho(); + int Id(); + int MyRank(); + void Marshall(char** pmarshalled_dataset); + int MarshallSize(); + void Demarshall(char** pmarshalled_dataset); + int Enum(); + Object* copy(); + /*}}}*/ +}; + +#endif /* _SEGMENT_H_ */ Index: /issm/trunk/src/c/objects/objects.h =================================================================== --- /issm/trunk/src/c/objects/objects.h (revision 4772) +++ /issm/trunk/src/c/objects/objects.h (revision 4773) @@ -14,4 +14,5 @@ #include "./Vertex.h" #include "./Node.h" +#include "./Segment.h" #include "./NodeSets.h" #include "./IoModel.h" Index: /issm/trunk/src/m/classes/@model/model.m =================================================================== --- /issm/trunk/src/m/classes/@model/model.m (revision 4772) +++ /issm/trunk/src/m/classes/@model/model.m (revision 4773) @@ -300,7 +300,8 @@ md.variabledescriptors=NaN; md.responsedescriptors=NaN; - md.qmu_mass_flux_profile=NaN; + md.qmu_mass_flux_num_profiles=NaN; + md.qmu_mass_flux_profile_directory=NaN; + md.qmu_mass_flux_profiles=NaN; md.qmu_mass_flux_segments=NaN; - md.qmu_mass_flux_num_segments=NaN; %partitioner: Index: /issm/trunk/src/m/qmu/MassFluxProcessProfile.m =================================================================== --- /issm/trunk/src/m/qmu/MassFluxProcessProfile.m (revision 4772) +++ /issm/trunk/src/m/qmu/MassFluxProcessProfile.m (revision 4773) @@ -1,3 +1,3 @@ -function segments=MassFluxProcessProfile(md); +function segments=MassFluxProcessProfile(md,directory,profilename); %MASSFLUXPROCESSPROFILE: process an argus domain outlien profile into a list of segments. % @@ -9,5 +9,5 @@ %first read the profile points. -profile=expread(['../' md.qmu_mass_flux_profile],1); %go to ../ because the exp file is in the root directory. +profile=expread([directory '/' profilename]); %project this profile onto mesh. Index: /issm/trunk/src/m/qmu/process_qmu_response_data.m =================================================================== --- /issm/trunk/src/m/qmu/process_qmu_response_data.m (revision 4772) +++ /issm/trunk/src/m/qmu/process_qmu_response_data.m (revision 4773) @@ -7,53 +7,42 @@ +%preliminary data +process_mass_flux_profiles=0; + +%loop through response descriptors, and act accordingly for i=1:numel(md.responsedescriptors), - if strncmpi(md.responsedescriptors{i},'mass_flux',9), - indx=str2int(md.responsedescriptors{i}(10:end)); - if isempty(indx) || ~indx - indx=1; - end - %we need a profile of points on which to compute the mass_flux, is it here? - if isnans(md.qmu_mass_flux_profile), - error('process_qmu_response_data error message: could not find a mass_flux exp profile!'); - end - if ~(ischar(md.qmu_mass_flux_profile) && (indx == 1)) && ... - ~(iscell(md.qmu_mass_flux_profile) && length(md.qmu_mass_flux_profile) >= indx && ischar(md.qmu_mass_flux_profile{indx})), - error('process_qmu_response_data error message: mass_flux exp profile should be a domain outline name'); - end - - %ok, process the qmu_mass_flux_profile to build a list of segments: - if ischar(md.qmu_mass_flux_profile), - md.qmu_mass_flux_segments=MassFluxProcessProfile(md); - md.qmu_mass_flux_num_segments=size(md.qmu_mass_flux_segments,1); - else - if isempty(md.qmu_mass_flux_segments) || ~iscell(md.qmu_mass_flux_segments) - md.qmu_mass_flux_segments=cell(numel(md.qmu_mass_flux_profile),1); - end - if isempty(md.qmu_mass_flux_num_segments) || ... - (numel(md.qmu_mass_flux_num_segments) == 1 && isnan(md.qmu_mass_flux_num_segments)) - md.qmu_mass_flux_num_segments=zeros(numel(md.qmu_mass_flux_profile),1); - end - save=md.qmu_mass_flux_profile; - md.qmu_mass_flux_profile=md.qmu_mass_flux_profile{indx}; - md.qmu_mass_flux_segments{indx,1}=MassFluxProcessProfile(md); - md.qmu_mass_flux_num_segments(indx,1)=size(md.qmu_mass_flux_segments{indx},1); - md.qmu_mass_flux_profile=save; - clear save - end + %Do we have to process mass flux profiles? + if strncmpi(md.responsedescriptors{i},'MassFlux',8), + process_mass_flux_profiles=1; end end -% for now, concatenate all segments from cells to make a double array for parallel -if iscell(md.qmu_mass_flux_segments) - segments=zeros(0,5); - for i=1:length(md.qmu_mass_flux_num_segments) - if md.qmu_mass_flux_num_segments(i) - segments=[segments; md.qmu_mass_flux_segments{i}]; - end + +%deal with mass flux profiles +if process_mass_flux_profiles, + + %we need a profile of points on which to compute the mass_flux, is it here? + if isnans(md.qmu_mass_flux_profiles), + error('process_qmu_response_data error message: could not find a mass_flux exp profile!'); end - md.qmu_mass_flux_segments=segments; - clear segments + + if ~iscell(md.qmu_mass_flux_profiles), + error('process_qmu_response_data error message: qmu_mass_flux_profiles field should be a cell array of domain outline names'); + end + + if isempty(md.qmu_mass_flux_profiles), + error('process_qmu_response_data error message: qmu_mass_flux_profiles cannot be empty!'); + end + + %ok, process the domains named in qmu_mass_flux_profiles, to build a list of segments + md.qmu_mass_flux_segments=cell(numel(md.qmu_mass_flux_profiles),1); + + for i=1:numel(md.qmu_mass_flux_profiles), + md.qmu_mass_flux_segments{i}=MassFluxProcessProfile(md,md.qmu_mass_flux_profile_directory,md.qmu_mass_flux_profiles{i}); + end end +error; + Index: /issm/trunk/src/m/utils/Mesh/ProfileProjectOntoMesh.m =================================================================== --- /issm/trunk/src/m/utils/Mesh/ProfileProjectOntoMesh.m (revision 4772) +++ /issm/trunk/src/m/utils/Mesh/ProfileProjectOntoMesh.m (revision 4773) @@ -14,4 +14,5 @@ %[x0,y0] = intersections(profile.x,profile.y,x,y,1); [x0,y0,indices,j] = intersections(profile.x,profile.y,x,y); + % sort intersections to create segments in order and continuous along profile [indices,isort]=sort(indices); Index: /issm/trunk/src/mex/Makefile.am =================================================================== --- /issm/trunk/src/mex/Makefile.am (revision 4772) +++ /issm/trunk/src/mex/Makefile.am (revision 4773) @@ -44,4 +44,5 @@ Mergesolutionfromftog\ MeshPartition\ + MeshProfileIntersection\ Misfit \ ModelProcessor \ @@ -228,4 +229,7 @@ MeshPartition/MeshPartition.h +MeshProfileIntersection_SOURCES = MeshProfileIntersection/MeshProfileIntersection.cpp\ + MeshProfileIntersection/MeshProfileIntersection.h + ModelProcessor_SOURCES = ModelProcessor/ModelProcessor.cpp \ ModelProcessor/ModelProcessor.h Index: /issm/trunk/src/mex/MeshProfileIntersection/MeshProfileIntersection.cpp =================================================================== --- /issm/trunk/src/mex/MeshProfileIntersection/MeshProfileIntersection.cpp (revision 4773) +++ /issm/trunk/src/mex/MeshProfileIntersection/MeshProfileIntersection.cpp (revision 4773) @@ -0,0 +1,117 @@ +/*! \file MeshProfileIntersection.cpp + \brief: takes a .exp file (made of several profiles), and figures out its intersection + with a mesh. + + usage: + [segments]=MeshProfileIntersection(index,x,y,filename); + + where: + input: + index,x,y is a triangulation + filename: name of Argus style .exp file containing the segments (can be groups of disconnected segments) + output: + segments: array made of x1,y1,x2,y2,element_id lines (x1,y1) and (x2,y2) are segment extremitis for a segment + belonging to the elemnt_id element. there are as many lines in segments as there are segments intersecting the + mesh. +*/ + +#include "./MeshProfileIntersection.h" + +void mexFunction( int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) { + + int i,j; + + /* required input: */ + //mesh + double* double_index=NULL; + int* index=NULL; + int nel; + double* x=NULL; + double* y=NULL; + int nods; + int dummy; + + //contours + mxArray* matlabstructure=NULL; + Contour** contours=NULL; + int numcontours; + + + /*intermediate: */ + Contour* contouri=NULL; + + /* output: */ + double* segments=NULL; + int numsegs; + + /*Boot module: */ + MODULEBOOT(); + + /*checks on arguments on the matlab side: */ + CheckNumMatlabArguments(nlhs,NLHS,nrhs,NRHS,__FUNCT__,&MeshProfileIntersectionUsage); + + /*First, call expread on filename to build a contour array in the matlab workspace: */ + mexCallMATLAB( 1, &matlabstructure, 1, (mxArray**)&FILENAME, "expread"); + + /*Fetch inputs: */ + //index + FetchData(&double_index,&nel,&dummy,INDEX); + if(dummy!=3)ISSMERROR(" element triangulation should be of 3 column width!"); + index=(int*)xmalloc(nel*3*sizeof(int)); + for(i=0;inods=(int)mxGetScalar(mxGetField(matlabstructure,i,"nods")); + //set pointers. + contouri->x=mxGetPr(mxGetField(matlabstructure,i,"x")); + contouri->y=mxGetPr(mxGetField(matlabstructure,i,"y")); + *(contours+i)=contouri; + } + + /* Debugging of contours :{{{1*/ + /*for(i=0;inods); + for (j=0;jnods;j++){ + printf(" %lf %lf\n",*(contouri->x+j),*(contouri->y+j)); + } + }*/ + /*}}}*/ + + /*Run interpolation routine: */ + MeshProfileIntersectionx(&segments,&numsegs,index,x,y,nel,nods,contours,numcontours); + + /* output: */ + WriteData(SEGMENTS,segments,numsegs,5); + + /*end module: */ + MODULEEND(); + +} + +void MeshProfileIntersectionUsage(void){ + printf(" usage:\n"); + printf(" [segments]=MeshProfileIntersection(index,x,y,filename);\n"); + printf(" where:\n"); + printf(" input:\n"); + printf(" index,x,y is a triangulation\n"); + printf(" filename: name of Argus style .exp file containing the segments (can be groups of disconnected segments)\n"); + printf(" output:\n"); + printf(" segments: array made of x1,y1,x2,y2,element_id lines (x1,y1) and (x2,y2) are segment extremitis for a segment \n"); + printf(" belonging to the elemnt_id element. there are as many lines in segments as there are segments intersecting the \n"); + printf(" mesh.\n"); +} Index: /issm/trunk/src/mex/MeshProfileIntersection/MeshProfileIntersection.h =================================================================== --- /issm/trunk/src/mex/MeshProfileIntersection/MeshProfileIntersection.h (revision 4773) +++ /issm/trunk/src/mex/MeshProfileIntersection/MeshProfileIntersection.h (revision 4773) @@ -0,0 +1,36 @@ + +/* + MeshProfileIntersection.h +*/ + + +#ifndef _MESHPROFILEINTERSECTION_H +#define _MESHPROFILEINTERSECTION_H + +/* local prototypes: */ +void MeshProfileIntersectionUsage(void); + +#include "../../c/modules/modules.h" +#include "../../c/Container/Container.h" +#include "../../c/shared/shared.h" + +#undef __FUNCT__ +#define __FUNCT__ "MeshProfileIntersection" + +/* input macros: */ +#define INDEX prhs[0] +#define X prhs[1] +#define Y prhs[2] +#define FILENAME prhs[3] + +/* serial output macros: */ +#define SEGMENTS (mxArray**)&plhs[0] + +/* serial arg counts: */ +#undef NLHS +#define NLHS 1 +#undef NRHS +#define NRHS 4 + +#endif /* _MESHPROFILEINTERSECTION_H */ + Index: /issm/trunk/src/mex/ModelProcessor/ModelProcessor.cpp =================================================================== --- /issm/trunk/src/mex/ModelProcessor/ModelProcessor.cpp (revision 4772) +++ /issm/trunk/src/mex/ModelProcessor/ModelProcessor.cpp (revision 4773) @@ -35,5 +35,5 @@ FetchData(&double_analyses,&numanalyses,ANALYSES); analyses=(int*)xmalloc(numanalyses*sizeof(int)); - for(i=0;i