Index: /issm/trunk-jpl/src/c/analyses/DakotaSpawnCore.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/DakotaSpawnCore.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/analyses/DakotaSpawnCore.cpp (revision 16149) @@ -41,5 +41,4 @@ char **responses_descriptors = NULL; //these are our! there are only numresponsedescriptors of them, not d_numresponses!!! int numresponsedescriptors; - char *string = NULL; int solution_type; bool control_analysis = false; Index: /issm/trunk-jpl/src/c/analyses/ResetBoundaryConditions.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/ResetBoundaryConditions.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/analyses/ResetBoundaryConditions.cpp (revision 16149) @@ -18,7 +18,4 @@ femmodel->SetCurrentConfiguration(analysis_type); - /*recover nodes: */ - nodes=femmodel->nodes; - /*retrieve boundary conditions from element inputs :*/ GetSolutionFromInputsx(&yg,femmodel); Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 16149) @@ -8220,20 +8220,12 @@ /*Intermediaries */ - int i,j,approximation; + int i,j; IssmDouble Jdet,viscosity,FSreconditioning,diameter,rigidity; IssmDouble xyz_list[NUMVERTICES][3]; IssmDouble epsilon[6]; /* epsilon=[exx,eyy,ezz,exy,exz,eyz];*/ - IssmDouble B[8][24]; - IssmDouble B_prime[8][24]; - IssmDouble B_stab[3][NUMVERTICES]; - IssmDouble D_scalar,D_scalar_stab; - IssmDouble D[8][8]={0.0}; - IssmDouble D_stab[3][3]={0.0}; - IssmDouble Ke_temp[24][24]={0.0}; //for the six nodes - IssmDouble Ke_temp_stab[6][6]={0.0}; //for the six nodes GaussPenta *gauss=NULL; /*Stabilization*/ - bool stabilization = true; + IssmDouble D_scalar; IssmDouble dbasis[3][6]; IssmDouble dmu[3]; @@ -8331,5 +8323,5 @@ /*Intermediaries */ - int i,j,approximation; + int i,approximation; IssmDouble Jdet,viscosity,FSreconditioning,D_scalar; IssmDouble xyz_list[NUMVERTICES][3]; @@ -10940,5 +10932,4 @@ int* pdoflist=NULL; IssmDouble FSreconditioning; - GaussPenta *gauss; /*Fetch number of nodes and dof for this finite element*/ @@ -11276,5 +11267,4 @@ IssmDouble h[NUMVERTICES],s[NUMVERTICES],b[NUMVERTICES],r[NUMVERTICES]; IssmDouble melting[NUMVERTICES],phi[NUMVERTICES]; - bool grounded[NUMVERTICES],floating[NUMVERTICES]; if(!IsOnBed()) return; @@ -11301,6 +11291,4 @@ b[i] = r[i]; s[i] = b[i]+h[i]; - floating[i] = false; - grounded[i] = true; } } @@ -11314,12 +11302,8 @@ s[i] = (1-density)*h[i]; b[i] = -density*h[i]; - floating[i] = true; - grounded[i] = false; } else if(migration_style==SoftMigrationEnum && phi_ungrounding[vertices[i]->Pid()]<0.){ s[i] = (1-density)*h[i]; b[i] = -density*h[i]; - floating[i] = true; - grounded[i] = false; } else{ Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 16149) @@ -1777,5 +1777,6 @@ switch(type){ case VertexPIdEnum: - for (int i=0;i(NUMVERTICES); + for(int i=0;ivertices[i]->Pid()]; } @@ -1787,12 +1788,12 @@ this->inputs->AddInput(new TriaInput(name,values,P1Enum)); } - return; case VertexSIdEnum: - for (int i=0;i(NUMVERTICES); + for(int i=0;ivertices[i]->Sid()]; } /*update input*/ - if (name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsDamageDEnum || name==MaterialsDamageDbarEnum){ + if(name==MaterialsRheologyBbarEnum || name==MaterialsRheologyBEnum || name==MaterialsDamageDEnum || name==MaterialsDamageDbarEnum){ material->inputs->AddInput(new TriaInput(name,values,P1Enum)); } @@ -1800,5 +1801,4 @@ this->inputs->AddInput(new TriaInput(name,values,P1Enum)); } - return; case NodesEnum: @@ -1809,16 +1809,9 @@ GetDofList(&doflist,NoneApproximationEnum,GsetEnum); - /*Use the dof list to index into the vector: */ for(int i=0;i(values[i])) _error_("NaN found in vector"); } - - /*Add input to the element: */ this->inputs->AddInput(new TriaInput(name,values,this->element_type)); - - /*Free ressources:*/ - xDelete(doflist); - return; case NodeSIdEnum: @@ -1832,7 +1825,5 @@ if(xIsNan(values[i])) _error_("NaN found in vector"); } - /*Add input to the element: */ this->inputs->AddInput(new TriaInput(name,values,this->element_type)); - return; default: @@ -7618,5 +7609,4 @@ IssmDouble melting[NUMVERTICES],phi[NUMVERTICES];; IssmDouble h[NUMVERTICES],s[NUMVERTICES],b[NUMVERTICES],r[NUMVERTICES]; - bool grounded[NUMVERTICES],floating[NUMVERTICES]; /*Recover info at the vertices: */ @@ -7641,6 +7631,4 @@ b[i] = r[i]; s[i] = b[i]+h[i]; - floating[i] = false; - grounded[i] = true; } } @@ -7654,12 +7642,8 @@ s[i] = (1-density)*h[i]; b[i] = -density*h[i]; - floating[i] = true; - grounded[i] = false; } else if(migration_style==SoftMigrationEnum && phi_ungrounding[vertices[i]->Pid()]<0.){ s[i] = (1-density)*h[i]; b[i] = -density*h[i]; - floating[i] = true; - grounded[i] = false; } else{ Index: /issm/trunk-jpl/src/c/classes/ExternalResults/GenericExternalResult.h =================================================================== --- /issm/trunk-jpl/src/c/classes/ExternalResults/GenericExternalResult.h (revision 16148) +++ /issm/trunk-jpl/src/c/classes/ExternalResults/GenericExternalResult.h (revision 16149) @@ -22,10 +22,10 @@ private: - int id; - int enum_type; + int id; + int enum_type; ResultType value; - int M; - int N; - int step; + int M; + int N; + int step; IssmDouble time; @@ -60,21 +60,27 @@ /*GenericExternalResult constructors and destructors*/ GenericExternalResult(){ /*{{{*/ - id = 0; + id = 0; enum_type = NoneEnum; - M=0; - N=0; - step=0; - time=0; + M = 0; + N = 0; + step = 0; + time = 0; } /*}}}*/ GenericExternalResult(int in_id, int in_enum_type,ResultType in_values, int in_M,int in_N,int in_step,IssmDouble in_time){/*{{{*/ + id = 0; + enum_type = NoneEnum; + M = 0; + N = 0; + step = 0; + time = 0; _error_("template GenericExternalResult(int in_id, int in_enum_type,double* in_values, int in_M,int in_N,int in_step,IssmDouble in_time) not implemented for this ResultType\n"); } /*}}}*/ GenericExternalResult(int in_id, int in_enum_type,ResultType in_value,int in_step, IssmDouble in_time){ /*{{{*/ - id=in_id; - enum_type=in_enum_type; - value=in_value; - step=in_step; - time=in_time; + id = in_id; + enum_type = in_enum_type; + value = in_value; + step = in_step; + time = in_time; } /*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp (revision 16149) @@ -396,13 +396,13 @@ /*}}}*/ /*FUNCTION Inputs::AXPY{{{*/ -void Inputs::AXPY(int MeshYEnum, IssmDouble scalar, int MeshXEnum){ +void Inputs::AXPY(int inputy_enum, IssmDouble scalar, int inputx_enum){ /*Find x and y inputs: */ - Input* xinput=dynamic_cast(this->GetInput(MeshXEnum)); - Input* yinput=dynamic_cast(this->GetInput(MeshYEnum)); + Input* xinput=dynamic_cast(this->GetInput(inputx_enum)); + Input* yinput=dynamic_cast(this->GetInput(inputy_enum)); /*some checks: */ - if(!xinput) _error_("input " << EnumToStringx(MeshXEnum) << " could not be found!"); - if(!yinput) _error_("input " << EnumToStringx(MeshYEnum) << " could not be found!"); + if(!xinput) _error_("input " << EnumToStringx(inputx_enum) << " could not be found!"); + if(!yinput) _error_("input " << EnumToStringx(inputy_enum) << " could not be found!"); /*Apply AXPY: */ Index: /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h (revision 16148) +++ /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h (revision 16149) @@ -22,12 +22,12 @@ /*numerics*/ - int AddInput(Input* in_input); - void ChangeEnum(int enumtype,int new_enumtype); - void ConstrainMin(int constrain_enum, IssmDouble minimum); - int DeleteInput(int enum_type); - void DuplicateInput(int original_enum,int new_enum); - Input* GetInput(int enum_name); - Inputs* SpawnTriaInputs(int position); - void AXPY(int MeshYEnum, IssmDouble scalar, int MeshXEnum); + int AddInput(Input* in_input); + void ChangeEnum(int enumtype,int new_enumtype); + void ConstrainMin(int constrain_enum, IssmDouble minimum); + int DeleteInput(int enum_type); + void DuplicateInput(int original_enum,int new_enum); + Input* GetInput(int enum_name); + Inputs* SpawnTriaInputs(int position); + void AXPY(int inputy_enum, IssmDouble scalar, int inputx_enum); IssmDouble InfinityNorm(int enumtype); IssmDouble Max(int enumtype); @@ -35,9 +35,9 @@ IssmDouble Min(int enumtype); IssmDouble MinAbs(int enumtype); - void GetInputAverage(IssmDouble* pvalue, int enum_type); - void GetInputValue(bool* pvalue,int enum_type); - void GetInputValue(int* pvalue,int enum_type); - void GetInputValue(IssmDouble* pvalue,int enum_type); - void Configure(Parameters* parameters); + void GetInputAverage(IssmDouble* pvalue, int enum_type); + void GetInputValue(bool* pvalue,int enum_type); + void GetInputValue(int* pvalue,int enum_type); + void GetInputValue(IssmDouble* pvalue,int enum_type); + void Configure(Parameters* parameters); }; Index: /issm/trunk-jpl/src/c/classes/Params/DoubleMatParam.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Params/DoubleMatParam.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/classes/Params/DoubleMatParam.cpp (revision 16149) @@ -58,5 +58,5 @@ _printf_(" matrix size: " << this->M << "x" << this->N << "\n"); for(i=0;iM;i++){ - for(i=0;iN;i++){ + for(j=0;jN;j++){ _printf_(i << " " << j << " " << *(this->value+N*i+j) << "\n"); } Index: /issm/trunk-jpl/src/c/classes/Params/IntMatParam.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Params/IntMatParam.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/classes/Params/IntMatParam.cpp (revision 16149) @@ -58,5 +58,5 @@ _printf_(" matrix size: " << this->M << "x" << this->N << "\n"); for(i=0;iM;i++){ - for(i=0;iN;i++){ + for(j=0;jN;j++){ _printf_("(" << i << "," << j << ") " << *(this->value+N*i+j) << "\n"); } Index: /issm/trunk-jpl/src/c/modules/Shp2Expx/Shp2Expx.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/Shp2Expx/Shp2Expx.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/modules/Shp2Expx/Shp2Expx.cpp (revision 16149) @@ -12,6 +12,5 @@ #ifdef _HAVE_SHAPELIB_ //only works if Shapelib library has been compiled in. - return(Shp2Expx(filshp,filexp, - 0)); + return(Shp2Expx(filshp,filexp,0)); #else //ifdef _HAVE_SHAPELIB_ @@ -20,16 +19,12 @@ } -int Shp2Expx(char* filshp,char* filexp, - int sgn){ +int Shp2Expx(char* filshp,char* filexp,int sgn){ #ifdef _HAVE_SHAPELIB_ //only works if Shapelib library has been compiled in. double cm,sp; - - if (sgn) - Xy2lldef(&cm,&sp,sgn); - - return(Shp2Expx(filshp,filexp, - sgn,cm,sp)); + Xy2lldef(&cm,&sp,sgn); + + return(Shp2Expx(filshp,filexp,sgn,cm,sp)); #else //ifdef _HAVE_SHAPELIB_ @@ -38,6 +33,5 @@ } -int Shp2Expx(char* filshp,char* filexp, - int sgn,double cm,double sp){ +int Shp2Expx(char* filshp,char* filexp,int sgn,double cm,double sp){ #ifdef _HAVE_SHAPELIB_ //only works if Shapelib library has been compiled in. Index: /issm/trunk-jpl/src/c/modules/Shp2Kmlx/Shp2Kmlx.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/Shp2Kmlx/Shp2Kmlx.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/modules/Shp2Kmlx/Shp2Kmlx.cpp (revision 16149) @@ -8,16 +8,12 @@ #include "../../kml/kmlobjects.h" -int Shp2Kmlx(char* filshp,char* filkml, - int sgn){ +int Shp2Kmlx(char* filshp,char* filkml,int sgn){ #ifdef _HAVE_SHAPELIB_ //only works if Shapelib library has been compiled in. double cm,sp; - - if (sgn) - Xy2lldef(&cm,&sp,sgn); - - return(Shp2Kmlx(filshp,filkml, - sgn,cm,sp)); + Xy2lldef(&cm,&sp,sgn); + + return(Shp2Kmlx(filshp,filkml,sgn,cm,sp)); #else //ifdef _HAVE_SHAPELIB_ @@ -26,6 +22,5 @@ } -int Shp2Kmlx(char* filshp,char* filkml, - int sgn,double cm,double sp){ +int Shp2Kmlx(char* filshp,char* filkml,int sgn,double cm,double sp){ #ifdef _HAVE_SHAPELIB_ //only works if Shapelib library has been compiled in. @@ -103,9 +98,6 @@ pshapm=xNew(nshape); -/* loop over the list of shapes */ - - for( i = 0; i < nEntities; i++ ) - { -// int j; + /* loop over the list of shapes */ + for(i=0;i @@ -21,8 +21,7 @@ /*Mumps header files: */ -#include "dmumps_c.h" -/*}}}*/ - -void MumpsInit(DMUMPS_STRUC_C &theMumpsStruc) { +#include + +void MumpsInit(DMUMPS_STRUC_C &theMumpsStruc){ theMumpsStruc.par = 1; theMumpsStruc.sym = 0; @@ -71,18 +70,11 @@ } -void MumpsSolve(int n, - int nnz, - int local_nnz, - int* irn_loc, - int* jcn_loc, - IssmPDouble *a_loc, - IssmPDouble *rhs, - Parameters* parameters=0 /*unused here*/) { - /*Initialize mumps: {{{*/ +void MumpsSolve(int n,int nnz,int local_nnz,int* irn_loc,int* jcn_loc, IssmPDouble *a_loc, IssmPDouble *rhs, Parameters* parameters=0 /*unused here*/){ + /*Initialize mumps*/ DMUMPS_STRUC_C theMumpsStruc; MumpsInit(theMumpsStruc); MumpsSettings(theMumpsStruc); - /*}}}*/ - // now setup the rest of theMumpsStruc + + /*now setup the rest of theMumpsStruc */ theMumpsStruc.n=n; theMumpsStruc.nz=nnz; @@ -94,9 +86,9 @@ theMumpsStruc.nrhs=1; theMumpsStruc.lrhs=1; - /*Solve system: {{{*/ + + /*Solve system*/ MumpsAnalyze(theMumpsStruc); MumpsFactorize(theMumpsStruc); MumpsBacksubstitute(theMumpsStruc); - /*}}}*/ MumpsFinalize(theMumpsStruc); } @@ -116,34 +108,31 @@ void MpiDenseMumpsSolve( /*output: */ IssmDouble* uf, int uf_M, int uf_m, /*matrix input: */ IssmDouble* Kff, int Kff_M, int Kff_N, int Kff_m, /*right hand side vector: */ IssmDouble* pf, int pf_M, int pf_m, Parameters* parameters){ /*{{{*/ - /*Variables: {{{*/ - - ISSM_MPI_Comm comm; - int my_rank; - int num_procs; - int i; - int j; - int nnz ,local_nnz; - int *irn_loc = NULL; - int *jcn_loc = NULL; - IssmDouble *a_loc = NULL; - int count; - int lower_row; - int upper_row; - IssmDouble* rhs=NULL; - int* recvcounts=NULL; - int* displs=NULL; - /*}}}*/ - /*Communicator info:{{{ */ - my_rank=IssmComm::GetRank(); - num_procs=IssmComm::GetSize(); - comm=IssmComm::GetComm(); - /*}}}*/ - /*First, some checks:{{{ */ + /*Variables*/ + ISSM_MPI_Comm comm; + int my_rank; + int num_procs; + int i,j; + int nnz,local_nnz; + int *irn_loc = NULL; + int *jcn_loc = NULL; + IssmDouble *a_loc = NULL; + int count; + int lower_row; + int upper_row; + IssmDouble *rhs = NULL; + int *recvcounts = NULL; + int *displs = NULL; + + /*Communicator info */ + my_rank = IssmComm::GetRank(); + num_procs = IssmComm::GetSize(); + comm = IssmComm::GetComm(); + + /*First, some checks*/ if (Kff_M!=Kff_N)_error_("stiffness matrix Kff should be square"); if (uf_M!=Kff_M | uf_M!=pf_M)_error_("solution vector should be the same size as stiffness matrix Kff and load vector pf"); if (uf_m!=Kff_m | uf_m!=pf_m)_error_("solution vector should be locally the same size as stiffness matrix Kff and load vector pf"); - /*}}}*/ - /*Initialize matrix:{{{ */ - + + /*Initialize matrix */ /*figure out number of non-zero entries: */ local_nnz=0; @@ -192,18 +181,10 @@ ISSM_MPI_Gatherv(pf, pf_m, ISSM_MPI_DOUBLE, rhs, recvcounts, displs, ISSM_MPI_DOUBLE,0,comm); - MumpsSolve(Kff_M, - nnz, - local_nnz, - irn_loc, - jcn_loc, - a_loc, - rhs, - parameters); - - /*Now scatter from cpu 0 to all other cpus: {{{*/ + MumpsSolve(Kff_M,nnz,local_nnz,irn_loc,jcn_loc,a_loc,rhs,parameters); + + /*Now scatter from cpu 0 to all other cpus*/ ISSM_MPI_Scatterv( rhs, recvcounts, displs, ISSM_MPI_DOUBLE, uf, uf_m, ISSM_MPI_DOUBLE, 0, comm); - /*}}}*/ - /*Cleanup: {{{*/ + /*Cleanup*/ xDelete(irn_loc); xDelete(jcn_loc); @@ -212,5 +193,4 @@ xDelete(recvcounts); xDelete(displs); - /*}}}*/ } /*}}}*/ @@ -245,12 +225,5 @@ } -void MumpsSolve(int n, - int nnz, - int local_nnz, - int* irn_loc, - int* jcn_loc, - IssmDouble *a_loc, - IssmDouble *rhs, - Parameters* parameters) { +void MumpsSolve(int n,int nnz,int local_nnz,int* irn_loc,int* jcn_loc,IssmDouble *a_loc,IssmDouble *rhs,Parameters* parameters){ int packedDimsSparseArrLength=1+1+1+local_nnz+local_nnz; int *packedDimsSparseArr=xNew(packedDimsSparseArrLength); Index: /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscSolver.cpp =================================================================== --- /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscSolver.cpp (revision 16148) +++ /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscSolver.cpp (revision 16149) @@ -170,5 +170,4 @@ IssmDouble* df_local=NULL; int df_local_size; - int i; int* pressure_indices=NULL; @@ -188,5 +187,5 @@ pressure_num=0; velocity_num=0; - for(i=0;i