From 48a38f3cf516ed5596dfb7052d461cb87b4f0948 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 9 Jan 2024 14:43:13 +0100 Subject: [PATCH] harmonize status reporting --- src/grid/grid_mech_FEM.f90 | 14 +++++++------- src/grid/grid_mech_spectral_basic.f90 | 13 +++++++------ src/grid/grid_mech_spectral_polarization.f90 | 15 ++++++++------- src/grid/grid_mech_utilities.f90 | 6 ++---- 4 files changed, 24 insertions(+), 24 deletions(-) diff --git a/src/grid/grid_mech_FEM.f90 b/src/grid/grid_mech_FEM.f90 index 140b22817..258a1b511 100644 --- a/src/grid/grid_mech_FEM.f90 +++ b/src/grid/grid_mech_FEM.f90 @@ -28,7 +28,7 @@ module grid_mechanical_FEM use homogenization use discretization use discretization_grid - + use constants #if (PETSC_VERSION_MAJOR==3 && PETSC_VERSION_MINOR>14) && !defined(PETSC_HAVE_MPI_F90MODULE_VISIBILITY) implicit none(type,external) @@ -84,7 +84,7 @@ module grid_mechanical_FEM err_BC !< deviation from stress BC integer :: totalIter = 0 !< total iteration in current increment - logical :: broken + integer(kind(STATUS_OK)) :: status public :: & grid_mechanical_FEM_init, & @@ -271,7 +271,7 @@ subroutine grid_mechanical_FEM_init(num_grid) end if restartRead call utilities_updateCoords(F) - call utilities_constitutiveResponse(broken,P_current,P_av,C_volAvg,devNull, & ! stress field, stress avg, global average of stiffness and (min+max)/2 + call utilities_constitutiveResponse(status,P_current,P_av,C_volAvg,devNull, & ! stress field, stress avg, global average of stiffness and (min+max)/2 F, & ! target F 0.0_pREAL) ! time increment call DMDAVecRestoreArrayF90(DM_mech,u_PETSc,u,err_PETSc) @@ -325,7 +325,7 @@ function grid_mechanical_FEM_solution(incInfoIn) result(solution) solution%converged = reason > 0 solution%iterationsNeeded = totalIter - solution%termIll = broken + solution%termIll = STATUS_OK /= status P_aim = merge(P_av,P_aim,params%stress_mask) end function grid_mechanical_FEM_solution @@ -492,7 +492,7 @@ subroutine converged(snes_local,PETScIter,devNull1,devNull2,fnorm,reason,dummy,e BCTol = max(maxval(abs(P_av))*num%eps_stress_rtol, num%eps_stress_atol) if ((totalIter >= num%itmin .and. all([err_div/divTol, err_BC/BCTol] < 1.0_pREAL)) & - .or. broken) then + .or. STATUS_OK /= status) then reason = 1 elseif (totalIter >= num%itmax) then reason = -1 @@ -570,10 +570,10 @@ subroutine formResidual(da_local,x_local, & !-------------------------------------------------------------------------------------------------- ! evaluate constitutive response - call utilities_constitutiveResponse(broken,P_current,& + call utilities_constitutiveResponse(status,P_current,& P_av,C_volAvg,devNull, & F,params%Delta_t,params%rotation_BC) - call MPI_Allreduce(MPI_IN_PLACE,broken,1_MPI_INTEGER_KIND,MPI_LOGICAL,MPI_LOR,MPI_COMM_WORLD,err_MPI) + call MPI_Allreduce(MPI_IN_PLACE,status,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_SUM,MPI_COMM_WORLD,err_MPI) call parallelization_chkerr(err_MPI) !-------------------------------------------------------------------------------------------------- diff --git a/src/grid/grid_mech_spectral_basic.f90 b/src/grid/grid_mech_spectral_basic.f90 index cc92d4e13..2653be0d6 100644 --- a/src/grid/grid_mech_spectral_basic.f90 +++ b/src/grid/grid_mech_spectral_basic.f90 @@ -26,6 +26,7 @@ module grid_mechanical_spectral_basic use grid_mech_utilities use homogenization use discretization_grid + use constants #if (PETSC_VERSION_MAJOR==3 && PETSC_VERSION_MINOR>14) && !defined(PETSC_HAVE_MPI_F90MODULE_VISIBILITY) implicit none(type,external) @@ -84,7 +85,7 @@ module grid_mechanical_spectral_basic err_div !< RMS of div of P integer :: totalIter = 0 !< total iteration in current increment - logical :: broken + integer(kind(STATUS_OK)) :: status public :: & grid_mechanical_spectral_basic_init, & @@ -228,7 +229,7 @@ subroutine grid_mechanical_spectral_basic_init(num_grid) end if restartRead call utilities_updateCoords(reshape(F,shape(F_lastInc))) - call utilities_constitutiveResponse(broken,P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 + call utilities_constitutiveResponse(status,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(DM_mech,F_PETSc,F,err_PETSc) ! deassociate pointer @@ -287,7 +288,7 @@ function grid_mechanical_spectral_basic_solution(incInfoIn) result(solution) solution%converged = reason > 0 solution%iterationsNeeded = totalIter - solution%termIll = broken + solution%termIll = STATUS_OK /= status P_aim = merge(P_av,P_aim,params%stress_mask) end function grid_mechanical_spectral_basic_solution @@ -452,7 +453,7 @@ subroutine converged(snes_local,PETScIter,devNull1,devNull2,devNull3,reason,dumm BCTol = max(maxval(abs(P_av))*num%eps_stress_rtol, num%eps_stress_atol) if ((totalIter >= num%itmin .and. all([err_div/divTol, err_BC/BCTol] < 1.0_pREAL)) & - .or. broken) then + .or. STATUS_OK /= status) then reason = 1 elseif (totalIter >= num%itmax) then reason = -1 @@ -514,10 +515,10 @@ subroutine formResidual(residual_subdomain, F, & end if newIteration associate (P => r) - call utilities_constitutiveResponse(broken,P, & + call utilities_constitutiveResponse(status,P, & P_av,C_volAvg,C_minMaxAvg, & F,params%Delta_t,params%rotation_BC) - call MPI_Allreduce(MPI_IN_PLACE,broken,1_MPI_INTEGER_KIND,MPI_LOGICAL,MPI_LOR,MPI_COMM_WORLD,err_MPI) + call MPI_Allreduce(MPI_IN_PLACE,status,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_SUM,MPI_COMM_WORLD,err_MPI) call parallelization_chkerr(err_MPI) err_div = utilities_divergenceRMS(P) end associate diff --git a/src/grid/grid_mech_spectral_polarization.f90 b/src/grid/grid_mech_spectral_polarization.f90 index 4328e591f..93e3cfc4a 100644 --- a/src/grid/grid_mech_spectral_polarization.f90 +++ b/src/grid/grid_mech_spectral_polarization.f90 @@ -27,6 +27,7 @@ module grid_mechanical_spectral_polarization use config use homogenization use discretization_grid + use constants #if (PETSC_VERSION_MAJOR==3 && PETSC_VERSION_MINOR>14) && !defined(PETSC_HAVE_MPI_F90MODULE_VISIBILITY) implicit none(type,external) @@ -96,7 +97,7 @@ module grid_mechanical_spectral_polarization err_div !< RMS of div of P integer :: totalIter = 0 !< total iteration in current increment - logical :: broken + integer(kind(STATUS_OK)) :: status public :: & grid_mechanical_spectral_polarization_init, & @@ -257,7 +258,7 @@ subroutine grid_mechanical_spectral_polarization_init(num_grid) end if restartRead call utilities_updateCoords(reshape(F,shape(F_lastInc))) - call utilities_constitutiveResponse(broken,P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 + call utilities_constitutiveResponse(status,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(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc) ! deassociate pointer @@ -322,7 +323,7 @@ function grid_mechanical_spectral_polarization_solution(incInfoIn) result(soluti solution%converged = reason > 0 solution%iterationsNeeded = totalIter - solution%termIll = broken + solution%termIll = STATUS_OK /= status P_aim = merge(P_av,P_aim,params%stress_mask) end function grid_mechanical_spectral_polarization_solution @@ -516,7 +517,7 @@ subroutine converged(snes_local,PETScIter,devNull1,devNull2,devNull3,reason,dumm BCTol = max(maxval(abs(P_av))*num%eps_stress_rtol, num%eps_stress_atol) if ((totalIter >= num%itmin .and. all([err_div/divTol, err_curl/curlTol, err_BC/BCTol] < 1.0_pREAL)) & - .or. broken) then + .or. STATUS_OK /= status) then reason = 1 elseif (totalIter >= num%itmax) then reason = -1 @@ -604,14 +605,14 @@ subroutine formResidual(residual_subdomain, FandF_tau, & err_curl = utilities_curlRMS(F) #ifdef __GFORTRAN__ - call utilities_constitutiveResponse(broken,r_F, & + call utilities_constitutiveResponse(status,r_F, & #else associate (P => r_F) - call utilities_constitutiveResponse(broken, P, & + call utilities_constitutiveResponse(status, P, & #endif P_av,C_volAvg,C_minMaxAvg, & F - r_F_tau/num%beta,params%Delta_t,params%rotation_BC) - call MPI_Allreduce(MPI_IN_PLACE,broken,1_MPI_INTEGER_KIND,MPI_LOGICAL,MPI_LOR,MPI_COMM_WORLD,err_MPI) + call MPI_Allreduce(MPI_IN_PLACE,status,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_SUM,MPI_COMM_WORLD,err_MPI) #ifdef __GFORTRAN__ err_div = utilities_divergenceRMS(r_F) #else diff --git a/src/grid/grid_mech_utilities.f90 b/src/grid/grid_mech_utilities.f90 index d05d4a51f..a4080efc2 100644 --- a/src/grid/grid_mech_utilities.f90 +++ b/src/grid/grid_mech_utilities.f90 @@ -113,10 +113,10 @@ end function utilities_maskedCompliance !-------------------------------------------------------------------------------------------------- !> @brief Calculate constitutive response. !-------------------------------------------------------------------------------------------------- -subroutine utilities_constitutiveResponse(broken, P,P_av,C_volAvg,C_minmaxAvg,& +subroutine utilities_constitutiveResponse(status, P,P_av,C_volAvg,C_minmaxAvg,& F,Delta_t,rotation_BC) - logical, intent(out) :: broken + integer(kind(STATUS_OK)), intent(out) :: status real(pREAL), intent(out), dimension(3,3,3,3) :: C_volAvg, C_minmaxAvg !< average stiffness real(pREAL), intent(out), dimension(3,3) :: P_av !< average PK stress real(pREAL), intent(out), dimension(3,3,cells(1),cells(2),cells3) :: P !< PK stress @@ -129,7 +129,6 @@ subroutine utilities_constitutiveResponse(broken, P,P_av,C_volAvg,C_minmaxAvg,& real(pREAL), dimension(3,3,3,3) :: dPdF_max, dPdF_min real(pREAL) :: dPdF_norm_max, dPdF_norm_min real(pREAL), dimension(2) :: valueAndRank !< pair of min/max norm of dPdF to synchronize min/max of dPdF - integer(kind(STATUS_OK)) :: status print'(/,1x,a)', '... evaluating constitutive response ......................................' @@ -138,7 +137,6 @@ subroutine utilities_constitutiveResponse(broken, P,P_av,C_volAvg,C_minmaxAvg,& homogenization_F = reshape(F,[3,3,product(cells(1:2))*cells3]) ! set materialpoint target F to estimated field call homogenization_mechanical_response(status,Delta_t,1,product(cells(1:2))*cells3) ! calculate P field - broken = STATUS_OK /= status P = reshape(homogenization_P, [3,3,cells(1),cells(2),cells3]) P_av = sum(sum(sum(P,dim=5),dim=4),dim=3) * wgt