1d8606c27SBarry Smith! 2d8606c27SBarry Smith! This program is modified from a user's contribution. 3d8606c27SBarry Smith! It illustrates how to PetscCallA(MUMPS's LU solver 4d8606c27SBarry Smith! 5d8606c27SBarry Smith 6d8606c27SBarry Smith program main 7d8606c27SBarry Smith#include <petsc/finclude/petscmat.h> 8d8606c27SBarry Smith use petscmat 9d8606c27SBarry Smith implicit none 10d8606c27SBarry Smith 11d8606c27SBarry Smith Vec x,b,u 12d8606c27SBarry Smith Mat A, fact 13d8606c27SBarry Smith PetscInt i,j,II,JJ,m 14d8606c27SBarry Smith PetscInt Istart, Iend 15d8606c27SBarry Smith PetscInt ione, ifive 16d8606c27SBarry Smith PetscBool wmumps 17d8606c27SBarry Smith PetscBool flg 18d8606c27SBarry Smith PetscScalar one, v 19d8606c27SBarry Smith IS perm,iperm 20d8606c27SBarry Smith PetscErrorCode ierr 21*ce78bad3SBarry Smith MatFactorInfo info 22d8606c27SBarry Smith 23d8606c27SBarry Smith PetscCallA(PetscInitialize(PETSC_NULL_CHARACTER, ierr)) 24d8606c27SBarry Smith m = 10 25d8606c27SBarry Smith one = 1.0 26d8606c27SBarry Smith ione = 1 27d8606c27SBarry Smith ifive = 5 28d8606c27SBarry Smith 29d8606c27SBarry Smith wmumps = PETSC_FALSE 30d8606c27SBarry Smith 31d8606c27SBarry Smith PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-m',m,flg, ierr)) 32d8606c27SBarry Smith PetscCallA(PetscOptionsGetBool(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-use_mumps',wmumps,flg,ierr)) 33d8606c27SBarry Smith 34d8606c27SBarry Smith PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr)) 35d8606c27SBarry Smith PetscCallA(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, m*m, m*m, ierr)) 36d8606c27SBarry Smith PetscCallA(MatSetType(A, MATAIJ, ierr)) 37d8606c27SBarry Smith PetscCallA(MatSetFromOptions(A, ierr)) 385d83a8b1SBarry Smith PetscCallA(MatSeqAIJSetPreallocation(A,ifive, PETSC_NULL_INTEGER_ARRAY, ierr)) 395d83a8b1SBarry Smith PetscCallA(MatMPIAIJSetPreallocation(A,ifive,PETSC_NULL_INTEGER_ARRAY,ifive,PETSC_NULL_INTEGER_ARRAY,ierr)) 40d8606c27SBarry Smith 41d8606c27SBarry Smith PetscCallA(MatGetOwnershipRange(A,Istart,Iend,ierr)) 42d8606c27SBarry Smith 43d8606c27SBarry Smith do 10, II=Istart,Iend - 1 44d8606c27SBarry Smith v = -1.0 45d8606c27SBarry Smith i = II/m 46d8606c27SBarry Smith j = II - i*m 47d8606c27SBarry Smith if (i.gt.0) then 48d8606c27SBarry Smith JJ = II - m 495d83a8b1SBarry Smith PetscCallA(MatSetValues(A,ione,[II],ione,[JJ],[v],INSERT_VALUES,ierr)) 50d8606c27SBarry Smith endif 51d8606c27SBarry Smith if (i.lt.m-1) then 52d8606c27SBarry Smith JJ = II + m 535d83a8b1SBarry Smith PetscCallA(MatSetValues(A,ione,[II],ione,[JJ],[v],INSERT_VALUES,ierr)) 54d8606c27SBarry Smith endif 55d8606c27SBarry Smith if (j.gt.0) then 56d8606c27SBarry Smith JJ = II - 1 575d83a8b1SBarry Smith PetscCallA(MatSetValues(A,ione,[II],ione,[JJ],[v],INSERT_VALUES,ierr)) 58d8606c27SBarry Smith endif 59d8606c27SBarry Smith if (j.lt.m-1) then 60d8606c27SBarry Smith JJ = II + 1 615d83a8b1SBarry Smith PetscCallA(MatSetValues(A,ione,[II],ione,[JJ],[v],INSERT_VALUES,ierr)) 62d8606c27SBarry Smith endif 63d8606c27SBarry Smith v = 4.0 645d83a8b1SBarry Smith PetscCallA(MatSetValues(A,ione,[II],ione,[II],[v],INSERT_VALUES,ierr)) 65d8606c27SBarry Smith 10 continue 66d8606c27SBarry Smith 67d8606c27SBarry Smith PetscCallA(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY, ierr)) 68d8606c27SBarry Smith PetscCallA(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY, ierr)) 69d8606c27SBarry Smith 70d8606c27SBarry Smith PetscCallA(VecCreate(PETSC_COMM_WORLD, u, ierr)) 71d8606c27SBarry Smith PetscCallA(VecSetSizes(u, PETSC_DECIDE, m*m, ierr)) 72d8606c27SBarry Smith PetscCallA(VecSetFromOptions(u, ierr)) 73d8606c27SBarry Smith PetscCallA(VecDuplicate(u,b,ierr)) 74d8606c27SBarry Smith PetscCallA(VecDuplicate(b,x,ierr)) 75d8606c27SBarry Smith PetscCallA(VecSet(u, one, ierr)) 76d8606c27SBarry Smith PetscCallA(MatMult(A, u, b, ierr)) 77d8606c27SBarry Smith 78d8606c27SBarry Smith PetscCallA(MatFactorInfoInitialize(info,ierr)) 79d8606c27SBarry Smith PetscCallA(MatGetOrdering(A,MATORDERINGNATURAL,perm,iperm,ierr)) 80d8606c27SBarry Smith if (wmumps) then 81d8606c27SBarry Smith write(*,*) 'use MUMPS LU...' 82d8606c27SBarry Smith PetscCallA(MatGetFactor(A,MATSOLVERMUMPS,MAT_FACTOR_LU,fact,ierr)) 83d8606c27SBarry Smith else 84d8606c27SBarry Smith write(*,*) 'use PETSc LU...' 85d8606c27SBarry Smith PetscCallA(MatGetFactor(A,MATSOLVERPETSC,MAT_FACTOR_LU,fact,ierr)) 86d8606c27SBarry Smith endif 87d8606c27SBarry Smith PetscCallA(MatLUFactorSymbolic(fact, A, perm, iperm, info, ierr)) 88d8606c27SBarry Smith PetscCallA(ISDestroy(perm,ierr)) 89d8606c27SBarry Smith PetscCallA(ISDestroy(iperm,ierr)) 90d8606c27SBarry Smith 91d8606c27SBarry Smith PetscCallA(MatLUFactorNumeric(fact, A, info, ierr)) 92d8606c27SBarry Smith PetscCallA(MatSolve(fact, b, x,ierr)) 93d8606c27SBarry Smith PetscCallA(MatDestroy(fact, ierr)) 94d8606c27SBarry Smith 95d8606c27SBarry Smith PetscCallA(MatDestroy(A, ierr)) 96d8606c27SBarry Smith PetscCallA(VecDestroy(u, ierr)) 97d8606c27SBarry Smith PetscCallA(VecDestroy(x, ierr)) 98d8606c27SBarry Smith PetscCallA(VecDestroy(b, ierr)) 99d8606c27SBarry Smith 100d8606c27SBarry Smith PetscCallA(PetscFinalize(ierr)) 101d8606c27SBarry Smith end 102d8606c27SBarry Smith 103d8606c27SBarry Smith!/*TEST 104d8606c27SBarry Smith! 105d8606c27SBarry Smith! test: 106d8606c27SBarry Smith! 107d8606c27SBarry Smith! test: 108d8606c27SBarry Smith! suffix: 2 109d8606c27SBarry Smith! args: -use_mumps 110d8606c27SBarry Smith! requires: mumps 111d8606c27SBarry Smith! 112d8606c27SBarry Smith!TEST*/ 113