I am attempting to plot the determinant of a Jacobian matrix that is defined by deriving scalar functions and combining these derivatives into a matrix function (the Jacobian) and thentaking its determinant.

I define shape functions Ni which are functions of xi and eta. For each shape function the derivative with respect to xi and eta is required. These derivatives are then assembled into a matrix GN that is multiplied by a matrix of coordinates XY. The Jacobian is then J=GN.XY and I want to plot Det[J].

My attempt is shown below

`(* Define Ni *) N1[xi_, eta_] := 1/2 (1 - xi) (1 - eta); N2[xi_, eta_] := 1/2 (1 - xi) (1 + eta); N3[xi_, eta_] := 1/2 (1 + xi) (1 + eta); N4[xi_, eta_] := 1/2 (1 - xi) (1 + eta); (* Compute derivatives *) dN1dxi = D[N1[xi, eta], xi]; dN1deta = D[N1[xi, eta], eta]; dN2dxi = D[N2[xi, eta], xi]; dN2deta = D[N2[xi, eta], eta]; dN3dxi = D[N3[xi, eta], xi]; dN3deta = D[N3[xi, eta], eta]; dN4dxi = D[N4[xi, eta], xi]; dN4deta = D[N4[xi, eta], eta]; (* Assemble GN *) GN = {{dN1dxi, dN2dxi, dN3dxi, dN4dxi}, {dN1deta, dN2deta, dN3deta, dN4deta}}; (* Define coordinates and assemble XY *) X1 = {-1, 1}; X2 = {1, -1}; X3 = {1, 1}; X4 = {-1, 1}; XY = {X1, X2, X3, X4}; (* Compute J and Det[J] *) J = GN.XY; Print[Simplify[J]]; (* For interest *) detJ = Det[J]; (* Plot *) DensityPlot[detJ[xi, eta], {xi, -1, 1}, {eta, -1, 1}] `

However, the plot is blank. See below.

Any help would be appreciated!