# Using predefined expressions as free variables in a pure function [duplicate]

I want to have a function value of an expression where some variables are solutions to some set of equations, with some values of parameters. I had an idea to use pure functions for that.

However, since both the expression and the equations are lengthy I'd like to place them in a separate expressions defined earlier on. Here is Minimal Working Example:

MWE1 = b x^2 + a y^2 + 3 + δ;
Dx = D[MWE1, x];
Dy = D[MWE1, y];
MWF[a_, b_] :=
Function[{δ}, (Evaluate[MWE1]) /.
FindRoot[
{Evaluate[Dx] == 0, Evaluate[Dy] == 0}, {{x, -3, 3}, {y, -3, 3}},
AccuracyGoal -> Infinity]]

And here is an attempt to test it:

test3 = MWF[1, 1];
test3[1]
3 + b x^2 + a y^2 + δ /.
FindRoot[
{Evaluate[Dx] == 0, Evaluate[Dy] == 0}, {{x, -3, 3}, {y, -3, 3}},
AccuracyGoal -> ∞]

However, if I paste the expressions MWE1, Dx, Dy directly to the definition of MWF, the function works like charm.

Is it possible to use predefined expressions inside the definition? Do so would make the code much clearer.

• Have a look at With – Sascha Aug 22 '16 at 13:43
• @Sascha Yes but that tends to mess things up due to how scoping is implemented using variable renaming. Replace is better, but can also end up with surprising results, see my answer. I don't have a good solution! :'( – Szabolcs Aug 22 '16 at 13:48
• @Szabolcs very good point, I haven't though of that at all – Sascha Aug 22 '16 at 13:52
• Related: (10086), (11461), (13757), (20766), (46751), (69590) – Mr.Wizard Aug 22 '16 at 21:38

Your issue is quite simple and fundamental: you try to use lexical scoping (SetDelayed and Function) as it would be dynamic scoping (Block). So your problem can be solved easily by outsourcing the dynamic scoping to Block with minimal modification of your original code:

MWE1 = b x^2 + a y^2 + 3 + δ;
Dx = D[MWE1, x];
Dy = D[MWE1, y];
MWF[aa_, bb_] =
Function[{δδ},
Block[{δ = δδ, a = aa, b = bb},
MWE1 /. FindRoot[{Dx == 0, Dy == 0}, {{x, -3, 3}, {y, -3, 3}},
AccuracyGoal -> Infinity]]]

test3 = MWF[1, 1];
test3[1]

NDSolve`ScaledVectorNorm::urange: The scaling vector has a zero component, so the zero absolute tolerance leads to an unbounded value.

4.

### Update

I just realized that it is possible to achieve what you want by means of anonymous pure function # &:

MWE1 = b x^2 + a y^2 + 3 + δ;
Dx = D[MWE1, x];
Dy = D[MWE1, y];
MWF[a_, b_] = Function[{δ},
#1 /. FindRoot[{#2 == 0, #3 == 0}, {{x, -3, 3}, {y, -3, 3}}, AccuracyGoal -> Infinity]
] &[MWE1, Dx, Dy]

test3 = MWF[1, 1];
test3[1]

This method utilizes the fact that anonymous pure functions do not rename variables in the nested scoping constructs.

• Alexey, would you please review my Related: links and tell me if you believe this question is a duplicate? – Mr.Wizard Aug 22 '16 at 21:42
• @Alexey Popkov I tried that method, defining MWF[aa_, bb_] := Function[{δδ}, Module[{a = aa, b = bb, δ = δδ}, MWE1 /. FindRoot[{Dx == 0, Dy == 0}, {{x, -3, 3}, {y, -3, 3}}, AccuracyGoal -> Infinity]]] but it does not work. What version of mathematica are you running? I've recently upgraded, may that be source of that difference? – Picek Aug 23 '16 at 10:00
• @Picek Module performs lexical scoping, so it can't work in place of Block where exactly dynamic scoping is necessary! Replace Module with Block and this will work. – Alexey Popkov Aug 23 '16 at 10:02
• @Picek I have added another solution which doesn't use dynamic scoping and despite this doesn't imply any artificial tricks (like fooling enclosing scoping construct in order to avoid renaming of variables). – Alexey Popkov Aug 23 '16 at 12:14
• @Mr.Wizard Taking into account that 3 people already voted for (10086), I've voted in favor of (69590). In my opinion none of the questions you linked to is an exact duplicate of this one but (10086) and (69590) together cover the scope of this question. – Alexey Popkov Aug 24 '16 at 7:31

### Why doesn't it work?

Evaluate only has effect at the first level in a held expression.

{Hold[1 + 1], Hold[Evaluate[1 + 1]], Hold[f[1 + 1]]}
(* {Hold[1 + 1], Hold[2], Hold[f[1 + 1]]} *)

### Why does Evaluate only work at the first level?

Because the evaluator does not look for exceptions beyond the first level inside an expression having a head with HoldAll. If it did, it may potentially spend a lot of time scanning large expressions.

### What's the solution?

Typically one has to use Replace or ReplaceAll.

Function[x, y] /. y -> x^2
(* Function[x, x^2] *)

All this falls under the label of , a difficult topic where things may go wrong all too easily due to how Mathematica implements scoping using variable renaming.

Example of how things may go wrong:

Can't we just use With instead of Replace? Let's try:

With[{y = x^2}, Function[x, y]]

(* Function[x$, x^2] *) Oops! The function argument got renamed from x to x$ to prevent a name collision. At the same time it "broke" our function. This is the reason why Replace and ReplaceAll are better than With for such purposes.

But is Replace bulletproof? Not at all, unfortunately. When doing such a transformation, we might want to localize y. Let's see what happens.

Module[{y},
Function[x, y] /. y -> x^2
]

(* Function[x\$, x^2] *)

Oops! The renaming kicked in again. (Note: Block won't rename here but it's not always a good replacement to Module.)

Doing such things can also lead to situations where SetSystemOptions["StrictLexicalScoping" -> True] may make a significant difference (may fix or break things). This is why Stan Wagon's FindRoots2D function won't work (as posted on this site) if strict lexical scoping is turned on.

Sometimes we can use a trick like this to avoid variable renaming:

Module[{y}, Function @@ Hold[x, y] /. y -> x^2]
(* Function[x, x^2] *)

Writing Function @@ Hold[...] instead of Function[...] prevents the renaming (for better or worse).

To sum up, code generation seems to be quite difficult in Mathematica. Caution is advised. Now if someone could prove me wrong, I would be very happy ... Looking forward to other answers!

• Szabolcs, would you please review my Related: links and tell me if you believe this question is a duplicate? – Mr.Wizard Aug 22 '16 at 21:40
• Would you vote for 13757 as a 2nd duplicate question? I have already voted for 10086. – m_goldberg Aug 23 '16 at 12:52

It can be done, but I suspect you will find it a bit more tricky than you expected. Here is one way to do it. Note that I don't introduce a variable named δ. It is easier to work with Slot ( # ).

{a, b, x, y} = Range[4]; (* to ensure scope conflicts assert themselves *)

MWF[a_, b_] =
Block[{MWE1, Dx, Dy, a, b, x, y},
MWE1 = b x^2 + a y^2 + 3 + #;
Dx = D[MWE1, x];
Dy = D[MWE1, y];
With[{f = MWE1, dx = Dx, dy = Dy},
Function[f /.
FindRoot[{dx == 0, dy == 0}, {{x, -3, 3}, {y, -3, 3}},
AccuracyGoal -> ∞]]]]
3 + b x^2 + a y^2 + #1 /.
FindRoot[{2 b x == 0, 2 a y == 0}, {{x, -3, 3}, {y, -3, 3}},
AccuracyGoal -> ∞] &
testF = MWF[1, 1]
3 + 1 x^2 + 1 y^2 + #1 /.
FindRoot[{2 x == 0, 2 y == 0}, {{x, -3, 3}, {y, -3, 3}},
AccuracyGoal -> ∞] &
Block[{x, y}, testF[1]]

4.

• m_goldberg, would you please review my Related: links and tell me if you believe this question is a duplicate? – Mr.Wizard Aug 22 '16 at 21:42
• @Mr.Wizard. I would vote for the aggregate of 10086, 13757, and 46751 as constituting a duplicate. Can we close that way? If three people each vote for a different one, would that do it? – m_goldberg Aug 23 '16 at 0:56
• Yes, that would do it. – Mr.Wizard Aug 23 '16 at 12:39
• @Mr.Wizard. I will vote for 10086 to get things started. – m_goldberg Aug 23 '16 at 12:49

Another method would be to use delayed evaluation to work around the issue of variable renaming.

MWE1[a_, b_][x_, y_][δ_] := b x^2 + a  y^2 + 3 + δ;
Dx[a_, b_][x_, y_][δ_] := D[MWE1[a, b][x, y][δ], x];
Dy[a_, b_][x_, y_][δ_] := D[MWE1[a, b][x, y][δ], y];

Note that one has to introduce a Module to properly scope x and y.

MWF[a_, b_] := Module[{x, y},
Function[{δ}, (MWE1[a, b][x, y][δ]) /.
FindRoot[{
Dx[a, b][x, y][δ] == 0,
Dy[a, b][x, y][δ] == 0},
{{x, -3, 3}, {y, -3, 3}},
AccuracyGoal -> Infinity]]
]

testing confirms that all symbols are properly scoped and that the function can be called

a = "test"; b = "test"; x = "test"; y = "test"; δ = "test";;
MWF[1, 1][1]

• That works for sure, but it does not have the desired result since the whole point is that in my original calculation MWE1 is much more complicated and executing derivatives on every allocation of MWF must be avoided – Picek Aug 23 '16 at 9:37