From 52b5ef57004fcdd28552f6be4ad79eaac7df8c2a Mon Sep 17 00:00:00 2001 From: Dave Whipp <dwhipp@dal.ca> Date: Mon, 14 Feb 2011 18:44:14 +0000 Subject: [PATCH] Fixed dilatation rate calculation, added option of using 2D strain rate invariants and removed some commented lines --- src/do_leaf_measurements.f90 | 15 +++++---------- 1 file changed, 5 insertions(+), 10 deletions(-) diff --git a/src/do_leaf_measurements.f90 b/src/do_leaf_measurements.f90 index 657efb31..e9e31efa 100644 --- a/src/do_leaf_measurements.f90 +++ b/src/do_leaf_measurements.f90 @@ -150,12 +150,12 @@ do i=1+iproc,nleaves,nproc exy=exy+(dhdx(k)*vy(k)+dhdy(k)*vx(k))/2.d0 eyz=eyz+(dhdy(k)*vz(k)+dhdz(k)*vy(k))/2.d0 ezx=ezx+(dhdz(k)*vx(k)+dhdx(k)*vz(k))/2.d0 -! exy=exy+(dhdx(k)*vy(k)-dhdy(k)*vx(k))/2.d0 -! eyz=eyz+(dhdy(k)*vz(k)-dhdz(k)*vy(k))/2.d0 -! ezx=ezx+(dhdz(k)*vx(k)-dhdx(k)*vz(k))/2.d0 enddo - call deviatoric (exx,eyy,ezz,exy,eyz,ezx) + ! Dilatation rate + dilatrtemp(i)=(exx+eyy+ezz) + + call deviatoric (params%invariants_2d,exx,eyy,ezz,exy,eyz,ezx) ! if (.not.osolve%is_plastic(i)) then ! qtemp(i)=0.d0 @@ -163,14 +163,9 @@ do i=1+iproc,nleaves,nproc qtemp(i)=2.d0*osolve%eviscosity(i)*sqrt(second_invariant(exx,eyy,ezz,exy,eyz,ezx)) ! end if - ! Dilatation rate - dilatrtemp(i)=(exx+eyy+ezz) - J2d=second_invariant(exx,eyy,ezz,exy,eyz,ezx) J3d=third_invariant (exx,eyy,ezz,exy,eyz,ezx) - !e2dtemp(i)=J2d - ! Set e2d as effective stress (sqrt of 2nd invariant) - dwhipp 12/09 e2dtemp(i)=sqrt(J2d) if (J2d/=0.d0) lodetemp(i)=lode_angle (J2d,J3d)/pi @@ -223,4 +218,4 @@ return end !------------------------------------------------------------------------------| -!------------------------------------------------------------------------------| +!------------------------------------------------------------------------------| \ No newline at end of file -- GitLab