From adb68ff7927d6dc8a43604ab0e96fcab15a823b9 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 14 Jan 2019 12:45:07 +0100 Subject: [PATCH] avoid converstions and use of global variables --- PRIVATE | 2 +- src/CPFEM.f90 | 8 ++++---- src/DAMASK_abaqus.f | 8 +++----- src/DAMASK_marc.f90 | 11 ++++------- src/math.f90 | 12 ++++++------ 5 files changed, 18 insertions(+), 23 deletions(-) diff --git a/PRIVATE b/PRIVATE index 5ed6a1f60..6e7550042 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 5ed6a1f60b412eb46ff6820cf03b684095ff1f75 +Subproject commit 6e7550042259f46992329d202f5804df48ce99ff diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index 515b28ab0..11bfd8be7 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -594,7 +594,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt ! translate from P to CS Kirchhoff = math_mul33x33(materialpoint_P(1:3,1:3,ip,elCP), transpose(materialpoint_F(1:3,1:3,ip,elCP))) J_inverse = 1.0_pReal / math_det33(materialpoint_F(1:3,1:3,ip,elCP)) - CPFEM_cs(1:6,ip,elCP) = math_33to6(J_inverse * Kirchhoff) + CPFEM_cs(1:6,ip,elCP) = math_33to6(J_inverse * Kirchhoff,weighted=.false.) ! translate from dP/dF to dCS/dE H = 0.0_pReal @@ -610,7 +610,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt forall(i=1:3, j=1:3,k=1:3,l=1:3) & H_sym(i,j,k,l) = 0.25_pReal * (H(i,j,k,l) + H(j,i,k,l) + H(i,j,l,k) + H(j,i,l,k)) - CPFEM_dcsde(1:6,1:6,ip,elCP) = math_3333to66(J_inverse * H_sym) + CPFEM_dcsde(1:6,1:6,ip,elCP) = math_3333to66(J_inverse * H_sym,weighted=.false.) endif terminalIllness endif validCalculation @@ -637,7 +637,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt !*** remember extreme values of stress ... - cauchyStress33 = math_6to33(CPFEM_cs(1:6,ip,elCP)) + cauchyStress33 = math_6to33(CPFEM_cs(1:6,ip,elCP),weighted=.false.) if (maxval(cauchyStress33) > debug_stressMax) then debug_stressMaxLocation = [elCP, ip] debug_stressMax = maxval(cauchyStress33) @@ -647,7 +647,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt debug_stressMin = minval(cauchyStress33) endif !*** ... and Jacobian - jacobian3333 = math_66to3333(CPFEM_dcsdE(1:6,1:6,ip,elCP)) + jacobian3333 = math_66to3333(CPFEM_dcsdE(1:6,1:6,ip,elCP),weighted=.false.) if (maxval(jacobian3333) > debug_jacobianMax) then debug_jacobianMaxLocation = [elCP, ip] debug_jacobianMax = maxval(jacobian3333) diff --git a/src/DAMASK_abaqus.f b/src/DAMASK_abaqus.f index 6c6434e4a..9072de95d 100644 --- a/src/DAMASK_abaqus.f +++ b/src/DAMASK_abaqus.f @@ -102,8 +102,6 @@ subroutine UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,& calcMode, & terminallyIll, & symmetricSolver - use math, only: & - invnrmMandel use debug, only: & debug_info, & debug_reset, & @@ -305,9 +303,9 @@ subroutine UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,& ! ABAQUS implicit: 11, 22, 33, 12, 13, 23 ! ABAQUS implicit: 11, 22, 33, 12 - forall(i=1:ntens) ddsdde(1:ntens,i) = invnrmMandel(i)*ddsdde_h(1:ntens,i)*invnrmMandel(1:ntens) - stress(1:ntens) = stress_h(1:ntens)*invnrmMandel(1:ntens) - if(symmetricSolver) ddsdde(1:ntens,1:ntens) = 0.5_pReal*(ddsdde(1:ntens,1:ntens) + transpose(ddsdde(1:ntens,1:ntens))) + ddsdde = ddsdde_h(1:ntens,1:ntens) + stress = stress_h(1:ntens) + if(symmetricSolver) ddsdde = 0.5_pReal*(ddsdde + transpose(ddsdde)) if(ntens == 6) then stress_h = stress stress(5) = stress_h(6) diff --git a/src/DAMASK_marc.f90 b/src/DAMASK_marc.f90 index f3130c5cd..0c7d1adeb 100644 --- a/src/DAMASK_marc.f90 +++ b/src/DAMASK_marc.f90 @@ -127,9 +127,6 @@ subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,nshear,disp, & calcMode, & terminallyIll, & symmetricSolver - use math, only: & - math_transpose33,& - invnrmMandel use debug, only: & debug_level, & debug_LEVELBASIC, & @@ -235,9 +232,9 @@ subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,nshear,disp, & write(6,'(a,i12)') ' Nodes: ', nnode write(6,'(a,i1)') ' Deformation gradient: ', itel write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') ' Deformation gradient at t=n:', & - math_transpose33(ffn) + transpose(ffn) write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') ' Deformation gradient at t=n+1:', & - math_transpose33(ffn1) + transpose(ffn1) endif !$ defaultNumThreadsInt = omp_get_num_threads() ! remember number of threads set by Marc @@ -357,8 +354,8 @@ subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,nshear,disp, & ! Marc: 11, 22, 33, 12, 23, 13 ! Marc: 11, 22, 33, 12 - forall(i=1:ngens) d(1:ngens,i) = invnrmMandel(i)*ddsdde(1:ngens,i)*invnrmMandel(1:ngens) - s(1:ndi+nshear) = stress(1:ndi+nshear)*invnrmMandel(1:ndi+nshear) + d = ddsdde(1:ngens,1:ngens) + s = stress(1:ndi+nshear) g = 0.0_pReal if(symmetricSolver) d = 0.5_pReal*(d+transpose(d)) diff --git a/src/math.f90 b/src/math.f90 index 3d60eba0c..1b782dce2 100644 --- a/src/math.f90 +++ b/src/math.f90 @@ -29,7 +29,7 @@ module math 1.0_pReal, 1.0_pReal, 1.0_pReal, & sqrt(2.0_pReal), sqrt(2.0_pReal), sqrt(2.0_pReal) ] !< weighting for Mandel notation (forward) - real(pReal), dimension(6), parameter , public :: & + real(pReal), dimension(6), parameter , private :: & invnrmMandel = [& 1.0_pReal, 1.0_pReal, 1.0_pReal, & 1.0_pReal/sqrt(2.0_pReal), 1.0_pReal/sqrt(2.0_pReal), 1.0_pReal/sqrt(2.0_pReal) ] !< weighting for Mandel notation (backward) @@ -1180,7 +1180,7 @@ end function math_6to33 !-------------------------------------------------------------------------------------------------- -!> @brief convert 3333 matrix into vector 99 +!> @brief convert 3333 matrix into 99 matrix !-------------------------------------------------------------------------------------------------- pure function math_3333to99(m3333) @@ -1197,7 +1197,7 @@ end function math_3333to99 !-------------------------------------------------------------------------------------------------- -!> @brief convert 99 vector into 3333 matrix +!> @brief convert 99 matrix into 3333 matrix !-------------------------------------------------------------------------------------------------- pure function math_99to3333(m99) @@ -1214,7 +1214,7 @@ end function math_99to3333 !-------------------------------------------------------------------------------------------------- -!> @brief convert symmetric 3333 matrix into 66 vector +!> @brief convert symmetric 3333 matrix into 66 matrix !> @details Weighted conversion (default) rearranges according to Nye and weights shear ! components according to Mandel. Advisable for matrix operations. ! Unweighted conversion only changes order according to Nye @@ -1242,7 +1242,7 @@ end function math_3333to66 !-------------------------------------------------------------------------------------------------- -!> @brief convert 66 vector into symmetric 3333 matrix +!> @brief convert 66 matrix into symmetric 3333 matrix !> @details Weighted conversion (default) rearranges according to Nye and weights shear ! components according to Mandel. Advisable for matrix operations. ! Unweighted conversion only changes order according to Nye @@ -1274,7 +1274,7 @@ end function math_66to3333 !-------------------------------------------------------------------------------------------------- -!> @brief convert Voigt matrix 66 back to symmetric 3333 tensor +!> @brief convert 66 Voigt matrix into symmetric 3333 matrix !-------------------------------------------------------------------------------------------------- pure function math_Voigt66to3333(m66)