use Fortran internals for error handling

This commit is contained in:
Martin Diehl 2020-11-11 15:57:05 +01:00
parent 3dd5eaf1c1
commit cd7ada0da9
5 changed files with 8 additions and 13 deletions

View File

@ -463,10 +463,6 @@ subroutine IO_error(error_ID,el,ip,g,instance,ext_msg)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! math errors ! math errors
case (400)
msg = 'matrix inversion error'
case (401)
msg = 'error in Eigenvalue calculation'
case (402) case (402)
msg = 'invalid orientation specified' msg = 'invalid orientation specified'

View File

@ -713,7 +713,6 @@ function utilities_maskedCompliance(rot_BC,mask_stress,C)
allocate(s_reduced,mold = c_reduced) allocate(s_reduced,mold = c_reduced)
call math_invert(s_reduced, errmatinv, c_reduced) ! invert reduced stiffness call math_invert(s_reduced, errmatinv, c_reduced) ! invert reduced stiffness
if (any(IEEE_is_NaN(s_reduced))) errmatinv = .true. if (any(IEEE_is_NaN(s_reduced))) errmatinv = .true.
if (errmatinv) call IO_error(error_ID=400,ext_msg='utilities_maskedCompliance')
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! check if inversion was successful ! check if inversion was successful
@ -725,7 +724,7 @@ function utilities_maskedCompliance(rot_BC,mask_stress,C)
write(IO_STDOUT,trim(formatString),advance='no') ' C * S (load) ', & write(IO_STDOUT,trim(formatString),advance='no') ' C * S (load) ', &
transpose(matmul(c_reduced,s_reduced)) transpose(matmul(c_reduced,s_reduced))
write(IO_STDOUT,trim(formatString),advance='no') ' S (load) ', transpose(s_reduced) write(IO_STDOUT,trim(formatString),advance='no') ' S (load) ', transpose(s_reduced)
if(errmatinv) call IO_error(error_ID=400,ext_msg='utilities_maskedCompliance') if(errmatinv) error stop 'matrix inversion error'
endif endif
temp99_real = reshape(unpack(reshape(s_reduced,[size_reduced**2]),reshape(mask,[81]),0.0_pReal),[9,9]) temp99_real = reshape(unpack(reshape(s_reduced,[size_reduced**2]),reshape(mask,[81]),0.0_pReal),[9,9])
else else

View File

@ -2209,11 +2209,11 @@ function equivalent_nu(C,assumption) result(nu)
/ 9.0_pReal / 9.0_pReal
elseif(IO_lc(assumption) == 'reuss') then elseif(IO_lc(assumption) == 'reuss') then
call math_invert(S,error,C) call math_invert(S,error,C)
if(error) call IO_error(0) if(error) error stop 'matrix inversion failed'
K = 1.0_pReal & K = 1.0_pReal &
/ (S(1,1)+S(2,2)+S(3,3) +2.0_pReal*(S(1,2)+S(2,3)+S(1,3))) / (S(1,1)+S(2,2)+S(3,3) +2.0_pReal*(S(1,2)+S(2,3)+S(1,3)))
else else
call IO_error(0) error stop 'invalid assumption'
K = 0.0_pReal K = 0.0_pReal
endif endif
@ -2241,11 +2241,11 @@ function equivalent_mu(C,assumption) result(mu)
/ 15.0_pReal / 15.0_pReal
elseif(IO_lc(assumption) == 'reuss') then elseif(IO_lc(assumption) == 'reuss') then
call math_invert(S,error,C) call math_invert(S,error,C)
if(error) call IO_error(0) if(error) error stop 'matrix inversion failed'
mu = 15.0_pReal & mu = 15.0_pReal &
/ (4.0_pReal*(S(1,1)+S(2,2)+S(3,3)) -4.0_pReal*(S(1,2)+S(2,3)+S(1,3)) +3.0_pReal*(S(4,4)+S(5,5)+S(6,6))) / (4.0_pReal*(S(1,1)+S(2,2)+S(3,3)) -4.0_pReal*(S(1,2)+S(2,3)+S(1,3)) +3.0_pReal*(S(4,4)+S(5,5)+S(6,6)))
else else
call IO_error(0) error stop 'invalid assumption'
mu = 0.0_pReal mu = 0.0_pReal
endif endif

