following structure of mech_spectral_basic
This commit is contained in:
parent
7ccb5c1fa3
commit
a06be13d49
|
@ -64,6 +64,7 @@ module grid_mech_spectral_basic
|
|||
grid_mech_spectral_basic_solution, &
|
||||
grid_mech_spectral_basic_forward
|
||||
private :: &
|
||||
converged, &
|
||||
formResidual
|
||||
|
||||
contains
|
||||
|
@ -151,7 +152,7 @@ subroutine grid_mech_spectral_basic_init
|
|||
call DMcreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr) ! global solution vector (grid x 9, i.e. every def grad tensor)
|
||||
call DMDASNESsetFunctionLocal(da,INSERT_VALUES,formResidual,PETSC_NULL_SNES,ierr) ! residual vector of same shape as solution vector
|
||||
CHKERRQ(ierr)
|
||||
call SNESsetConvergenceTest(snes,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,ierr)! specify custom convergence check function "_converged"
|
||||
call SNESsetConvergenceTest(snes,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,ierr)! specify custom convergence check function "converged"
|
||||
CHKERRQ(ierr)
|
||||
call SNESsetFromOptions(snes,ierr); CHKERRQ(ierr) ! pull it all together with additional CLI arguments
|
||||
|
||||
|
@ -267,7 +268,6 @@ function grid_mech_spectral_basic_solution(incInfoIn,timeinc,timeinc_old,stress_
|
|||
solution%termIll = terminallyIll
|
||||
terminallyIll = .false.
|
||||
|
||||
|
||||
end function grid_mech_spectral_basic_solution
|
||||
|
||||
|
||||
|
|
|
@ -20,9 +20,6 @@ module grid_mech_spectral_polarisation
|
|||
implicit none
|
||||
private
|
||||
|
||||
character (len=*), parameter, public :: &
|
||||
GRID_MECH_SPECTRAL_POLARISATION_LABEL = 'polarisation'
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! derived types
|
||||
type(tSolutionParams), private :: params
|
||||
|
@ -74,6 +71,7 @@ module grid_mech_spectral_polarisation
|
|||
grid_mech_spectral_polarisation_solution, &
|
||||
grid_mech_spectral_polarisation_forward
|
||||
private :: &
|
||||
converged, &
|
||||
formResidual
|
||||
|
||||
contains
|
||||
|
@ -86,10 +84,6 @@ subroutine grid_mech_spectral_polarisation_init
|
|||
IO_intOut, &
|
||||
IO_error, &
|
||||
IO_open_jobFile_binary
|
||||
use debug, only: &
|
||||
debug_level, &
|
||||
debug_spectral, &
|
||||
debug_spectralRestart
|
||||
use FEsolving, only: &
|
||||
restartInc
|
||||
use numerics, only: &
|
||||
|
@ -166,7 +160,7 @@ subroutine grid_mech_spectral_polarisation_init
|
|||
call DMcreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr) ! global solution vector (grid x 18, i.e. every def grad tensor)
|
||||
call DMDASNESsetFunctionLocal(da,INSERT_VALUES,formResidual,PETSC_NULL_SNES,ierr) ! residual vector of same shape as solution vector
|
||||
CHKERRQ(ierr)
|
||||
call SNESsetConvergenceTest(snes,grid_mech_spectral_polarisation_converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,ierr) ! specify custom convergence check function "_converged"
|
||||
call SNESsetConvergenceTest(snes,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,ierr) ! specify custom convergence check function "converged"
|
||||
CHKERRQ(ierr)
|
||||
call SNESsetFromOptions(snes,ierr); CHKERRQ(ierr) ! pull it all together with additional CLI arguments
|
||||
|
||||
|
@ -175,12 +169,9 @@ subroutine grid_mech_spectral_polarisation_init
|
|||
call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr) ! places pointer on PETSc data
|
||||
F => FandF_tau( 0: 8,:,:,:)
|
||||
F_tau => FandF_tau( 9:17,:,:,:)
|
||||
|
||||
restart: if (restartInc > 0) then
|
||||
if (iand(debug_level(debug_spectral),debug_spectralRestart) /= 0) then
|
||||
write(6,'(/,a,'//IO_intOut(restartInc)//',a)') &
|
||||
'reading values of increment ', restartInc, ' from file'
|
||||
flush(6)
|
||||
endif
|
||||
write(6,'(/,a,'//IO_intOut(restartInc)//',a)') 'reading values of increment ', restartInc, ' from file'
|
||||
|
||||
fileUnit = IO_open_jobFile_binary('F_aimDot')
|
||||
read(fileUnit) F_aimDot; close(fileUnit)
|
||||
|
@ -218,10 +209,8 @@ subroutine grid_mech_spectral_polarisation_init
|
|||
call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr) ! deassociate pointer
|
||||
|
||||
restartRead: if (restartInc > 0) then
|
||||
if (iand(debug_level(debug_spectral),debug_spectralRestart)/= 0 .and. worldrank == 0) &
|
||||
write(6,'(/,a,'//IO_intOut(restartInc)//',a)') &
|
||||
'reading more values of increment ', restartInc, ' from file'
|
||||
flush(6)
|
||||
write(6,'(/,a,'//IO_intOut(restartInc)//',a)') 'reading more values of increment ', restartInc, ' from file'
|
||||
|
||||
fileUnit = IO_open_jobFile_binary('C_volAvg')
|
||||
read(fileUnit) C_volAvg; close(fileUnit)
|
||||
fileUnit = IO_open_jobFile_binary('C_volAvgLastInv')
|
||||
|
@ -309,77 +298,6 @@ function grid_mech_spectral_polarisation_solution(incInfoIn,timeinc,timeinc_old,
|
|||
end function grid_mech_spectral_polarisation_solution
|
||||
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
!> @brief convergence check
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
subroutine grid_mech_spectral_polarisation_converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dummy,ierr)
|
||||
use numerics, only: &
|
||||
itmax, &
|
||||
itmin, &
|
||||
err_div_tolRel, &
|
||||
err_div_tolAbs, &
|
||||
err_curl_tolRel, &
|
||||
err_curl_tolAbs, &
|
||||
err_stress_tolRel, &
|
||||
err_stress_tolAbs
|
||||
use math, only: &
|
||||
math_mul3333xx33
|
||||
use FEsolving, only: &
|
||||
terminallyIll
|
||||
|
||||
implicit none
|
||||
SNES :: snes_local
|
||||
PetscInt :: PETScIter
|
||||
PetscReal :: &
|
||||
xnorm, & ! not used
|
||||
snorm, & ! not used
|
||||
fnorm ! not used
|
||||
SNESConvergedReason :: reason
|
||||
PetscObject :: dummy
|
||||
PetscErrorCode :: ierr
|
||||
real(pReal) :: &
|
||||
curlTol, &
|
||||
divTol, &
|
||||
BCTol
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! stress BC handling
|
||||
F_aim = F_aim - math_mul3333xx33(S, ((P_av - params%stress_BC))) ! S = 0.0 for no bc
|
||||
err_BC = maxval(abs((1.0_pReal-params%stress_mask) * math_mul3333xx33(C_scale,F_aim-F_av) + &
|
||||
params%stress_mask * (P_av-params%stress_BC))) ! mask = 0.0 for no bc
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! error calculation
|
||||
curlTol = max(maxval(abs(F_aim-math_I3))*err_curl_tolRel ,err_curl_tolAbs)
|
||||
divTol = max(maxval(abs(P_av)) *err_div_tolRel ,err_div_tolAbs)
|
||||
BCTol = max(maxval(abs(P_av)) *err_stress_tolRel,err_stress_tolAbs)
|
||||
|
||||
converged: if ((totalIter >= itmin .and. &
|
||||
all([ err_div /divTol, &
|
||||
err_curl/curlTol, &
|
||||
err_BC /BCTol ] < 1.0_pReal)) &
|
||||
.or. terminallyIll) then
|
||||
reason = 1
|
||||
elseif (totalIter >= itmax) then converged
|
||||
reason = -1
|
||||
else converged
|
||||
reason = 0
|
||||
endif converged
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! report
|
||||
write(6,'(1/,a)') ' ... reporting .............................................................'
|
||||
write(6,'(/,a,f12.2,a,es8.2,a,es9.2,a)') ' error divergence = ', &
|
||||
err_div/divTol, ' (',err_div, ' / m, tol = ',divTol,')'
|
||||
write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error curl = ', &
|
||||
err_curl/curlTol,' (',err_curl,' -, tol = ',curlTol,')'
|
||||
write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error BC = ', &
|
||||
err_BC/BCTol, ' (',err_BC, ' Pa, tol = ',BCTol,')'
|
||||
write(6,'(/,a)') ' ==========================================================================='
|
||||
flush(6)
|
||||
|
||||
end subroutine grid_mech_spectral_polarisation_converged
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
!> @brief forwarding routine
|
||||
!> @details find new boundary conditions and best F estimate for end of current timestep
|
||||
|
@ -467,7 +385,7 @@ subroutine grid_mech_spectral_polarisation_forward(guess,timeinc,timeinc_old,loa
|
|||
write(fileUnit) F_tau_lastInc; close (fileUnit)
|
||||
endif
|
||||
|
||||
call CPFEM_age() ! age state and kinematics
|
||||
call CPFEM_age ! age state and kinematics
|
||||
call utilities_updateIPcoords(F)
|
||||
|
||||
C_volAvgLastInc = C_volAvg
|
||||
|
@ -504,7 +422,6 @@ subroutine grid_mech_spectral_polarisation_forward(guess,timeinc,timeinc_old,loa
|
|||
!--------------------------------------------------------------------------------------------------
|
||||
! update average and local deformation gradients
|
||||
F_aim = F_aim_lastInc + F_aimDot * timeinc
|
||||
|
||||
F = reshape(Utilities_forwardField(timeinc,F_lastInc,Fdot, & ! estimate of F at end of time+timeinc that matches rotated F_aim on average
|
||||
math_rotate_backward33(F_aim,rotation_BC)),&
|
||||
[9,grid(1),grid(2),grid3])
|
||||
|
@ -530,6 +447,76 @@ subroutine grid_mech_spectral_polarisation_forward(guess,timeinc,timeinc_old,loa
|
|||
end subroutine grid_mech_spectral_polarisation_forward
|
||||
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
!> @brief convergence check
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
subroutine converged(snes_local,PETScIter,xnorm,snorm,fnorm,reason,dummy,ierr)
|
||||
use numerics, only: &
|
||||
itmax, &
|
||||
itmin, &
|
||||
err_div_tolRel, &
|
||||
err_div_tolAbs, &
|
||||
err_curl_tolRel, &
|
||||
err_curl_tolAbs, &
|
||||
err_stress_tolRel, &
|
||||
err_stress_tolAbs
|
||||
use math, only: &
|
||||
math_mul3333xx33
|
||||
use FEsolving, only: &
|
||||
terminallyIll
|
||||
|
||||
implicit none
|
||||
SNES :: snes_local
|
||||
PetscInt :: PETScIter
|
||||
PetscReal :: &
|
||||
xnorm, & ! not used
|
||||
snorm, & ! not used
|
||||
fnorm ! not used
|
||||
SNESConvergedReason :: reason
|
||||
PetscObject :: dummy
|
||||
PetscErrorCode :: ierr
|
||||
real(pReal) :: &
|
||||
curlTol, &
|
||||
divTol, &
|
||||
BCTol
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! stress BC handling
|
||||
F_aim = F_aim - math_mul3333xx33(S, ((P_av - params%stress_BC))) ! S = 0.0 for no bc
|
||||
err_BC = maxval(abs((1.0_pReal-params%stress_mask) * math_mul3333xx33(C_scale,F_aim-F_av) + &
|
||||
params%stress_mask * (P_av-params%stress_BC))) ! mask = 0.0 for no bc
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! error calculation
|
||||
curlTol = max(maxval(abs(F_aim-math_I3))*err_curl_tolRel ,err_curl_tolAbs)
|
||||
divTol = max(maxval(abs(P_av)) *err_div_tolRel ,err_div_tolAbs)
|
||||
BCTol = max(maxval(abs(P_av)) *err_stress_tolRel,err_stress_tolAbs)
|
||||
|
||||
if ((totalIter >= itmin .and. &
|
||||
all([ err_div /divTol, &
|
||||
err_curl/curlTol, &
|
||||
err_BC /BCTol ] < 1.0_pReal)) &
|
||||
.or. terminallyIll) then
|
||||
reason = 1
|
||||
elseif (totalIter >= itmax) then
|
||||
reason = -1
|
||||
else
|
||||
reason = 0
|
||||
endif
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! report
|
||||
write(6,'(1/,a)') ' ... reporting .............................................................'
|
||||
write(6,'(/,a,f12.2,a,es8.2,a,es9.2,a)') ' error divergence = ', &
|
||||
err_div/divTol, ' (',err_div, ' / m, tol = ',divTol,')'
|
||||
write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error curl = ', &
|
||||
err_curl/curlTol,' (',err_curl,' -, tol = ',curlTol,')'
|
||||
write(6, '(a,f12.2,a,es8.2,a,es9.2,a)') ' error BC = ', &
|
||||
err_BC/BCTol, ' (',err_BC, ' Pa, tol = ',BCTol,')'
|
||||
write(6,'(/,a)') ' ==========================================================================='
|
||||
flush(6)
|
||||
|
||||
end subroutine converged
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
!> @brief forms the polarisation residual vector
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
|
@ -570,9 +557,9 @@ subroutine formResidual(in, FandF_tau, &
|
|||
|
||||
implicit none
|
||||
DMDALocalInfo, dimension(DMDA_LOCAL_INFO_SIZE) :: in !< DMDA info (needs to be named "in" for macros like XRANGE to work)
|
||||
PetscScalar, dimension(3,3,2, XG_RANGE,YG_RANGE,ZG_RANGE), &
|
||||
PetscScalar, dimension(3,3,2,XG_RANGE,YG_RANGE,ZG_RANGE), &
|
||||
target, intent(in) :: FandF_tau
|
||||
PetscScalar, dimension(3,3,2, X_RANGE, Y_RANGE, Z_RANGE),&
|
||||
PetscScalar, dimension(3,3,2,X_RANGE,Y_RANGE,Z_RANGE),&
|
||||
target, intent(out) :: residuum !< residuum field
|
||||
PetscScalar, pointer, dimension(:,:,:,:,:) :: &
|
||||
F, &
|
||||
|
@ -668,7 +655,7 @@ subroutine formResidual(in, FandF_tau, &
|
|||
! calculating curl
|
||||
tensorField_real = 0.0_pReal
|
||||
tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3) = F
|
||||
call utilities_FFTtensorForward()
|
||||
call utilities_FFTtensorForward
|
||||
err_curl = Utilities_curlRMS()
|
||||
|
||||
end subroutine formResidual
|
||||
|
|
Loading…
Reference in New Issue