Changeset 23044
- Timestamp:
- 08/02/18 16:17:04 (7 years ago)
- Location:
- issm/trunk-jpl/src/c/toolkits
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk-jpl/src/c/toolkits/issm/IssmDenseMat.h
r22737 r23044 272 272 pf=(IssmSeqVec<IssmDouble>*)pfin; 273 273 274 #ifdef _HAVE_MUMPS_ 275 /*Assume we have a sequential vec, downcast*/ 276 uf=((IssmSeqVec<IssmDouble>*)pfin)->Duplicate(); 277 SeqDenseMumpsSolve(uf->vector,uf->M,uf->M, /*stiffness matrix:*/ this->matrix,this->M,this->N,this->M, /*right hand side load vector: */ pf->vector,pf->M,pf->M,parameters); 278 return uf; 279 #endif 280 274 281 #ifdef _HAVE_GSL_ 275 282 DenseGslSolve(/*output*/ &x,/*stiffness matrix:*/ this->matrix,this->M,this->N, /*right hand side load vector: */ pf->vector,pf->M,parameters); … … 277 284 uf=new IssmSeqVec<IssmDouble>(x,this->N); xDelete(x); 278 285 return uf; 279 #else280 _error_("GSL support not compiled in!");281 286 #endif 287 288 _error_("No solver available"); 289 return NULL; 282 290 283 291 }/*}}}*/ -
issm/trunk-jpl/src/c/toolkits/issm/IssmMat.h
r16415 r23044 17 17 #include "../ToolkitOptions.h" 18 18 #include "./IssmToolkitUtils.h" 19 19 #include "../mumps/mumpsincludes.h" 20 20 /*}}}*/ 21 21 -
issm/trunk-jpl/src/c/toolkits/mumps/MumpsSolve.cpp
r21615 r23044 27 27 theMumpsStruc.par = 1; 28 28 theMumpsStruc.sym = 0; 29 #ifdef _HAVE_MPI_ 29 30 theMumpsStruc.comm_fortran = MPI_Comm_c2f(IssmComm::GetComm()); 31 #endif 30 32 theMumpsStruc.job = -1; 31 33 dmumps_c(&theMumpsStruc); … … 35 37 void MumpsSettings(DMUMPS_STRUC_C &theMumpsStruc) { 36 38 /*Control statements:{{{ */ 37 theMumpsStruc.icntl[1-1] = 6; //error verbose: default 6, 0 or negative -> suppressed39 theMumpsStruc.icntl[1-1] = 0; //error verbose: default 6, 0 or negative -> suppressed 38 40 theMumpsStruc.icntl[2-1] = 0; //std verbose: default 1, 0 or negative -> suppressed 39 41 theMumpsStruc.icntl[3-1] = 0; //global information verbose: default 6, 0 or negative -> suppressed … … 108 110 #endif 109 111 112 #ifdef _HAVE_MPI_ 110 113 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){ /*{{{*/ 111 114 … … 201 204 xDelete<int>(displs); 202 205 } /*}}}*/ 203 204 206 void MpiSparseMumpsSolve( /*output: */ IssmDouble* uf, int uf_M, int uf_m, /*matrix input: */ SparseRow<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){ /*{{{*/ 205 207 … … 287 289 xDelete<int>(displs); 288 290 } /*}}}*/ 291 #endif 292 void SeqDenseMumpsSolve(IssmDouble* uf,int uf_M,int uf_m, IssmDouble* Kff,int Kff_M, int Kff_N, int Kff_m, IssmDouble* pf, int pf_M, int pf_m, Parameters* parameters){/*{{{*/ 293 294 /*Variables*/ 295 int *irn_loc = NULL; 296 int *jcn_loc = NULL; 297 IssmDouble *a_loc = NULL; 298 IssmDouble *rhs = NULL; 299 300 /*First, some checks*/ 301 if (Kff_M!=Kff_N)_error_("stiffness matrix Kff should be square"); 302 if (Kff_M!=Kff_m)_error_("stiffness matrix Kff is not serial"); 303 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"); 304 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"); 305 306 /*Initialize matrix */ 307 /*figure out number of non-zero entries: */ 308 int nnz = 0; 309 for(int i=0;i<Kff_M;i++){ 310 for(int j=0;j<Kff_N;j++){ 311 if(Kff[i*Kff_N+j]!=0)nnz++; 312 } 313 } 314 315 /*Allocate: */ 316 if(nnz){ 317 irn_loc = xNew<int>(nnz); 318 jcn_loc = xNew<int>(nnz); 319 a_loc = xNew<IssmDouble>(nnz); 320 } 321 322 /*Populate the triplets: */ 323 int count=0; 324 for(int i=0;i<Kff_M;i++){ 325 for(int j=0;j<Kff_N;j++){ 326 if(Kff[i*Kff_N+j]!=0){ 327 irn_loc[count] = i+1; //fortran indexing 328 jcn_loc[count] = j+1; //fortran indexing 329 a_loc[count] = Kff[i*Kff_N+j]; 330 count++; 331 } 332 } 333 } 334 335 /*Deal with right hand side. We need to ISSM_MPI_Gather it onto cpu 0: */ 336 // AD performance is sensitive to calls to ensurecontiguous. 337 // Providing "t" will cause ensurecontiguous to be called. 338 #ifdef _HAVE_AD_ 339 rhs=xNew<IssmDouble>(pf_M,"t"); 340 #else 341 rhs=xNew<IssmDouble>(pf_M); 342 #endif 343 xMemCpy<IssmDouble>(rhs,pf,Kff_M); 344 345 MumpsSolve(Kff_M,nnz,nnz,irn_loc,jcn_loc,a_loc,rhs,parameters); 346 347 /*Now scatter from cpu 0 to all other cpus*/ 348 xMemCpy<IssmDouble>(uf,rhs,Kff_M); 349 350 /*Cleanup*/ 351 xDelete<int>(irn_loc); 352 xDelete<int>(jcn_loc); 353 xDelete<IssmDouble>(a_loc); 354 xDelete<IssmDouble>(rhs); 355 }/*}}}*/ 289 356 290 357 #ifdef _HAVE_ADOLC_ -
issm/trunk-jpl/src/c/toolkits/mumps/mumpsincludes.h
r16415 r23044 19 19 template <class doubletype> class SparseRow; 20 20 21 #ifdef _HAVE_MPI_ 21 22 void MpiDenseMumpsSolve(IssmDouble* uf,int uf_M,int uf_n, IssmDouble* Kff,int Kff_M, int Kff_N, int Kff_m, IssmDouble* pf, int pf_M, int pf_m, Parameters* parameters); 22 23 void MpiSparseMumpsSolve(IssmDouble* uf,int uf_M,int uf_n, SparseRow<IssmDouble>** Kff,int Kff_M, int Kff_N, int Kff_m, IssmDouble* pf, int pf_M, int pf_m, Parameters* parameters); 24 #endif 25 void SeqDenseMumpsSolve(IssmDouble* uf,int uf_M,int uf_n, IssmDouble* Kff,int Kff_M, int Kff_N, int Kff_m, IssmDouble* pf, int pf_M, int pf_m, Parameters* parameters); 23 26 24 27 #if defined(_HAVE_ADOLC_) && !defined(_WRAPPERS_)
Note:
See TracChangeset
for help on using the changeset viewer.