2
$\begingroup$

As the result of a symbolic integration, I have a very long array with symbolic variables u[i] which I would like to substitute for real number values. Moreover, I would like to do this using a compiled function for efficiency.

For simplicity say we have the array

aLong := {u[1]/u[2], u[1]*u[2]}

I want to compile a function whose job is to evaluate this expression.

However, if I try:

f=Compile[{{w,_Real,1}},aLong/.{u[1]->w[[1]],u[2]->w[[2]]}]
f[{1,2}]

it gives errors:

CompiledFunction::cfse: Compiled expression {0.5,2.} should be a machine-size real number.

CompiledFunction::cfex: Could not complete external evaluation at instruction 1; proceeding with uncompiled evaluation.

If instead I try Evaluate inside of the Compile, then I get the following error:

f = Compile[{{w, _Real, 1}}, Evaluate[aLong /. {u[1] -> w[[1]], u[2] -> w[[2]]}]];
f[{1, 1}]

Part::partd: Part specification w[[1]] is longer than depth of object.

Part::partd: Part specification w[[2]] is longer than depth of object.

Any suggestions that can help me fix this? (there must be something very basic I have not understood) Thank you!

Improve this question
$\endgroup$

2 Answers 2

Reset to default
3
$\begingroup$

It's not really an error, only a warning message, if w is atomic. It's not an error in the sense that w[[1]] is the result of Part evaluating before w has a value and that this result is what is intended. Nonetheless, the error message is irritating.

One way:

f = Block[{w}, (* in case w has a value *)
  Quiet[
   Compile[{{w, _Real, 1}}, 
    Evaluate[aLong /. {u[1] -> w[[1]], u[2] -> w[[2]]}]],
   Part::partd] (* turn off the specific message *)
   ];
f[{1, 1}]
(*  {1., 1.}  *)

An alternative that avoids "errors":

f = Block[{Compile, Part, w}, (* suspend evaluation until safe *)
   Compile[{{w, _Real, 1}}, 
    aLong /. {u[1] -> w[[1]], u[2] -> w[[2]]}]
   ];
f[{1, 1}]
(*  {1., 1.}  *)

Update. Third way (Block[] not needed):

f = Hold@Compile[{{w, _Real, 1}}, aLong] /.
     OwnValues@aLong /. 
    {u[1] :> w[[1]], u[2] :> w[[2]]} //
   ReleaseHold;
f[{1, 1}]
(*  {1., 1.}  *)
Improve this answer
$\endgroup$
2
  • $\begingroup$ Ok, so let's see if I understood, the idea is: 1. Use Hold to maintain the whole expression unevaluated. 2. Use Replace + HoldPattern to replace the string aLong for the vector 3. Use Replace together with :> to substitute the u[i] for w[[i]] without the warning appearing because of the hold on :> 4. Finally, release the hold. This way of working with unevaluated expressions is new to me. Very useful though. Thank you! $\endgroup$
    – mmen
    Sep 20, 2022 at 13:52
  • $\begingroup$ @mmendez Yes (for the summary of the third alternative). The ability to generate code is a notable feature of Mma. See mathematica.stackexchange.com/… or code-generation for more. -- Here a shorter version: f = Hold@Compile[{{w, _Real, 1}}, aLong] /. OwnValues@aLong /. u[i_] :> w[[i]] // ReleaseHold $\endgroup$
    – Michael E2
    Sep 20, 2022 at 14:09
2
$\begingroup$

This would fix the issue (it might not be the cleanest solution though)

aLong := {u[1]/u[2], u[1]*u[2]}
expr = aLong /. {u[1] -> w[[1]], u[2] -> w[[2]]} // Quiet;

f = Compile[{{w, _Real, 1}},
  Evaluate@expr
];

f[{1,2}]
>> {0.5, 2.}

Here Quiet[] prevents errors from w[[1]] being undefined yet

Improve this answer
$\endgroup$
2
  • $\begingroup$ If I understand correctly you are only hiding the error messages. Right? Is this common practice in Mathematica? I was hoping for an approach that would not produce errors in the first place. $\endgroup$
    – mmen
    Sep 20, 2022 at 12:20
  • $\begingroup$ @mmendez I would not say it is common place to simply hide errors but in some cases it can be a good solution. Micheal produced a solution without errors but arguably it is more complicated $\endgroup$
    – AccidentalTaylorExpansion
    Sep 20, 2022 at 12:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.