Actual source code: ex26f90.F90
1: program ex62f90
2: #include "petsc/finclude/petscdmplex.h"
3: use petscdmplex
4: implicit none
5: #include "exodusII.inc"
7: ! Get the Fortran kind associated with PetscInt and PetscReal so that we can use literal constants.
8: PetscInt :: dummyPetscInt
9: PetscReal :: dummyPetscreal
10: integer, parameter :: kPI = kind(dummyPetscInt)
11: integer, parameter :: kPR = kind(dummyPetscReal)
13: type(tDM) :: dm, dmU, dmA, dmS, dmUA, dmUA2, pDM
14: type(tDM), dimension(:), pointer :: dmList
15: type(tVec) :: X, U, A, S, UA, UA2
16: type(tIS) :: isU, isA, isS, isUA
17: type(tPetscSection) :: section
18: PetscInt :: fieldU = 0
19: PetscInt :: fieldA = 2
20: PetscInt :: fieldS = 1
21: PetscInt, dimension(2) :: fieldUA = [0, 2]
22: character(len=PETSC_MAX_PATH_LEN) :: ifilename, ofilename, IOBuffer
23: integer :: exoid = -1
24: type(tIS) :: csIS
25: PetscInt, dimension(:), pointer :: csID
26: PetscInt, dimension(:), pointer :: pStartDepth, pEndDepth
27: PetscInt :: order = 1
28: Integer :: i
29: PetscInt :: sdim, d, pStart, pEnd, p, numCS, set, j
30: PetscMPIInt :: rank, numProc
31: PetscBool :: flg
32: PetscErrorCode :: ierr
33: MPI_Comm :: comm
34: type(tPetscViewer) :: viewer
36: Character(len=MXSTLN) :: sJunk
37: PetscInt :: numstep = 3, step
38: Integer :: numNodalVar, numZonalVar
39: character(len=MXNAME), dimension(4) :: nodalVarName2D = ["U_x ", &
40: "U_y ", &
41: "Alpha", &
42: "Beta "]
43: character(len=MXNAME), dimension(3) :: zonalVarName2D = ["Sigma_11", &
44: "Sigma_12", &
45: "Sigma_22"]
46: character(len=MXNAME), dimension(5) :: nodalVarName3D = ["U_x ", &
47: "U_y ", &
48: "U_z ", &
49: "Alpha", &
50: "Beta "]
51: character(len=MXNAME), dimension(6) :: zonalVarName3D = ["Sigma_11", &
52: "Sigma_22", &
53: "Sigma_33", &
54: "Sigma_23", &
55: "Sigma_13", &
56: "Sigma_12"]
57: logical, dimension(:, :), pointer :: truthtable
59: type(tIS) :: cellIS
60: PetscInt, dimension(:), pointer :: cellID
61: PetscInt :: numCells, cell, closureSize
62: PetscInt, dimension(:), pointer :: closureA, closure
64: type(tPetscSection) :: sectionUA, coordSection
65: type(tVec) :: UALoc, coord
66: PetscScalar, dimension(:), pointer :: cval, xyz
67: PetscInt :: dofUA, offUA, c
69: ! dof layout ordered by increasing height in the DAG: cell, face, edge, vertex
70: PetscInt, dimension(3), target :: dofS2D = [0, 0, 3]
71: PetscInt, dimension(3), target :: dofUP1Tri = [2, 0, 0]
72: PetscInt, dimension(3), target :: dofAP1Tri = [1, 0, 0]
73: PetscInt, dimension(3), target :: dofUP2Tri = [2, 2, 0]
74: PetscInt, dimension(3), target :: dofAP2Tri = [1, 1, 0]
75: PetscInt, dimension(3), target :: dofUP1Quad = [2, 0, 0]
76: PetscInt, dimension(3), target :: dofAP1Quad = [1, 0, 0]
77: PetscInt, dimension(3), target :: dofUP2Quad = [2, 2, 2]
78: PetscInt, dimension(3), target :: dofAP2Quad = [1, 1, 1]
79: PetscInt, dimension(4), target :: dofS3D = [0, 0, 0, 6]
80: PetscInt, dimension(4), target :: dofUP1Tet = [3, 0, 0, 0]
81: PetscInt, dimension(4), target :: dofAP1Tet = [1, 0, 0, 0]
82: PetscInt, dimension(4), target :: dofUP2Tet = [3, 3, 0, 0]
83: PetscInt, dimension(4), target :: dofAP2Tet = [1, 1, 0, 0]
84: PetscInt, dimension(4), target :: dofUP1Hex = [3, 0, 0, 0]
85: PetscInt, dimension(4), target :: dofAP1Hex = [1, 0, 0, 0]
86: PetscInt, dimension(4), target :: dofUP2Hex = [3, 3, 3, 3]
87: PetscInt, dimension(4), target :: dofAP2Hex = [1, 1, 1, 1]
88: PetscInt, dimension(:), pointer :: dofU, dofA, dofS
90: type(tPetscSF) :: migrationSF
91: PetscPartitioner :: part
92: PetscLayout :: layout
94: type(tVec) :: tmpVec
95: PetscReal :: norm
96: PetscReal :: time = 1.234_kPR
98: PetscCallA(PetscInitialize(PETSC_NULL_CHARACTER, ierr))
99: if (ierr /= 0) then
100: print *, 'Unable to initialize PETSc'
101: stop
102: end if
104: PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))
105: PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD, numProc, ierr))
106: PetscCallA(PetscOptionsGetString(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-i', ifilename, flg, ierr))
107: PetscCheckA(flg, PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, 'missing input file name -i ')
108: PetscCallA(PetscOptionsGetString(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-o', ofilename, flg, ierr))
109: PetscCheckA(flg, PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, 'missing output file name -o