From a92b945e3f826add6402e8cf665073c99da4aefd Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 15 Sep 2018 11:40:07 +0200 Subject: [PATCH] does not crash anymore xx is nc, the number of components needs to be 3 in the current case --- src/FEM_mech.f90 | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/src/FEM_mech.f90 b/src/FEM_mech.f90 index 1163f470a..52e8a859c 100644 --- a/src/FEM_mech.f90 +++ b/src/FEM_mech.f90 @@ -5,6 +5,7 @@ !> @brief FEM PETSc solver !-------------------------------------------------------------------------------------------------- module FEM_mech +use PETSC #include #include #include @@ -99,7 +100,7 @@ subroutine FEM_mech_init(fieldBC) DMLabel :: BCLabel PetscInt, allocatable, target :: numComp(:), numDoF(:), bcField(:) PetscInt, pointer :: pNumComp(:), pNumDof(:), pBcField(:), pBcPoint(:) - PetscInt :: numBC, bcSize, xx + PetscInt :: numBC, bcSize, nc IS :: bcPoint IS, allocatable, target :: bcComps(:), bcPoints(:) IS, pointer :: pBcComps(:), pBcPoints(:) @@ -137,8 +138,9 @@ subroutine FEM_mech_init(fieldBC) qWeightsP => qWeights call PetscQuadratureCreate(PETSC_COMM_SELF,mechQuad,ierr); CHKERRQ(ierr) CHKERRQ(ierr) - call PetscQuadratureSetData(mechQuad,dimPlex,xx,nQuadrature,qPointsP,qWeightsP,ierr) - ! description for xx not on http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DT/PetscQuadratureSetData.html + ! what is the number of components, nc? + nc = 3 + call PetscQuadratureSetData(mechQuad,dimPlex,nc,nQuadrature,qPointsP,qWeightsP,ierr) CHKERRQ(ierr) call PetscFECreateDefault(mech_mesh,dimPlex,dimPlex,PETSC_TRUE,prefix, & integrationOrder,mechFE,ierr); CHKERRQ(ierr) @@ -253,8 +255,7 @@ subroutine FEM_mech_init(fieldBC) do basis = 0, nBasis-1 call PetscDualSpaceGetFunctional(mechDualSpace,basis,functional,ierr) CHKERRQ(ierr) - call PetscQuadratureGetData(functional,dimPlex,xx,nNodalPoints,nodalPointsP,nodalWeightsP,ierr) - ! description for xx not on http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DT/PetscQuadratureSetData.html + call PetscQuadratureGetData(functional,dimPlex,nc,nNodalPoints,nodalPointsP,nodalWeightsP,ierr) CHKERRQ(ierr) x_scal(basis*dimPlex+1:(basis+1)*dimPlex) = pV0 + matmul(transpose(cellJMat),nodalPointsP + 1.0) enddo