Multiply list of pure functions by list of scalars?

First, to evaluate the functions embedded in calcShapeFunctions[3], you could do the following:

calcShapeFunctions[3] /. f_Function :> f[x, y]

{1/2 - x/2, x/2 - y/2, 1/2 + y/2}

But to find the derivatives, you need to alter how calcShapeFunctions function is defined so that the bodies of your Functions are evaluated and the part of the list of formulas is extracted (line 4):

calcShapeFunctions[nnpe_] := 
  Module[{m, xiCoord, etaCoord, (*p,*) c}, 
   m = Array[Function[{xi, eta}, 
         Evaluate@{1, xi, eta, xi*eta, xi^2, eta^2, xi^2*eta, xi*eta^2, xi^2*eta^2}[[#]]] &,
        nnpe];
   xiCoord = {-1, 1, 1, -1, 0, 1, 0, -1, 0};
   etaCoord = {-1, -1, 1, 1, -1, 0, 1, 0, 0};
   (*p = m[[1 ;; nnpe]];*) (* same as p = m *)
   c = Array[Through[m[xiCoord[[#]], etaCoord[[#]]]] &, nnpe];
   Return[m.Inverse[c]];];

Then the derivatives can be done as follows:

calcShapeFunctions[3] /. f_Function :> Derivative[1, 0][f][x, y]

{-(1/2), 1/2, 0}

calcShapeFunctions[3] /. f_Function :> Derivative[0, 1][f][x, y]

{0, -(1/2), 1/2}

(Note: In calcShapeFunctions, the variable p was the same as m, so I replaced it.)

Response to comment

In response to the OP's request as to how I would modify Mr. Wizard's answer, I offer the following. As I suggested in a comment, I would Evaluate the function body inside With:

calcFn2[nnpe_] := 
 Module[{m, xiCoord, etaCoord, p, c, xi, eta}, 
  m = {1, #, #2, # #2, #^2, #2^2, #^2 #2, # #2^2, #^2 #2^2} &;
  m = m[[{1}, ;; nnpe]];
  xiCoord = {-1, 1, 1, -1, 0, 1, 0, -1, 0};
  etaCoord = {-1, -1, 1, 1, -1, 0, 1, 0, 0};
  c = Array[m[xiCoord[[#]], etaCoord[[#]]] &, nnpe];
  With[{ic = Inverse[c], mev = m}, Evaluate[mev[#1, #2].ic] &]]

calcFn2[4]
(* {1/4 - #1/4 - #2/4 + (#1 #2)/4, 1/4 + #1/4 - #2/4 - (#1 #2)/4, 
    1/4 + #1/4 + #2/4 + (#1 #2)/4, 1/4 - #1/4 + #2/4 - (#1 #2)/4} & *)

One can take the Derivative of the Function:

Derivative[0, 1][calcFn2[4]]
(* {-(1/4) + #1/4, -(1/4) - #1/4, 1/4 + #1/4, 1/4 - #1/4} & *)

And plug in expressions:

Derivative[0, 1][calcFn2[4]][x, y]
(* {-(1/4) + x/4, -(1/4) - x/4, 1/4 + x/4, 1/4 - x/4} *)

Or one can apply D to the expression calcFn2[3][x, y]:

D[calcFn2[4][x, y], y]
(* {-(1/4) + x/4, -(1/4) - x/4, 1/4 + x/4, 1/4 - x/4} *)

As with the last question on which I believe this one is based, it may be simpler to use a single function rather than a list of functions.

To understand my code you will need to know that:

• You can manipulate a Function with Part just as you would any other expression
• Part always wraps element sequences with the original head of the expression

I chose to use Slot (# and #2) because it makes this manipulation much simpler, though it would be possible with named parameters with more work.

Other notes:

• I made output of the calc function to be a single function as well, simplifying both the code and the application of it
• I eliminated p as it appeared to be a copy of m
• With is use to evaluate certain parts of the output Function

The code:

calcFn2[nnpe_] := 
 Module[{m, xiCoord, etaCoord, p, c, xi, eta}, 
  m = {1, #, #2, # #2, #^2, #2^2, #^2 #2, # #2^2, #^2 #2^2} &;
  m = m[[{1}, ;; nnpe]];
  xiCoord = {-1, 1, 1, -1, 0, 1, 0, -1, 0};
  etaCoord = {-1, -1, 1, 1, -1, 0, 1, 0, 0};
  c = Array[m[xiCoord[[#]], etaCoord[[#]]] &, nnpe];
  With[{ic = Inverse[c], mev = m}, mev[##].ic &]
 ]

calcFn2[3][x, y]

{1/2 - x/2, x/2 - y/2, 1/2 + y/2}