Intel® Software Network

i posted this to some math forums, but got no response.  perhaps one of you knows the answer.  it is not a FORTRAN question. 


SHORT VERSION:

I've got some code that finds a determinant of a matrix F by taking the
square root of an expression of the form (e*g-f*f), where e,g, and f are
specified elements of (F^T)*F.

I've never seen this before.  Have you?  Can you provide me a reference?


LONG VERSION:

Suppose you have an elastic cube.  Suppose you gently deform it by pushing
and pulling on the corners.  Just bend it a little bit; don't jump on it
until it pops or anything like that.

Choose one of the sides, and a point on that side.
Use "X" to denote the original position of the point (before you deformed
it), and "x" to denote the new position.

Suppose you determine the deformation gradient F of the face mapping that
takes X into x.

If that last sentence (or anything else) didn't make sense to you, ignore
all of the above and suppose you have a symmetric, positive definite matrix
F.  If you're still confused, suppose you have a nice 3x3 matrix F.


Suppose you want the determinant of F.  What would you do?  This is what
someone else did:


if (  (face.eq.BOTTOM_FACE)
.or.(face.eq.TOP_FACE))
then
ee = sum(F(:,1)*F(:,1))
ff = sum(F(:,1)*F(:,2))
gg = sum(F(:,2)*F(:,2))

else if ((face.eq.INTERIOR_FACE)
.or.(face.eq.EXTERIOR_FACE))
then
ee = sum(F(:,1)*F(:,1))
ff = sum(F(:,1)*F(:,3))
gg = sum(F(:,3)*F(:,3))

else if ((face.eq.RIGHT_FACE)
.or.(face.eq.LEFT_FACE))
then
ee = sum(F(:,2)*F(:,2))
ff = sum(F(:,2)*F(:,3))
gg = sum(F(:,3)*F(:,3))
else
write(*,*) "ERROR: face not found"
endif

! Jacobian of face mapping
detJb = sqrt(ee*gg-ff*ff)


Anyone know what formula this is?  Can you provide a reference?

He then goes on to use this to find the normal:


if (face.eq.BOTTOM_FACE) then
nor(1) = F(2,2)*F(3,1) - F(3,2)*F(2,1)
nor(2) = F(3,2)*F(1,1) - F(1,2)*F(3,1)
nor(3) = F(1,2)*F(2,1) - F(2,2)*F(1,1)
else if(face.eq.INTERIOR_FACE) then
nor(1) = F(2,1)*F(3,3) - F(3,1)*F(2,3)
nor(2) = F(3,1)*F(1,3) - F(1,1)*F(3,3)
nor(3) = F(1,1)*F(2,3) - F(2,1)*F(1,3)
else if(face.eq.RIGHT_FACE) then
nor(1) = F(2,2)*F(3,3) - F(3,2)*F(2,3)
nor(2) = F(3,2)*F(1,3) - F(1,2)*F(3,3)
nor(3) = F(1,2)*F(2,3) - F(2,2)*F(1,3)
else if(face.eq.EXTERIOR_FACE) then
nor(1) = -F(2,1)*F(3,3) + F(3,1)*F(2,3)
nor(2) = -F(3,1)*F(1,3) + F(1,1)*F(3,3)
nor(3) = -F(1,1)*F(2,3) + F(2,1)*F(1,3)
else if(face.eq.LEFT_FACE) then
nor(1) = -F(2,2)*F(3,3) + F(3,2)*F(2,3)
nor(2) = -F(3,2)*F(1,3) + F(1,2)*F(3,3)
nor(3) = -F(1,2)*F(2,3) + F(2,2)*F(1,3)
else if(face.eq.TOP_FACE) then
nor(1) = -F(2,2)*F(3,1) + F(3,2)*F(2,1)
nor(2) = -F(3,2)*F(1,1) + F(1,2)*F(3,1)
nor(3) = -F(1,2)*F(2,1) + F(2,2)*F(1,1)
endif

! make the normal a unit vector
tmp = sqrt( sum(nor(:)*nor(:)) )
nor(:) = nor(:)/tmp


Reference?

Thanks bunches!
Nooj