Ice Sheet System Model  4.18
Code documentation
Functions
MumpsSolve.cpp File Reference

: solve dense matrix system with MUMPS More...

#include "../../shared/Numerics/types.h"
#include "../../shared/MemOps/MemOps.h"
#include "../../shared/Exceptions/exceptions.h"
#include "../../shared/io/Comm/IssmComm.h"
#include "../../classes/Params/GenericParam.h"
#include "../../classes/Params/Parameters.h"
#include "../mpi/issmmpi.h"
#include "../adolc/adolcincludes.h"
#include "../codipack/codipackincludes.h"
#include "../issm/SparseRow.h"
#include "./mumpsincludes.h"
#include <dmumps_c.h>

Go to the source code of this file.

Functions

void MumpsInit (DMUMPS_STRUC_C &theMumpsStruc)
 
void MumpsSettings (DMUMPS_STRUC_C &theMumpsStruc)
 
void MumpsAnalyze (DMUMPS_STRUC_C &theMumpsStruc)
 
void MumpsFactorize (DMUMPS_STRUC_C &theMumpsStruc)
 
void MumpsBacksubstitute (DMUMPS_STRUC_C &theMumpsStruc)
 
void MumpsFinalize (DMUMPS_STRUC_C &theMumpsStruc)
 
void MumpsSolve (int n, int nnz, int local_nnz, int *irn_loc, int *jcn_loc, IssmPDouble *a_loc, IssmPDouble *rhs, Parameters *parameters=0)
 
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)
 

Detailed Description

: solve dense matrix system with MUMPS

Definition in file MumpsSolve.cpp.

Function Documentation

◆ MumpsInit()

void MumpsInit ( DMUMPS_STRUC_C &  theMumpsStruc)

Definition at line 27 of file MumpsSolve.cpp.

27  {
28  theMumpsStruc.n = 0;
29  theMumpsStruc.nz = 0;
30  theMumpsStruc.a = NULL;
31  theMumpsStruc.jcn = NULL;
32  theMumpsStruc.irn = NULL;
33  theMumpsStruc.par = 1;
34  theMumpsStruc.sym = 0;
35  #ifdef _HAVE_MPI_
36  theMumpsStruc.comm_fortran = MPI_Comm_c2f(IssmComm::GetComm());
37  #endif
38  theMumpsStruc.job = -1;
39  dmumps_c(&theMumpsStruc);
40 }

◆ MumpsSettings()

void MumpsSettings ( DMUMPS_STRUC_C &  theMumpsStruc)

Definition at line 43 of file MumpsSolve.cpp.

43  {
44  /*Control statements:{{{ */
45  theMumpsStruc.icntl[1-1] = 0; //error verbose: default 6, 0 or negative -> suppressed
46  theMumpsStruc.icntl[2-1] = 0; //std verbose: default 1, 0 or negative -> suppressed
47  theMumpsStruc.icntl[3-1] = 0; //global information verbose: default 6, 0 or negative -> suppressed
48  theMumpsStruc.icntl[4-1] = 0; //verbose everything: default is 4
49  theMumpsStruc.icntl[5-1] = 0;
50  theMumpsStruc.icntl[18-1] = 3;
51 
52  theMumpsStruc.icntl[20-1] = 0;
53  theMumpsStruc.icntl[21-1] = 0;
54  theMumpsStruc.icntl[30-1] = 0;
55  /*}}}*/
56 }

◆ MumpsAnalyze()

void MumpsAnalyze ( DMUMPS_STRUC_C &  theMumpsStruc)

Definition at line 59 of file MumpsSolve.cpp.

59  {
60  theMumpsStruc.job = 1;
61  dmumps_c(&theMumpsStruc);
62 }

◆ MumpsFactorize()

void MumpsFactorize ( DMUMPS_STRUC_C &  theMumpsStruc)

Definition at line 65 of file MumpsSolve.cpp.

65  {
66  theMumpsStruc.job = 2;
67  dmumps_c(&theMumpsStruc);
68 }

◆ MumpsBacksubstitute()

void MumpsBacksubstitute ( DMUMPS_STRUC_C &  theMumpsStruc)

Definition at line 71 of file MumpsSolve.cpp.

71  {
72  theMumpsStruc.job = 3;
73  dmumps_c(&theMumpsStruc);
74 }

◆ MumpsFinalize()

void MumpsFinalize ( DMUMPS_STRUC_C &  theMumpsStruc)

Definition at line 77 of file MumpsSolve.cpp.

77  {
78  theMumpsStruc.job = -2;
79  dmumps_c(&theMumpsStruc);
80 }

◆ MumpsSolve()

void MumpsSolve ( int  n,
int  nnz,
int  local_nnz,
int *  irn_loc,
int *  jcn_loc,
IssmPDouble *  a_loc,
IssmPDouble *  rhs,
Parameters *  parameters = 0 
)

Definition at line 82 of file MumpsSolve.cpp.

