For the purpose of this minimal example, let's say we have a list of functions, like this:
f[y_?NumericQ] := {NIntegrate[z*y, {z, 0, 1}], a y}
I want to integrate an expression involving f, say
NDSolve[{y'[t] == f[y[t]][[1]], y[0] == 1}, y[t], {t, 0, 5}]
Now, the problem is that this doesn't return the expected result, because f[y[t]][[1]] evaluates to y[t] inside the NDSolve.
How can this be done correctly?
You can include the part extraction as an argument of your function, perhaps as a SubValues definition:
ClearAll[f]
f[y_?NumericQ][part_] := {NIntegrate[z*y, {z, 0, 1}], a y}[[part]]
NDSolve[{y'[t] == f[y[t]][1], y[0] == 1}, y[t], {t, 0, 5}]
{{y[t]->InterpolatingFunction[{{0.,5.}},<>][t]}}
Or, inside the primary body as an optional argument:
ClearAll[f]
f[y_?NumericQ, part_: All] := {NIntegrate[z*y, {z, 0, 1}], a y}[[part]]
NDSolve[{y'[t] == f[y[t], 1], y[0] == 1}, y[t], {t, 0, 5}]
This second method returns both values by default:
f[3.6]
{1.8, 3.6 a}
An alternative that comes to mind is to use a custom Part function that won't trigger when it should not, e.g.:
ClearAll[f]
listPart[x_List, part__] := x[[part]]
f[y_?NumericQ] := {NIntegrate[z*y, {z, 0, 1}], a y}
NDSolve[{y'[t] == f[y[t]] ~listPart~ 1, y[0] == 1}, y[t], {t, 0, 5}]
Mathematica 10 introduces something like listPart, with additional functionality, in Indexed:
Indexed can be used to indicate components of symbolic vectors, matrices, tensors, etc.
When expr is a list, Indexed[expr,i] gives expr[[i]].
When expr is a list, Indexed[expr,{i,j,...}] gives Indexed[expr[[i]],{j,...}].
Indexed can be used instead of Part for single-element extraction (at any level).
Applied to the case at hand:
f[y_?NumericQ] := {NIntegrate[z*y, {z, 0, 1}], a y}
NDSolve[{y'[t] == Indexed[f[y[t]], 1], y[0] == 1}, y[t], {t, 0, 5}]
{{y[t] -> InterpolatingFunction[{{0., 5.}}, <>][t]}}