| 1 | /*!\file: NewMat.cpp
|
|---|
| 2 | * \brief create matrix using the Petsc library
|
|---|
| 3 | */
|
|---|
| 4 |
|
|---|
| 5 | #ifdef HAVE_CONFIG_H
|
|---|
| 6 | #include <config.h>
|
|---|
| 7 | #else
|
|---|
| 8 | #error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
|
|---|
| 9 | #endif
|
|---|
| 10 |
|
|---|
| 11 | /*Petsc includes: */
|
|---|
| 12 | #include <petscversion.h>
|
|---|
| 13 | #include <petscmat.h>
|
|---|
| 14 | #include <petscvec.h>
|
|---|
| 15 | #include <petscksp.h>
|
|---|
| 16 |
|
|---|
| 17 | #include "./petscpatches.h"
|
|---|
| 18 | #include "../../../shared/shared.h"
|
|---|
| 19 | #include "../../mpi/issmmpi.h"
|
|---|
| 20 |
|
|---|
| 21 | /*NewMat(int M,int N){{{*/
|
|---|
| 22 | Mat NewMat(int M,int N,ISSM_MPI_Comm comm){
|
|---|
| 23 |
|
|---|
| 24 | /*output:*/
|
|---|
| 25 | Mat outmatrix=NULL;
|
|---|
| 26 |
|
|---|
| 27 | /*parameters: */
|
|---|
| 28 | double sparsity=0.001; //default
|
|---|
| 29 | int m,n;
|
|---|
| 30 | int d_nz,o_nz,nnz;
|
|---|
| 31 |
|
|---|
| 32 | /*Determine local sizes: */
|
|---|
| 33 | m=DetermineLocalSize(M,comm);
|
|---|
| 34 | n=DetermineLocalSize(N,comm);
|
|---|
| 35 |
|
|---|
| 36 | nnz=(int)((double)M*(double)N*sparsity); //number of non zeros.
|
|---|
| 37 | d_nz=(int)((double)nnz/(double)M/2.0); //number of non zeros per row/2
|
|---|
| 38 | o_nz=(int)((double)nnz/(double)M/2.0); //number of non zeros per row/2
|
|---|
| 39 |
|
|---|
| 40 | #if PETSC_VERSION_GT(3,2,0)
|
|---|
| 41 | MatCreateAIJ(comm,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix);
|
|---|
| 42 | #else
|
|---|
| 43 | MatCreateMPIAIJ(comm,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix);
|
|---|
| 44 | #endif
|
|---|
| 45 |
|
|---|
| 46 | return outmatrix;
|
|---|
| 47 | }
|
|---|
| 48 | /*}}}*/
|
|---|
| 49 | /*NewMat(int M,int N,double sparsity,ISSM_MPI_Comm comm){{{*/
|
|---|
| 50 | Mat NewMat(int M,int N,double sparsity,ISSM_MPI_Comm comm){
|
|---|
| 51 |
|
|---|
| 52 | /*output:*/
|
|---|
| 53 | Mat outmatrix=NULL;
|
|---|
| 54 |
|
|---|
| 55 | /*parameters: */
|
|---|
| 56 | int m,n;
|
|---|
| 57 | int d_nz,o_nz;
|
|---|
| 58 | int nnz;
|
|---|
| 59 |
|
|---|
| 60 | /*Determine local sizes: */
|
|---|
| 61 | m=DetermineLocalSize(M,comm);
|
|---|
| 62 | n=DetermineLocalSize(N,comm);
|
|---|
| 63 |
|
|---|
| 64 | nnz=(int)((double)M*(double)N*sparsity); //number of non zeros.
|
|---|
| 65 | d_nz=(int)((double)nnz/(double)M/2.0); //number of non zeros per row/2
|
|---|
| 66 | o_nz=(int)((double)nnz/(double)M/2.0); //number of non zeros per row/2
|
|---|
| 67 |
|
|---|
| 68 | #if PETSC_VERSION_GT(3,2,0)
|
|---|
| 69 | if(sparsity==1){
|
|---|
| 70 | MatCreateDense(comm,m,n,M,N,NULL,&outmatrix);
|
|---|
| 71 | }
|
|---|
| 72 | else{
|
|---|
| 73 | MatCreateAIJ(comm,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix);
|
|---|
| 74 | }
|
|---|
| 75 | #else
|
|---|
| 76 | MatCreateMPIAIJ(comm,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix);
|
|---|
| 77 | #endif
|
|---|
| 78 |
|
|---|
| 79 | return outmatrix;
|
|---|
| 80 | }
|
|---|
| 81 | /*}}}*/
|
|---|
| 82 | /*NewMat(int M,int N,int connectivity,int numberofdofspernode){{{*/
|
|---|
| 83 | Mat NewMat(int M,int N,int connectivity,int numberofdofspernode,ISSM_MPI_Comm comm){
|
|---|
| 84 |
|
|---|
| 85 | /*output:*/
|
|---|
| 86 | Mat outmatrix=NULL;
|
|---|
| 87 |
|
|---|
| 88 | /*parameters: */
|
|---|
| 89 | int m,n;
|
|---|
| 90 | int d_nz,o_nz;
|
|---|
| 91 |
|
|---|
| 92 | #if PETSC_VERSION_MAJOR >= 3
|
|---|
| 93 | #if defined(_HAVE_PETSCDEV_) || PETSC_VERSION_MINOR >=4
|
|---|
| 94 | MatType type;
|
|---|
| 95 | #else
|
|---|
| 96 | const MatType type;
|
|---|
| 97 | #endif
|
|---|
| 98 | #else
|
|---|
| 99 | MatType type;
|
|---|
| 100 | #endif
|
|---|
| 101 |
|
|---|
| 102 | /*Determine local sizes: */
|
|---|
| 103 | m=DetermineLocalSize(M,comm);
|
|---|
| 104 | n=DetermineLocalSize(N,comm);
|
|---|
| 105 |
|
|---|
| 106 | /*Figure out number of non zeros per row: */
|
|---|
| 107 | d_nz=(int)connectivity*numberofdofspernode/2;
|
|---|
| 108 | o_nz=(int)connectivity*numberofdofspernode/2;
|
|---|
| 109 |
|
|---|
| 110 | MatCreate(comm,&outmatrix);
|
|---|
| 111 | MatSetSizes(outmatrix,m,n,M,N);
|
|---|
| 112 | MatSetFromOptions(outmatrix);
|
|---|
| 113 |
|
|---|
| 114 | /*preallocation according to type: */
|
|---|
| 115 | MatGetType(outmatrix,&type);
|
|---|
| 116 |
|
|---|
| 117 | if((strcmp(type,"mpiaij")==0) || (strcmp(type,"mpidense")==0)){
|
|---|
| 118 | MatMPIAIJSetPreallocation(outmatrix,d_nz,NULL,o_nz,NULL);
|
|---|
| 119 | }
|
|---|
| 120 |
|
|---|
| 121 | return outmatrix;
|
|---|
| 122 | }
|
|---|
| 123 | /*}}}*/
|
|---|