82  {
83  /*Initialize mumps*/
84  DMUMPS_STRUC_C theMumpsStruc;
85  MumpsInit(theMumpsStruc);
86  MumpsSettings(theMumpsStruc);
87 
88  /*now setup the rest of theMumpsStruc */
89  theMumpsStruc.n=n;
90  theMumpsStruc.nz=nnz;
91  theMumpsStruc.nz_loc=local_nnz;
92  theMumpsStruc.irn_loc=irn_loc;
93  theMumpsStruc.jcn_loc=jcn_loc;
94  theMumpsStruc.a_loc=a_loc;
95  theMumpsStruc.rhs=rhs;
96  theMumpsStruc.nrhs=1;
97  theMumpsStruc.lrhs=1;
98 
99  /*Solve system*/
100  MumpsAnalyze(theMumpsStruc);
101  MumpsFactorize(theMumpsStruc);
102  MumpsBacksubstitute(theMumpsStruc);
103  MumpsFinalize(theMumpsStruc);
104 }

◆ SeqDenseMumpsSolve()

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 
)

Definition at line 296 of file MumpsSolve.cpp.

296  {/*{{{*/
297 
298  /*Variables*/
299  int *irn_loc = NULL;
300  int *jcn_loc = NULL;
301  IssmDouble *a_loc = NULL;
302  IssmDouble *rhs = NULL;
303 
304  /*First, some checks*/
305  if (Kff_M!=Kff_N)_error_("stiffness matrix Kff should be square");
306  if (Kff_M!=Kff_m)_error_("stiffness matrix Kff is not serial");
307  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");
308  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");
309 
310  /*Initialize matrix */
311  /*figure out number of non-zero entries: */
312  int nnz = 0;
313  for(int i=0;i<Kff_M;i++){
314  for(int j=0;j<Kff_N;j++){
315  if(Kff[i*Kff_N+j]!=0)nnz++;
316  }
317  }
318 
319  /*Allocate: */
320  if(nnz){
321  irn_loc = xNew<int>(nnz);
322  jcn_loc = xNew<int>(nnz);
323  a_loc = xNew<IssmDouble>(nnz);
324  }
325 
326  /*Populate the triplets: */
327  int count=0;
328  for(int i=0;i<Kff_M;i++){
329  for(int j=0;j<Kff_N;j++){
330  if(Kff[i*Kff_N+j]!=0){
331  irn_loc[count] = i+1; //fortran indexing
332  jcn_loc[count] = j+1; //fortran indexing
333  a_loc[count] = Kff[i*Kff_N+j];
334  count++;
335  }
336  }
337  }
338 
339  /*Deal with right hand side. We need to ISSM_MPI_Gather it onto cpu 0: */
340 // AD performance is sensitive to calls to ensurecontiguous.
341 // Providing "t" will cause ensurecontiguous to be called.
342 #ifdef _HAVE_AD_
343  rhs=xNew<IssmDouble>(pf_M,"t");
344 #else
345  rhs=xNew<IssmDouble>(pf_M);
346 #endif
347  xMemCpy<IssmDouble>(rhs,pf,Kff_M);
348 
349  MumpsSolve(Kff_M,nnz,nnz,irn_loc,jcn_loc,a_loc,rhs,parameters);
350 
351  /*Now scatter from cpu 0 to all other cpus*/
352  xMemCpy<IssmDouble>(uf,rhs,Kff_M);
353 
354  /*Cleanup*/
355  xDelete<int>(irn_loc);
356  xDelete<int>(jcn_loc);
357  xDelete<IssmDouble>(a_loc);
358  xDelete<IssmDouble>(rhs);
359 }/*}}}*/
IssmDouble
double IssmDouble
Definition: types.h:37
MumpsSolve
void MumpsSolve(int n, int nnz, int local_nnz, int *irn_loc, int *jcn_loc, IssmPDouble *a_loc, IssmPDouble *rhs, Parameters *parameters=0)
Definition: MumpsSolve.cpp:82
IssmComm::GetComm
static ISSM_MPI_Comm GetComm(void)
Definition: IssmComm.cpp:30
MumpsBacksubstitute
void MumpsBacksubstitute(DMUMPS_STRUC_C &theMumpsStruc)
Definition: MumpsSolve.cpp:71
MumpsInit
void MumpsInit(DMUMPS_STRUC_C &theMumpsStruc)
Definition: MumpsSolve.cpp:27
MumpsAnalyze
void MumpsAnalyze(DMUMPS_STRUC_C &theMumpsStruc)
Definition: MumpsSolve.cpp:59
_error_
#define _error_(StreamArgs)
Definition: exceptions.h:49
MumpsSettings
void MumpsSettings(DMUMPS_STRUC_C &theMumpsStruc)
Definition: MumpsSolve.cpp:43
MumpsFinalize
void MumpsFinalize(DMUMPS_STRUC_C &theMumpsStruc)
Definition: MumpsSolve.cpp:77
MumpsFactorize
void MumpsFactorize(DMUMPS_STRUC_C &theMumpsStruc)
Definition: MumpsSolve.cpp:65
Generated on Thu Jul 2 2020 08:09:25 for Ice Sheet System Model by   doxygen 1.8.19