Index: /issm/trunk/src/c/io/FetchData.cpp =================================================================== --- /issm/trunk/src/c/io/FetchData.cpp (revision 5697) +++ /issm/trunk/src/c/io/FetchData.cpp (revision 5698) @@ -14,7 +14,4 @@ #ifdef _SERIAL_ -/***************** **************** **************** **************** **************** **************** **************** **************** **************** **************** - Serial Fetch Data Routines, all overloaded. -**************** **************** **************** **************** **************** **************** **************** **************** **************** *****************/ #include /*FUNCTION FetchData(DataSet** pdataset,const mxArray* dataref){{{1*/ @@ -298,9 +295,10 @@ void FetchData(bool* pboolean,const mxArray* dataref){ - bool boolean; + bool* mxbool_ptr=NULL; if (mxIsClass(dataref,"logical")){ - /*Recover the double: */ - boolean=(bool)mxGetScalar(dataref); + if(mxGetM(dataref)!=1) ISSMERROR("input data is not of size 1x1"); + if(mxGetN(dataref)!=1) ISSMERROR("input data is not of size 1x1"); + mxbool_ptr=mxGetLogicals(dataref); } else{ @@ -308,6 +306,5 @@ } - /*Assign output pointers:*/ - *pboolean=boolean; + *pboolean=*mxbool_ptr; } /*}}}*/ @@ -315,7 +312,4 @@ #if defined(_PARALLEL_) && defined(_HAVE_PETSC_) -/***************** **************** **************** **************** **************** **************** **************** **************** **************** **************** - Parallel Fetch Data Routines, all overloaded. -**************** **************** **************** **************** **************** **************** **************** **************** **************** *****************/ /*FUNCTION FetchData(double** pmatrix, int* pM,int* pN,FILE* fid){{{1*/ void FetchData(double** pmatrix, int* pM,int* pN,FILE* fid){ Index: /issm/trunk/src/c/modules/Solverx/Solverx.cpp =================================================================== --- /issm/trunk/src/c/modules/Solverx/Solverx.cpp (revision 5697) +++ /issm/trunk/src/c/modules/Solverx/Solverx.cpp (revision 5698) @@ -98,5 +98,15 @@ } } - KSPSolve(ksp,pf,uf); + if(verbose){ + double max,min; + KSPSetComputeSingularValues(ksp,PETSC_TRUE); + KSPSetUp(ksp); + KSPSolve(ksp,pf,uf); + KSPComputeExtremeSingularValues(ksp,&max,&min); + printf("Condition number = %g %g\n",max,min); + } + else{ + KSPSolve(ksp,pf,uf); + } /*Check convergence*/ Index: /issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp =================================================================== --- /issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp (revision 5697) +++ /issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp (revision 5698) @@ -55,5 +55,5 @@ /*===================== DEBUGING ====================*/ - if (count==1) SystemMatricesx(&Kgg_nopenalty,&pg_nopenalty,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); + if (count==1) SystemMatricesx(&Kgg_nopenalty,&pg_nopenalty,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,true,true,false,false); MatDuplicate(Kgg_nopenalty,MAT_COPY_VALUES,&Kgg); VecDuplicatePatch(&pg,pg_nopenalty); Index: /issm/trunk/src/m/classes/public/solversettomumps.m =================================================================== --- /issm/trunk/src/m/classes/public/solversettomumps.m (revision 5697) +++ /issm/trunk/src/m/classes/public/solversettomumps.m (revision 5698) @@ -5,5 +5,5 @@ % md=solversettomumps(md) -md.solverstring=' -mat_type aijmumps -ksp_type preonly -pc_type lu -mat_mumps_icntl_14 120 -pc_factor_shift_positive_definite true'; +md.solverstring='-mat_type aijmumps -ksp_type preonly -pc_type lu -mat_mumps_icntl_14 120 -pc_factor_shift_positive_definite true'; %optional Index: /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m =================================================================== --- /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m (revision 5697) +++ /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m (revision 5698) @@ -21,5 +21,5 @@ % ================= DEBUG FOR NOW ==================== if count==1 - [K_gg_nopenalty,p_g_nopenalty,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); + [K_gg_nopenalty,p_g_nopenalty,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,true,true,false,false); end [K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); Index: /issm/trunk/src/mex/Solver/Solver.cpp =================================================================== --- /issm/trunk/src/mex/Solver/Solver.cpp (revision 5697) +++ /issm/trunk/src/mex/Solver/Solver.cpp (revision 5698) @@ -8,9 +8,11 @@ /*input datasets: */ - Mat Kff=NULL; - Vec pf=NULL; - Vec uf0=NULL; - char* solver_string=NULL; - Parameters* parameters=NULL; + Mat Kff = NULL; + Vec pf = NULL; + Vec uf0 = NULL; + Vec uf = NULL; + char *solver_string = NULL; + Parameters *parameters = NULL; + int verbose; /*Matlab solver: */ @@ -18,8 +20,4 @@ char* matlabstring="-ksp_type matlab"; - /* output datasets: */ - Vec uf=NULL; - - /*Boot module: */ MODULEBOOT(); @@ -30,4 +28,5 @@ FetchParams(¶meters,PARAMETERS); parameters->FindParam(&solver_string,SolverStringEnum); + int verbose; parameters->FindParam(&verbose,VerboseEnum); /*Fetch rest of data only if not running the matlab solver: */ @@ -44,4 +43,5 @@ else{ /*Matlab solver: */ + if (verbose) printf(" Solving\n"); array[0]=KFF; array[1]=PF; @@ -57,5 +57,4 @@ xfree((void**)&solver_string); - /*end module: */ MODULEEND(); } Index: /issm/trunk/src/mex/SystemMatrices/SystemMatrices.cpp =================================================================== --- /issm/trunk/src/mex/SystemMatrices/SystemMatrices.cpp (revision 5697) +++ /issm/trunk/src/mex/SystemMatrices/SystemMatrices.cpp (revision 5698) @@ -14,4 +14,5 @@ Materials *materials = NULL; Parameters *parameters = NULL; + bool kflag,pflag,penalty_kflag,penalty_pflag; /* output datasets: */ @@ -24,5 +25,8 @@ /*checks on arguments on the matlab side: */ - CheckNumMatlabArguments(nlhs,NLHS,nrhs,NRHS,__FUNCT__,&SystemMatricesUsage); + if((nlhs!=NLHS) || (nrhs!=6 && nrhs!=10)){ + SystemMatricesUsage(); + ISSMERROR(" usage. See above"); + } /*Input datasets: */ @@ -41,5 +45,13 @@ /*!Generate internal degree of freedom numbers: */ - SystemMatricesx(&Kgg,&pg,&kmax,elements,nodes,vertices,loads,materials,parameters); + if(nrhs==10){ + FetchData(&kflag,KFLAG); + FetchData(&pflag,PFLAG); + FetchData(&penalty_kflag,PENALTYKFLAG); + FetchData(&penalty_pflag,PENALTYPFLAG); + SystemMatricesx(&Kgg,&pg,&kmax,elements,nodes,vertices,loads,materials,parameters,kflag,pflag,penalty_kflag,penalty_pflag); + } + else + SystemMatricesx(&Kgg,&pg,&kmax,elements,nodes,vertices,loads,materials,parameters); /*write output datasets: */ @@ -66,4 +78,5 @@ _printf_("\n"); _printf_(" usage: [Kgg,pg,kmax] = %s(elements,nodes,vertices,loads,materials,parameters);\n",__FUNCT__); + _printf_(" usage: [Kgg,pg,kmax] = %s(elements,nodes,vertices,loads,materials,parameters,kflag,pflag,penalty_kflag,penalty_pflag);\n",__FUNCT__); _printf_("\n"); } Index: /issm/trunk/src/mex/SystemMatrices/SystemMatrices.h =================================================================== --- /issm/trunk/src/mex/SystemMatrices/SystemMatrices.h (revision 5697) +++ /issm/trunk/src/mex/SystemMatrices/SystemMatrices.h (revision 5698) @@ -18,10 +18,14 @@ /* serial input macros: */ -#define ELEMENTS (mxArray *)prhs[0] -#define NODES (mxArray *)prhs[1] -#define VERTICES (mxArray *)prhs[2] -#define LOADS (mxArray *)prhs[3] -#define MATERIALS (mxArray *)prhs[4] -#define PARAMETERS (mxArray *)prhs[5] +#define ELEMENTS (mxArray *)prhs[0] +#define NODES (mxArray *)prhs[1] +#define VERTICES (mxArray *)prhs[2] +#define LOADS (mxArray *)prhs[3] +#define MATERIALS (mxArray *)prhs[4] +#define PARAMETERS (mxArray *)prhs[5] +#define KFLAG (mxArray *)prhs[6] +#define PFLAG (mxArray *)prhs[7] +#define PENALTYKFLAG (mxArray *)prhs[8] +#define PENALTYPFLAG (mxArray *)prhs[9] /* serial output macros: */ @@ -34,5 +38,5 @@ #define NLHS 3 #undef NRHS -#define NRHS 6 +#define NRHS 10 #endif /* _SYSTEMMATRICES_H */