View File

@ -499,7 +499,7 @@ function math_invSym3333(A)
call dgetrf(6,6,temp66,6,ipiv6,ierr_i) call dgetrf(6,6,temp66,6,ipiv6,ierr_i)
call dgetri(6,temp66,6,ipiv6,work,size(work,1),ierr_f) call dgetri(6,temp66,6,ipiv6,work,size(work,1),ierr_f)
if (ierr_i /= 0 .or. ierr_f /= 0) then if (ierr_i /= 0 .or. ierr_f /= 0) then
call IO_error(400, ext_msg = 'math_invSym3333') error stop 'matrix inversion error'
else else
math_invSym3333 = math_66toSym3333(temp66) math_invSym3333 = math_66toSym3333(temp66)
endif endif

View File

@ -640,13 +640,13 @@ function om2ax(om) result(ax)
ax(1:3) = [ 0.0_pReal, 0.0_pReal, 1.0_pReal ] ax(1:3) = [ 0.0_pReal, 0.0_pReal, 1.0_pReal ]
else else
call dgeev('N','V',3,om_,3,Wr,Wi,devNull,3,VR,3,work,size(work,1),ierr) call dgeev('N','V',3,om_,3,Wr,Wi,devNull,3,VR,3,work,size(work,1),ierr)
if (ierr /= 0) call IO_error(401,ext_msg='Error in om2ax: DGEEV return not zero') if (ierr /= 0) error stop 'LAPACK error'
#if defined(__GFORTRAN__) && __GNUC__<9 || defined(__INTEL_COMPILER) && INTEL_COMPILER<1800 || defined(__PGI) #if defined(__GFORTRAN__) && __GNUC__<9 || defined(__INTEL_COMPILER) && INTEL_COMPILER<1800 || defined(__PGI)
i = maxloc(merge(1,0,cEq(cmplx(Wr,Wi,pReal),cmplx(1.0_pReal,0.0_pReal,pReal),tol=1.0e-14_pReal)),dim=1) i = maxloc(merge(1,0,cEq(cmplx(Wr,Wi,pReal),cmplx(1.0_pReal,0.0_pReal,pReal),tol=1.0e-14_pReal)),dim=1)
#else #else
i = findloc(cEq(cmplx(Wr,Wi,pReal),cmplx(1.0_pReal,0.0_pReal,pReal),tol=1.0e-14_pReal),.true.,dim=1) !find eigenvalue (1,0) i = findloc(cEq(cmplx(Wr,Wi,pReal),cmplx(1.0_pReal,0.0_pReal,pReal),tol=1.0e-14_pReal),.true.,dim=1) !find eigenvalue (1,0)
#endif #endif
if (i == 0) call IO_error(401,ext_msg='Error in om2ax Real: eigenvalue not found') if (i == 0) error stop 'om2ax conversion failed'
ax(1:3) = VR(1:3,i) ax(1:3) = VR(1:3,i)
where ( dNeq0([om(2,3)-om(3,2), om(3,1)-om(1,3), om(1,2)-om(2,1)])) & where ( dNeq0([om(2,3)-om(3,2), om(3,1)-om(1,3), om(1,2)-om(2,1)])) &
ax(1:3) = sign(ax(1:3),-P *[om(2,3)-om(3,2), om(3,1)-om(1,3), om(1,2)-om(2,1)]) ax(1:3) = sign(ax(1:3),-P *[om(2,3)-om(3,2), om(3,1)-om(1,3), om(1,2)-om(2,1)])