diff --git a/src/DAMASK_spectral.f90 b/src/DAMASK_spectral.f90 old mode 100755 new mode 100644 diff --git a/src/crystallite.f90 b/src/crystallite.f90 old mode 100755 new mode 100644 diff --git a/src/math.f90 b/src/math.f90 index ee228834a..7e35ca390 100644 --- a/src/math.f90 +++ b/src/math.f90 @@ -2755,8 +2755,7 @@ pure function math_rotate_forward33(tensor,rot_tensor) real(pReal), dimension(3,3) :: math_rotate_forward33 real(pReal), dimension(3,3), intent(in) :: tensor, rot_tensor - math_rotate_forward33 = math_mul33x33(rot_tensor,& - math_mul33x33(tensor,math_transpose33(rot_tensor))) + math_rotate_forward33 = math_mul33x33(rot_tensor,math_mul33x33(tensor,transpose(rot_tensor))) end function math_rotate_forward33 @@ -2770,8 +2769,7 @@ pure function math_rotate_backward33(tensor,rot_tensor) real(pReal), dimension(3,3) :: math_rotate_backward33 real(pReal), dimension(3,3), intent(in) :: tensor, rot_tensor - math_rotate_backward33 = math_mul33x33(math_transpose33(rot_tensor),& - math_mul33x33(tensor,rot_tensor)) + math_rotate_backward33 = math_mul33x33(transpose(rot_tensor),math_mul33x33(tensor,rot_tensor)) end function math_rotate_backward33 @@ -2792,8 +2790,8 @@ pure function math_rotate_forward3333(tensor,rot_tensor) do i = 1_pInt,3_pInt; do j = 1_pInt,3_pInt; do k = 1_pInt,3_pInt; do l = 1_pInt,3_pInt do m = 1_pInt,3_pInt; do n = 1_pInt,3_pInt; do o = 1_pInt,3_pInt; do p = 1_pInt,3_pInt math_rotate_forward3333(i,j,k,l) = math_rotate_forward3333(i,j,k,l) & - + rot_tensor(m,i) * rot_tensor(n,j) & - * rot_tensor(o,k) * rot_tensor(p,l) * tensor(m,n,o,p) + + rot_tensor(i,m) * rot_tensor(j,n) & + * rot_tensor(k,o) * rot_tensor(l,p) * tensor(m,n,o,p) enddo; enddo; enddo; enddo; enddo; enddo; enddo; enddo end function math_rotate_forward3333 diff --git a/src/spectral_mech_Basic.f90 b/src/spectral_mech_Basic.f90 index 2a5e1838d..171eeacad 100644 --- a/src/spectral_mech_Basic.f90 +++ b/src/spectral_mech_Basic.f90 @@ -301,7 +301,6 @@ subroutine BasicPETSC_formResidual(in,x_scal,f_scal,dummy,ierr) grid3 use math, only: & math_rotate_backward33, & - math_transpose33, & math_mul3333xx33 use debug, only: & debug_level, & @@ -349,9 +348,9 @@ subroutine BasicPETSC_formResidual(in,x_scal,f_scal,dummy,ierr) trim(incInfo), ' @ Iteration ', itmin, '≤',totalIter, '≤', itmax if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) & write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') & - ' deformation gradient aim (lab) =', math_transpose33(math_rotate_backward33(F_aim,params%rotation_BC)) + ' deformation gradient aim (lab) =', transpose(math_rotate_backward33(F_aim,params%rotation_BC)) write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') & - ' deformation gradient aim =', math_transpose33(F_aim) + ' deformation gradient aim =', transpose(F_aim) flush(6) endif newIteration @@ -516,13 +515,10 @@ subroutine BasicPETSc_forward(guess,timeinc,timeinc_old,loadCaseTime,deformation C_volAvgLastInc = C_volAvg C_minMaxAvgLastInc = C_minMaxAvg - - if (guess) then ! QUESTION: better with a = L ? x:y - F_aimDot = stress_BC%maskFloat * (F_aim - F_aim_lastInc)/timeinc_old ! initialize with correction based on last inc - else - F_aimDot = 0.0_pReal - endif ! components of deformation_BC%maskFloat always start out with zero + + F_aimDot = merge(stress_BC%maskFloat*(F_aim-F_aim_lastInc)/timeinc_old, 0.0_pReal, guess) F_aim_lastInc = F_aim + !-------------------------------------------------------------------------------------------------- ! calculate rate for aim if (deformation_BC%myType=='l') then ! calculate F_aimDot from given L and current F diff --git a/src/spectral_mech_Polarisation.f90 b/src/spectral_mech_Polarisation.f90 index 300b50142..acd713c70 100644 --- a/src/spectral_mech_Polarisation.f90 +++ b/src/spectral_mech_Polarisation.f90 @@ -333,7 +333,6 @@ subroutine Polarisation_formResidual(in,x_scal,f_scal,dummy,ierr) IO_intOut use math, only: & math_rotate_backward33, & - math_transpose33, & math_mul3333xx33, & math_invSym3333, & math_mul33x33 @@ -402,9 +401,9 @@ subroutine Polarisation_formResidual(in,x_scal,f_scal,dummy,ierr) trim(incInfo), ' @ Iteration ', itmin, '≤',totalIter, '≤', itmax if (iand(debug_level(debug_spectral),debug_spectralRotation) /= 0) & write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') & - ' deformation gradient aim (lab) =', math_transpose33(math_rotate_backward33(F_aim,params%rotation_BC)) + ' deformation gradient aim (lab) =', transpose(math_rotate_backward33(F_aim,params%rotation_BC)) write(6,'(/,a,/,3(3(f12.7,1x)/))',advance='no') & - ' deformation gradient aim =', math_transpose33(F_aim) + ' deformation gradient aim =', transpose(F_aim) flush(6) endif newIteration @@ -548,7 +547,6 @@ subroutine Polarisation_forward(guess,timeinc,timeinc_old,loadCaseTime,deformati use math, only: & math_mul33x33, & math_mul3333xx33, & - math_transpose33, & math_rotate_backward33 use numerics, only: & worldrank @@ -630,11 +628,7 @@ subroutine Polarisation_forward(guess,timeinc,timeinc_old,loadCaseTime,deformati C_volAvgLastInc = C_volAvg C_minMaxAvgLastInc = C_minMaxAvg - if (guess) then ! QUESTION: better with a = L ? x:y - F_aimDot = stress_BC%maskFloat * (F_aim - F_aim_lastInc)/timeinc_old ! initialize with correction based on last inc - else - F_aimDot = 0.0_pReal - endif + F_aimDot = merge(stress_BC%maskFloat*(F_aim-F_aim_lastInc)/timeinc_old, 0.0_pReal, guess) F_aim_lastInc = F_aim !-------------------------------------------------------------------------------------------------- ! calculate rate for aim @@ -676,16 +670,16 @@ subroutine Polarisation_forward(guess,timeinc,timeinc_old,loadCaseTime,deformati F_lambda33 = reshape(F_tau(1:9,i,j,k)-F(1:9,i,j,k),[3,3]) F_lambda33 = math_mul3333xx33(S_scale,math_mul33x33(F_lambda33, & math_mul3333xx33(C_scale,& - math_mul33x33(math_transpose33(F_lambda33),& + math_mul33x33(transpose(F_lambda33),& F_lambda33)-math_I3))*0.5_pReal)& + math_I3 F_tau(1:9,i,j,k) = reshape(F_lambda33,[9])+F(1:9,i,j,k) enddo; enddo; enddo endif - nullify(F) - nullify(F_tau) - call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr) + nullify(F) + nullify(F_tau) + call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr) end subroutine Polarisation_forward diff --git a/src/spectral_utilities.f90 b/src/spectral_utilities.f90 old mode 100755 new mode 100644 index 874558c57..e3383f3d1 --- a/src/spectral_utilities.f90 +++ b/src/spectral_utilities.f90 @@ -70,7 +70,6 @@ module spectral_utilities ! derived types type, public :: tSolutionState !< return type of solution from spectral solver variants logical :: converged = .true. - logical :: regrid = .false. logical :: stagConverged = .true. logical :: termIll = .false. integer(pInt) :: iterationsNeeded = 0_pInt @@ -782,7 +781,7 @@ function utilities_maskedCompliance(rot_BC,mask_stress,C) if(debugGeneral) then write(6,'(/,a)') ' ... updating masked compliance ............................................' write(6,'(/,a,/,9(9(2x,f12.7,1x)/))',advance='no') ' Stiffness C (load) / GPa =',& - transpose(temp99_Real)/1.e9_pReal + transpose(temp99_Real)*1.0e-9_pReal flush(6) endif k = 0_pInt ! calculate reduced stiffness @@ -837,7 +836,7 @@ function utilities_maskedCompliance(rot_BC,mask_stress,C) endif if(debugGeneral) then write(6,'(/,a,/,9(9(2x,f10.5,1x)/),/)',advance='no') & - ' Masked Compliance (load) / GPa =', transpose(temp99_Real*1.e-9_pReal) + ' Masked Compliance (load) * GPa =', transpose(temp99_Real)*1.0e9_pReal flush(6) endif utilities_maskedCompliance = math_Plain99to3333(temp99_Real) @@ -935,7 +934,6 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& debug_reset, & debug_info use math, only: & - math_transpose33, & math_rotate_forward33, & math_det33 use mesh, only: & @@ -994,10 +992,10 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& call MPI_Allreduce(MPI_IN_PLACE,P_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) if (debugRotation) & write(6,'(/,a,/,3(3(2x,f12.4,1x)/))',advance='no') ' Piola--Kirchhoff stress (lab) / MPa =',& - math_transpose33(P_av)*1.e-6_pReal + transpose(P_av)*1.e-6_pReal P_av = math_rotate_forward33(P_av,rotation_BC) write(6,'(/,a,/,3(3(2x,f12.4,1x)/))',advance='no') ' Piola--Kirchhoff stress / MPa =',& - math_transpose33(P_av)*1.e-6_pReal + transpose(P_av)*1.e-6_pReal flush(6) max_dPdF = 0.0_pReal