From 6bc57d4911460914c39f8585ecf77f5a316d331e Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 9 Mar 2020 14:17:41 +0100 Subject: [PATCH] mainly code duplication and not used --- src/math.f90 | 64 ---------------------------------------------------- 1 file changed, 64 deletions(-) diff --git a/src/math.f90 b/src/math.f90 index a6a3498e0..ac0974995 100644 --- a/src/math.f90 +++ b/src/math.f90 @@ -1022,70 +1022,6 @@ pure function math_eigenvectorBasisSym33(m) end function math_eigenvectorBasisSym33 -!-------------------------------------------------------------------------------------------------- -!> @brief logarithm eigenvector basis of symmetric 33 matrix m -!-------------------------------------------------------------------------------------------------- -pure function math_eigenvectorBasisSym33_log(m) - - real(pReal), dimension(3,3) :: math_eigenvectorBasisSym33_log - real(pReal), dimension(3) :: invariants, values - real(pReal), dimension(3,3), intent(in) :: m - real(pReal) :: P, Q, rho, phi - real(pReal), parameter :: TOL=1.e-14_pReal - real(pReal), dimension(3,3,3) :: N, EB - - invariants = math_invariantsSym33(m) - EB = 0.0_pReal - - P = invariants(2)-invariants(1)**2.0_pReal/3.0_pReal - Q = -2.0_pReal/27.0_pReal*invariants(1)**3.0_pReal+product(invariants(1:2))/3.0_pReal-invariants(3) - - threeSimilarEigenvalues: if(all(abs([P,Q]) < TOL)) then - values = invariants(1)/3.0_pReal - ! this is not really correct, but at least the basis is correct - EB(1,1,1)=1.0_pReal - EB(2,2,2)=1.0_pReal - EB(3,3,3)=1.0_pReal - else threeSimilarEigenvalues - rho=sqrt(-3.0_pReal*P**3.0_pReal)/9.0_pReal - phi=acos(math_clip(-Q/rho*0.5_pReal,-1.0_pReal,1.0_pReal)) - values = 2.0_pReal*rho**(1.0_pReal/3.0_pReal)* & - [cos(phi/3.0_pReal), & - cos((phi+2.0_pReal*PI)/3.0_pReal), & - cos((phi+4.0_pReal*PI)/3.0_pReal) & - ] + invariants(1)/3.0_pReal - N(1:3,1:3,1) = m-values(1)*math_I3 - N(1:3,1:3,2) = m-values(2)*math_I3 - N(1:3,1:3,3) = m-values(3)*math_I3 - twoSimilarEigenvalues: if(abs(values(1)-values(2)) < TOL) then - EB(1:3,1:3,3)=matmul(N(1:3,1:3,1),N(1:3,1:3,2))/ & - ((values(3)-values(1))*(values(3)-values(2))) - EB(1:3,1:3,1)=math_I3-EB(1:3,1:3,3) - elseif(abs(values(2)-values(3)) < TOL) then twoSimilarEigenvalues - EB(1:3,1:3,1)=matmul(N(1:3,1:3,2),N(1:3,1:3,3))/ & - ((values(1)-values(2))*(values(1)-values(3))) - EB(1:3,1:3,2)=math_I3-EB(1:3,1:3,1) - elseif(abs(values(3)-values(1)) < TOL) then twoSimilarEigenvalues - EB(1:3,1:3,2)=matmul(N(1:3,1:3,1),N(1:3,1:3,3))/ & - ((values(2)-values(1))*(values(2)-values(3))) - EB(1:3,1:3,1)=math_I3-EB(1:3,1:3,2) - else twoSimilarEigenvalues - EB(1:3,1:3,1)=matmul(N(1:3,1:3,2),N(1:3,1:3,3))/ & - ((values(1)-values(2))*(values(1)-values(3))) - EB(1:3,1:3,2)=matmul(N(1:3,1:3,1),N(1:3,1:3,3))/ & - ((values(2)-values(1))*(values(2)-values(3))) - EB(1:3,1:3,3)=matmul(N(1:3,1:3,1),N(1:3,1:3,2))/ & - ((values(3)-values(1))*(values(3)-values(2))) - endif twoSimilarEigenvalues - endif threeSimilarEigenvalues - - math_eigenvectorBasisSym33_log = log(sqrt(values(1))) * EB(1:3,1:3,1) & - + log(sqrt(values(2))) * EB(1:3,1:3,2) & - + log(sqrt(values(3))) * EB(1:3,1:3,3) - -end function math_eigenvectorBasisSym33_log - - !-------------------------------------------------------------------------------------------------- !> @brief rotational part from polar decomposition of 33 tensor m !--------------------------------------------------------------------------------------------------