diff --git a/src/grid/grid_damage_spectral.f90 b/src/grid/grid_damage_spectral.f90 index 48994ec15..96b72dbcc 100644 --- a/src/grid/grid_damage_spectral.f90 +++ b/src/grid/grid_damage_spectral.f90 @@ -112,8 +112,10 @@ subroutine grid_damage_spectral_init() !-------------------------------------------------------------------------------------------------- ! initialize solver specific parts of PETSc - call SNESCreate(PETSC_COMM_WORLD,SNES_damage,err_PETSc); CHKERRQ(err_PETSc) - call SNESSetOptionsPrefix(SNES_damage,'damage_',err_PETSc);CHKERRQ(err_PETSc) + call SNESCreate(PETSC_COMM_WORLD,SNES_damage,err_PETSc) + CHKERRQ(err_PETSc) + call SNESSetOptionsPrefix(SNES_damage,'damage_',err_PETSc) + CHKERRQ(err_PETSc) localK = 0_pPetscInt localK(worldrank) = int(grid3,pPetscInt) call MPI_Allreduce(MPI_IN_PLACE,localK,worldsize,MPI_INTEGER,MPI_SUM,MPI_COMM_WORLD,err_MPI) @@ -127,25 +129,39 @@ subroutine grid_damage_spectral_init() [int(grid(1),pPetscInt)],[int(grid(2),pPetscInt)],localK, & ! local grid damage_grid,err_PETSc) ! handle, error CHKERRQ(err_PETSc) - call DMsetFromOptions(damage_grid,err_PETSc); CHKERRQ(err_PETSc) - call DMsetUp(damage_grid,err_PETSc); CHKERRQ(err_PETSc) - call DMCreateGlobalVector(damage_grid,solution_vec,err_PETSc); CHKERRQ(err_PETSc) ! global solution vector (grid x 1, i.e. every def grad tensor) + call DMsetFromOptions(damage_grid,err_PETSc) + CHKERRQ(err_PETSc) + call DMsetUp(damage_grid,err_PETSc) + CHKERRQ(err_PETSc) + call DMCreateGlobalVector(damage_grid,solution_vec,err_PETSc) ! global solution vector (grid x 1, i.e. every def grad tensor) + CHKERRQ(err_PETSc) call DMDASNESSetFunctionLocal(damage_grid,INSERT_VALUES,formResidual,PETSC_NULL_SNES,err_PETSc) ! residual vector of same shape as solution vector CHKERRQ(err_PETSc) - call SNESSetDM(SNES_damage,damage_grid,err_PETSc); CHKERRQ(err_PETSc) - call SNESSetFromOptions(SNES_damage,err_PETSc); CHKERRQ(err_PETSc) ! pull it all together with additional CLI arguments - call SNESGetType(SNES_damage,snes_type,err_PETSc); CHKERRQ(err_PETSc) + call SNESSetDM(SNES_damage,damage_grid,err_PETSc) + CHKERRQ(err_PETSc) + call SNESSetFromOptions(SNES_damage,err_PETSc) ! pull it all together with additional CLI arguments + CHKERRQ(err_PETSc) + call SNESGetType(SNES_damage,snes_type,err_PETSc) + CHKERRQ(err_PETSc) if (trim(snes_type) == 'vinewtonrsls' .or. & trim(snes_type) == 'vinewtonssls') then - call DMGetGlobalVector(damage_grid,lBound,err_PETSc); CHKERRQ(err_PETSc) - call DMGetGlobalVector(damage_grid,uBound,err_PETSc); CHKERRQ(err_PETSc) - call VecSet(lBound,0.0_pReal,err_PETSc); CHKERRQ(err_PETSc) - call VecSet(uBound,1.0_pReal,err_PETSc); CHKERRQ(err_PETSc) + call DMGetGlobalVector(damage_grid,lBound,err_PETSc) + CHKERRQ(err_PETSc) + call DMGetGlobalVector(damage_grid,uBound,err_PETSc) + CHKERRQ(err_PETSc) + call VecSet(lBound,0.0_pReal,err_PETSc) + CHKERRQ(err_PETSc) + call VecSet(uBound,1.0_pReal,err_PETSc) + CHKERRQ(err_PETSc) call SNESVISetVariableBounds(SNES_damage,lBound,uBound,err_PETSc) ! variable bounds for variational inequalities - call DMRestoreGlobalVector(damage_grid,lBound,err_PETSc); CHKERRQ(err_PETSc) - call DMRestoreGlobalVector(damage_grid,uBound,err_PETSc); CHKERRQ(err_PETSc) + CHKERRQ(err_PETSc) + call DMRestoreGlobalVector(damage_grid,lBound,err_PETSc) + CHKERRQ(err_PETSc) + call DMRestoreGlobalVector(damage_grid,uBound,err_PETSc) + CHKERRQ(err_PETSc) end if - call VecSet(solution_vec,1.0_pReal,err_PETSc); CHKERRQ(err_PETSc) + call VecSet(solution_vec,1.0_pReal,err_PETSc) + CHKERRQ(err_PETSc) call updateReference() @@ -202,8 +218,10 @@ function grid_damage_spectral_solution(Delta_t) result(solution) call homogenization_set_phi(phi_current(i,j,k),ce) end do; end do; end do - call VecMin(solution_vec,devNull,phi_min,err_PETSc); CHKERRQ(err_PETSc) - call VecMax(solution_vec,devNull,phi_max,err_PETSc); CHKERRQ(err_PETSc) + call VecMin(solution_vec,devNull,phi_min,err_PETSc) + CHKERRQ(err_PETSc) + call VecMax(solution_vec,devNull,phi_max,err_PETSc) + CHKERRQ(err_PETSc) if (solution%converged) & print'(/,1x,a)', '... nonlocal damage converged .....................................' print'(/,1x,a,f8.6,2x,f8.6,2x,e11.4)', 'Minimum|Maximum|Delta Damage = ', phi_min, phi_max, stagNorm @@ -229,7 +247,8 @@ subroutine grid_damage_spectral_forward(cutBack) phi_stagInc = phi_lastInc !-------------------------------------------------------------------------------------------------- ! reverting damage field state - call SNESGetDM(SNES_damage,dm_local,err_PETSc); CHKERRQ(err_PETSc) + call SNESGetDM(SNES_damage,dm_local,err_PETSc) + CHKERRQ(err_PETSc) call DMDAVecGetArrayF90(dm_local,solution_vec,phi_PETSc,err_PETSc) !< get the data out of PETSc to work with CHKERRQ(err_PETSc) phi_PETSc = phi_current diff --git a/src/grid/grid_mech_FEM.f90 b/src/grid/grid_mech_FEM.f90 index fb3ec85f6..1dc421f09 100644 --- a/src/grid/grid_mech_FEM.f90 +++ b/src/grid/grid_mech_FEM.f90 @@ -203,9 +203,12 @@ subroutine grid_mechanical_FEM_init !-------------------------------------------------------------------------------------------------- ! init fields - call VecSet(solution_current,0.0_pReal,err_PETSc);CHKERRQ(err_PETSc) - call VecSet(solution_lastInc,0.0_pReal,err_PETSc);CHKERRQ(err_PETSc) - call VecSet(solution_rate ,0.0_pReal,err_PETSc);CHKERRQ(err_PETSc) + call VecSet(solution_current,0.0_pReal,err_PETSc) + CHKERRQ(err_PETSc) + call VecSet(solution_lastInc,0.0_pReal,err_PETSc) + CHKERRQ(err_PETSc) + call VecSet(solution_rate ,0.0_pReal,err_PETSc) + CHKERRQ(err_PETSc) call DMDAVecGetArrayF90(mechanical_grid,solution_current,u_current,err_PETSc) CHKERRQ(err_PETSc) call DMDAVecGetArrayF90(mechanical_grid,solution_lastInc,u_lastInc,err_PETSc) @@ -375,9 +378,11 @@ subroutine grid_mechanical_FEM_forward(cutBack,guess,Delta_t,Delta_t_old,t_remai call VecScale(solution_rate,1.0_pReal/Delta_t_old,err_PETSc) CHKERRQ(err_PETSc) else - call VecSet(solution_rate,0.0_pReal,err_PETSc); CHKERRQ(err_PETSc) + call VecSet(solution_rate,0.0_pReal,err_PETSc) + CHKERRQ(err_PETSc) endif - call VecCopy(solution_current,solution_lastInc,err_PETSc); CHKERRQ(err_PETSc) + call VecCopy(solution_current,solution_lastInc,err_PETSc) + CHKERRQ(err_PETSc) F_lastInc = F @@ -549,7 +554,8 @@ subroutine formResidual(da_local,x_local, & !-------------------------------------------------------------------------------------------------- ! get deformation gradient - call DMDAVecGetArrayF90(da_local,x_local,x_scal,err_PETSc);CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da_local,x_local,x_scal,err_PETSc) + CHKERRQ(err_PETSc) do k = grid3offset+1, grid3offset+grid3; do j = 1, grid(2); do i = 1, grid(1) ctr = 0 do kk = -1, 0; do jj = -1, 0; do ii = -1, 0 @@ -558,7 +564,8 @@ subroutine formResidual(da_local,x_local, & enddo; enddo; enddo F(1:3,1:3,i,j,k-grid3offset) = params%rotation_BC%rotate(F_aim,active=.true.) + transpose(matmul(BMat,x_elem)) enddo; enddo; enddo - call DMDAVecRestoreArrayF90(da_local,x_local,x_scal,err_PETSc);CHKERRQ(err_PETSc) + call DMDAVecRestoreArrayF90(da_local,x_local,x_scal,err_PETSc) + CHKERRQ(err_PETSc) !-------------------------------------------------------------------------------------------------- ! evaluate constitutive response @@ -575,9 +582,12 @@ subroutine formResidual(da_local,x_local, & !-------------------------------------------------------------------------------------------------- ! constructing residual - call VecSet(f_local,0.0_pReal,err_PETSc);CHKERRQ(err_PETSc) - call DMDAVecGetArrayF90(da_local,f_local,r,err_PETSc);CHKERRQ(err_PETSc) - call DMDAVecGetArrayF90(da_local,x_local,x_scal,err_PETSc);CHKERRQ(err_PETSc) + call VecSet(f_local,0.0_pReal,err_PETSc) + CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da_local,f_local,r,err_PETSc) + CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da_local,x_local,x_scal,err_PETSc) + CHKERRQ(err_PETSc) ele = 0 do k = grid3offset+1, grid3offset+grid3; do j = 1, grid(2); do i = 1, grid(1) ctr = 0 @@ -596,12 +606,15 @@ subroutine formResidual(da_local,x_local, & r(0:2,i+ii,j+jj,k+kk) = r(0:2,i+ii,j+jj,k+kk) + f_elem(ctr,1:3) enddo; enddo; enddo enddo; enddo; enddo - call DMDAVecRestoreArrayF90(da_local,x_local,x_scal,err_PETSc);CHKERRQ(err_PETSc) - call DMDAVecRestoreArrayF90(da_local,f_local,r,err_PETSc);CHKERRQ(err_PETSc) + call DMDAVecRestoreArrayF90(da_local,x_local,x_scal,err_PETSc) + CHKERRQ(err_PETSc) + call DMDAVecRestoreArrayF90(da_local,f_local,r,err_PETSc) + CHKERRQ(err_PETSc) !-------------------------------------------------------------------------------------------------- ! applying boundary conditions - call DMDAVecGetArrayF90(da_local,f_local,r,err_PETSc);CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da_local,f_local,r,err_PETSc) + CHKERRQ(err_PETSc) if (grid3offset == 0) then r(0:2,0, 0, 0) = 0.0_pReal r(0:2,grid(1),0, 0) = 0.0_pReal @@ -614,7 +627,8 @@ subroutine formResidual(da_local,x_local, & r(0:2,0, grid(2),grid(3)) = 0.0_pReal r(0:2,grid(1),grid(2),grid(3)) = 0.0_pReal end if - call DMDAVecRestoreArrayF90(da_local,f_local,r,err_PETSc);CHKERRQ(err_PETSc) + call DMDAVecRestoreArrayF90(da_local,f_local,r,err_PETSc) + CHKERRQ(err_PETSc) end subroutine formResidual @@ -646,7 +660,8 @@ subroutine formJacobian(da_local,x_local,Jac_pre,Jac,dummy,err_PETSc) CHKERRQ(err_PETSc) call MatSetOption(Jac,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE,err_PETSc) CHKERRQ(err_PETSc) - call MatZeroEntries(Jac,err_PETSc); CHKERRQ(err_PETSc) + call MatZeroEntries(Jac,err_PETSc) + CHKERRQ(err_PETSc) ce = 0 do k = grid3offset+1, grid3offset+grid3; do j = 1, grid(2); do i = 1, grid(1) ctr = 0 @@ -684,10 +699,14 @@ subroutine formJacobian(da_local,x_local,Jac_pre,Jac,dummy,err_PETSc) call MatSetValuesStencil(Jac,24_pPETScInt,row,24_pPetscInt,col,K_ele,ADD_VALUES,err_PETSc) CHKERRQ(err_PETSc) enddo; enddo; enddo - call MatAssemblyBegin(Jac,MAT_FINAL_ASSEMBLY,err_PETSc); CHKERRQ(err_PETSc) - call MatAssemblyEnd(Jac,MAT_FINAL_ASSEMBLY,err_PETSc); CHKERRQ(err_PETSc) - call MatAssemblyBegin(Jac_pre,MAT_FINAL_ASSEMBLY,err_PETSc); CHKERRQ(err_PETSc) - call MatAssemblyEnd(Jac_pre,MAT_FINAL_ASSEMBLY,err_PETSc); CHKERRQ(err_PETSc) + call MatAssemblyBegin(Jac,MAT_FINAL_ASSEMBLY,err_PETSc) + CHKERRQ(err_PETSc) + call MatAssemblyEnd(Jac,MAT_FINAL_ASSEMBLY,err_PETSc) + CHKERRQ(err_PETSc) + call MatAssemblyBegin(Jac_pre,MAT_FINAL_ASSEMBLY,err_PETSc) + CHKERRQ(err_PETSc) + call MatAssemblyEnd(Jac_pre,MAT_FINAL_ASSEMBLY,err_PETSc) + CHKERRQ(err_PETSc) !-------------------------------------------------------------------------------------------------- ! applying boundary conditions diff --git a/src/grid/grid_mech_spectral_basic.f90 b/src/grid/grid_mech_spectral_basic.f90 index 00ff756ec..c60c542d0 100644 --- a/src/grid/grid_mech_spectral_basic.f90 +++ b/src/grid/grid_mech_spectral_basic.f90 @@ -175,19 +175,25 @@ subroutine grid_mechanical_spectral_basic_init [int(grid(1),pPetscInt)],[int(grid(2),pPetscInt)],localK, & ! local grid da,err_PETSc) ! handle, error CHKERRQ(err_PETSc) - call DMsetFromOptions(da,err_PETSc); CHKERRQ(err_PETSc) - call DMsetUp(da,err_PETSc); CHKERRQ(err_PETSc) - call DMcreateGlobalVector(da,solution_vec,err_PETSc); CHKERRQ(err_PETSc) ! global solution vector (grid x 9, i.e. every def grad tensor) + call DMsetFromOptions(da,err_PETSc) + CHKERRQ(err_PETSc) + call DMsetUp(da,err_PETSc) + CHKERRQ(err_PETSc) + call DMcreateGlobalVector(da,solution_vec,err_PETSc) ! global solution vector (grid x 9, i.e. every def grad tensor) + CHKERRQ(err_PETSc) call DMDASNESsetFunctionLocal(da,INSERT_VALUES,formResidual,PETSC_NULL_SNES,err_PETSc) ! residual vector of same shape as solution vector CHKERRQ(err_PETSc) call SNESsetConvergenceTest(SNES_mechanical,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,err_PETSc) ! specify custom convergence check function "converged" CHKERRQ(err_PETSc) - call SNESSetDM(SNES_mechanical,da,err_PETSc); CHKERRQ(err_PETSc) - call SNESsetFromOptions(SNES_mechanical,err_PETSc); CHKERRQ(err_PETSc) ! pull it all together with additional CLI arguments + call SNESSetDM(SNES_mechanical,da,err_PETSc) + CHKERRQ(err_PETSc) + call SNESsetFromOptions(SNES_mechanical,err_PETSc) ! pull it all together with additional CLI arguments + CHKERRQ(err_PETSc) !-------------------------------------------------------------------------------------------------- ! init fields - call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc); CHKERRQ(err_PETSc) ! places pointer on PETSc data + call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc) ! places pointer on PETSc data + CHKERRQ(err_PETSc) restartRead: if (interface_restartInc > 0) then print'(/,1x,a,i0,a)', 'reading restart data of increment ', interface_restartInc, ' from file' @@ -220,7 +226,8 @@ subroutine grid_mechanical_spectral_basic_init call utilities_constitutiveResponse(P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 reshape(F,shape(F_lastInc)), & ! target F 0.0_pReal) ! time increment - call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc); CHKERRQ(err_PETSc) ! deassociate pointer + call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc) ! deassociate pointer + CHKERRQ(err_PETSc) restartRead2: if (interface_restartInc > 0) then print'(1x,a,i0,a)', 'reading more restart data of increment ', interface_restartInc, ' from file' @@ -309,7 +316,8 @@ subroutine grid_mechanical_spectral_basic_forward(cutBack,guess,Delta_t,Delta_t_ PetscScalar, pointer, dimension(:,:,:,:) :: F - call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc); CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc) + CHKERRQ(err_PETSc) if (cutBack) then C_volAvg = C_volAvgLastInc @@ -352,7 +360,8 @@ subroutine grid_mechanical_spectral_basic_forward(cutBack,guess,Delta_t,Delta_t_ F = reshape(utilities_forwardField(Delta_t,F_lastInc,Fdot, & ! estimate of F at end of time+Delta_t that matches rotated F_aim on average rotation_BC%rotate(F_aim,active=.true.)),[9,grid(1),grid(2),grid3]) - call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc); CHKERRQ(err_PETSc) + call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc) + CHKERRQ(err_PETSc) !-------------------------------------------------------------------------------------------------- ! set module wide available data @@ -371,9 +380,11 @@ subroutine grid_mechanical_spectral_basic_updateCoords PetscErrorCode :: err_PETSc PetscScalar, dimension(:,:,:,:), pointer :: F - call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc); CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc) + CHKERRQ(err_PETSc) call utilities_updateCoords(F) - call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc); CHKERRQ(err_PETSc) + call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc) + CHKERRQ(err_PETSc) end subroutine grid_mechanical_spectral_basic_updateCoords @@ -387,7 +398,8 @@ subroutine grid_mechanical_spectral_basic_restartWrite integer(HID_T) :: fileHandle, groupHandle PetscScalar, dimension(:,:,:,:), pointer :: F - call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc); CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc) + CHKERRQ(err_PETSc) print'(1x,a)', 'writing solver data required for restart to file'; flush(IO_STDOUT) @@ -414,7 +426,8 @@ subroutine grid_mechanical_spectral_basic_restartWrite if (num%update_gamma) call utilities_saveReferenceStiffness - call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc); CHKERRQ(err_PETSc) + call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc) + CHKERRQ(err_PETSc) end subroutine grid_mechanical_spectral_basic_restartWrite diff --git a/src/grid/grid_mech_spectral_polarisation.f90 b/src/grid/grid_mech_spectral_polarisation.f90 index 52c308ae9..810047403 100644 --- a/src/grid/grid_mech_spectral_polarisation.f90 +++ b/src/grid/grid_mech_spectral_polarisation.f90 @@ -178,7 +178,8 @@ subroutine grid_mechanical_spectral_polarisation_init !-------------------------------------------------------------------------------------------------- ! initialize solver specific parts of PETSc - call SNESCreate(PETSC_COMM_WORLD,SNES_mechanical,err_PETSc); CHKERRQ(err_PETSc) + call SNESCreate(PETSC_COMM_WORLD,SNES_mechanical,err_PETSc) + CHKERRQ(err_PETSc) call SNESSetOptionsPrefix(SNES_mechanical,'mechanical_',err_PETSc) CHKERRQ(err_PETSc) localK = 0_pPetscInt @@ -194,19 +195,25 @@ subroutine grid_mechanical_spectral_polarisation_init [int(grid(1),pPetscInt)],[int(grid(2),pPetscInt)],localK, & ! local grid da,err_PETSc) ! handle, error CHKERRQ(err_PETSc) - call DMsetFromOptions(da,err_PETSc); CHKERRQ(err_PETSc) - call DMsetUp(da,err_PETSc); CHKERRQ(err_PETSc) - call DMcreateGlobalVector(da,solution_vec,err_PETSc); CHKERRQ(err_PETSc) ! global solution vector (grid x 18, i.e. every def grad tensor) + call DMsetFromOptions(da,err_PETSc) + CHKERRQ(err_PETSc) + call DMsetUp(da,err_PETSc) + CHKERRQ(err_PETSc) + call DMcreateGlobalVector(da,solution_vec,err_PETSc) ! global solution vector (grid x 18, i.e. every def grad tensor) + CHKERRQ(err_PETSc) call DMDASNESsetFunctionLocal(da,INSERT_VALUES,formResidual,PETSC_NULL_SNES,err_PETSc) ! residual vector of same shape as solution vector CHKERRQ(err_PETSc) call SNESsetConvergenceTest(SNES_mechanical,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,err_PETSc) ! specify custom convergence check function "converged" CHKERRQ(err_PETSc) - call SNESSetDM(SNES_mechanical,da,err_PETSc); CHKERRQ(err_PETSc) - call SNESsetFromOptions(SNES_mechanical,err_PETSc); CHKERRQ(err_PETSc) ! pull it all together with additional CLI arguments + call SNESSetDM(SNES_mechanical,da,err_PETSc) + CHKERRQ(err_PETSc) + call SNESsetFromOptions(SNES_mechanical,err_PETSc) ! pull it all together with additional CLI arguments + CHKERRQ(err_PETSc) !-------------------------------------------------------------------------------------------------- ! init fields - call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc); CHKERRQ(err_PETSc) ! places pointer on PETSc data + call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc) ! places pointer on PETSc data + CHKERRQ(err_PETSc) F => FandF_tau(0: 8,:,:,:) F_tau => FandF_tau(9:17,:,:,:) @@ -343,7 +350,8 @@ subroutine grid_mechanical_spectral_polarisation_forward(cutBack,guess,Delta_t,D real(pReal), dimension(3,3) :: F_lambda33 - call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc); CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc) + CHKERRQ(err_PETSc) F => FandF_tau(0: 8,:,:,:) F_tau => FandF_tau(9:17,:,:,:) @@ -444,7 +452,8 @@ subroutine grid_mechanical_spectral_polarisation_restartWrite integer(HID_T) :: fileHandle, groupHandle PetscScalar, dimension(:,:,:,:), pointer :: FandF_tau, F, F_tau - call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc); CHKERRQ(err_PETSc) + call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc) + CHKERRQ(err_PETSc) F => FandF_tau(0: 8,:,:,:) F_tau => FandF_tau(9:17,:,:,:) diff --git a/src/grid/grid_thermal_spectral.f90 b/src/grid/grid_thermal_spectral.f90 index 9c6d10cca..f102cf2c1 100644 --- a/src/grid/grid_thermal_spectral.f90 +++ b/src/grid/grid_thermal_spectral.f90 @@ -113,8 +113,10 @@ subroutine grid_thermal_spectral_init(T_0) !-------------------------------------------------------------------------------------------------- ! initialize solver specific parts of PETSc - call SNESCreate(PETSC_COMM_WORLD,SNES_thermal,err_PETSc); CHKERRQ(err_PETSc) - call SNESSetOptionsPrefix(SNES_thermal,'thermal_',err_PETSc);CHKERRQ(err_PETSc) + call SNESCreate(PETSC_COMM_WORLD,SNES_thermal,err_PETSc) + CHKERRQ(err_PETSc) + call SNESSetOptionsPrefix(SNES_thermal,'thermal_',err_PETSc) + CHKERRQ(err_PETSc) localK = 0_pPetscInt localK(worldrank) = int(grid3,pPetscInt) call MPI_Allreduce(MPI_IN_PLACE,localK,worldsize,MPI_INTEGER,MPI_SUM,MPI_COMM_WORLD,err_MPI) @@ -128,14 +130,18 @@ subroutine grid_thermal_spectral_init(T_0) [int(grid(1),pPetscInt)],[int(grid(2),pPetscInt)],localK, & ! local grid thermal_grid,err_PETSc) ! handle, error CHKERRQ(err_PETSc) - call DMsetFromOptions(thermal_grid,err_PETSc); CHKERRQ(err_PETSc) - call DMsetUp(thermal_grid,err_PETSc); CHKERRQ(err_PETSc) + call DMsetFromOptions(thermal_grid,err_PETSc) + CHKERRQ(err_PETSc) + call DMsetUp(thermal_grid,err_PETSc) + CHKERRQ(err_PETSc) call DMCreateGlobalVector(thermal_grid,solution_vec,err_PETSc) ! global solution vector (grid x 1, i.e. every def grad tensor) CHKERRQ(err_PETSc) call DMDASNESSetFunctionLocal(thermal_grid,INSERT_VALUES,formResidual,PETSC_NULL_SNES,err_PETSc) ! residual vector of same shape as solution vector CHKERRQ(err_PETSc) - call SNESSetDM(SNES_thermal,thermal_grid,err_PETSc); CHKERRQ(err_PETSc) - call SNESSetFromOptions(SNES_thermal,err_PETSc); CHKERRQ(err_PETSc) ! pull it all together with additional CLI arguments + call SNESSetDM(SNES_thermal,thermal_grid,err_PETSc) + CHKERRQ(err_PETSc) + call SNESSetFromOptions(SNES_thermal,err_PETSc) ! pull it all together with additional CLI arguments + CHKERRQ(err_PETSc) call DMDAVecGetArrayF90(thermal_grid,solution_vec,T_PETSc,err_PETSc) CHKERRQ(err_PETSc) T_PETSc = T_current @@ -197,8 +203,10 @@ function grid_thermal_spectral_solution(Delta_t) result(solution) call homogenization_thermal_setField(T_current(i,j,k),(T_current(i,j,k)-T_lastInc(i,j,k))/params%Delta_t,ce) end do; end do; end do - call VecMin(solution_vec,devNull,T_min,err_PETSc); CHKERRQ(err_PETSc) - call VecMax(solution_vec,devNull,T_max,err_PETSc); CHKERRQ(err_PETSc) + call VecMin(solution_vec,devNull,T_min,err_PETSc) + CHKERRQ(err_PETSc) + call VecMax(solution_vec,devNull,T_max,err_PETSc) + CHKERRQ(err_PETSc) if (solution%converged) & print'(/,1x,a)', '... thermal conduction converged ..................................' print'(/,1x,a,f8.4,2x,f8.4,2x,f8.4)', 'Minimum|Maximum|Delta Temperature / K = ', T_min, T_max, stagNorm