Problem changing value of a fixed variable in a function

[clopen]

Hello, I am new to programming in Julia and currently working on transitioning code from Fortran to Julia.

In a file called codeA.jl, I have declared several variables, some of which are utilized as parameters within a function. For example:

a = 0.0
b = zeros(Float64, a)

function f(x::Int)
        c = b.*x .+ a
        return c
end 

This code from codeA.jl is then included in another file called codeB.jl.

In codeB.jl, the values of the a and b variables from codeA.jl are modified, and I call the function f within codeB.jl.

In an effort to reduce the number of allocations, following recommendations I found in this forum, I included the variables a and b as fixed inputs to f. This means that I modified f to be defined as:

a = 0.0
b = zeros(Float64, a)

function f(x::Int, av = a, bv = b)
        c = bv.*x .+ av
        return c
end 

After this modification, the number of allocations decreases almost to zero. However, if I change the values of a and b (from codeA.jl) in codeB.jl, the values of a and b are not updated when I call f in codeB.jl. This raises the following questions:

(1) Is this approach of fixing a and b as variables of f a natural solution to reduce the number of allocations? I need to call f several times in codeB.jl (inside a loop), so I prefer to define a and b outside the f function. Additionally, I prefer not to specify a and b when calling f (e.g., f(x, a, b)), as only x is meaningful in the problem modeled by f.

(2) In the way I have modeled this problem, how can I ensure that when I call the f function in codeB.jl, it recognizes the new values of a and b?

(3) Given the conditions outlined in (1), is there a better or more natural way to organize the code?

[PeterSimon]

Welcome to Julia! I also am a former Fortran programmer, and I now couldn’t dream of going back!

There is a problem with your code. zeros expects an integer as its second argument, to denote the length of the array to be constructed. So I really expect your second line of codeA.jl to generate a runtime error. Not sure what you intended here.

Also, I recommend that you remove the first two lines of codeA.jl and replace the function definition line with

function f(x::Int, av, bv)

You don’t need the default argument values—simply call the function with all three arguments set to the values that you want.

If you explain what it is that you are trying to accomplish with this style of programming, we can provide more guidance on achieving it in a more idiomatic way for Julia.

[Benny]

The typical advice to changing global variables to function inputs is for you to call f(3, a, b), not to use default values. That said, default values should be evaluated “freshly,” in other words calling f(3) is equivalent to f(3, a, b). It’s likely that there’s something about the include and evaluation order that isn’t said here.

[clopen]

Sorry, it should be a = 0 (I simplify my original code just for the example).

I would like to avoid have to call f by f(x,av,bv) and just call it by f(x), as only x it has a meaning in my problem (a and b can change, but are used as parameter).

[PeterSimon]

@clopen

have to call f by f(x,av,bv) and just call it by f(x), as only x it has a meaning in my problem (a and b can change

One way: In your second file, you can do this

a = ...
b = ...
f(x) = x -> f(x, a, b)

which will create a new method for f that can be called with a single argument, and which will use whatever the current values of a and b in your second file are. This is called a “closure”.

Another way:

function f(x)
    return f(x,a,b)
end

[Benny]

@PeterSimon

function f(x)
    return f(x,a,b)
end

Already happens with default values

julia> f(x, av = a, bv = b) = 0
f (generic function with 3 methods)

julia> @code_lowered f(2)
CodeInfo(
1 ─ %1 = (#self#)(x, Main.a, Main.b)
└──      return %1
)

Closures can be more performant when they capture unreassigned locals, and that sacrifices updates by reassignment. The closure in f(x) = x -> f(x, a, b) accesses globals, which doesn’t behave differently from the first version of f accessing a and b as globals. The argument x is also unused, and this method overwrites the f(x) method using the default values.

It’s still not clear how the issue is occurring. We need a MWE, a simplified codeA.jl and codeB.jl in their entirety throughout a step-by-step process where a and b are attempted to be updated, as well as the results. Here is an example:

# codeA.jl
a = 0.0
b = zeros(Float64, Int(round(a))) # fix the MethodError

function f(x::Int, av = a, bv = b)
  c = bv.*x .+ av
  return c
end

# codeB.jl
include("codeA.jl")
println("before changes: ", f(1))
a = 1.0
b = ones(Float64, Int(round(a))) # fix the MethodError
println("after changes: ", f(1))

Now in the REPL:

julia> include("codeB.jl")
before changes: Float64[]
after changes: [2.0]

And it’s evident the reassignments of the global a and b changed the results of f(1).

[uniment]

Most of your allocations and slowdowns are due to boxing and dynamic dispatch, which occur because Julia can’t know whether a or b will maintain the same type, let alone the same value, during the function’s execution.

Consider declaring a’s and b’s types in codeA.jl to enforce type stability. Something like this perhaps:

a::Float64 = 0.0
b::Vector{Float64} = fill(zero(a), n)

Or, if you want the types to be found automatically:

_a = 0.0
a::typeof(_a) = _a

_b = fill(zero(a), n)
b::typeof(_b) = _b

You can make a macro to automate this:

macro typed(var)
    v=gensym() # unique id
    :($v=$(var.args[2]); $(var.args[1])::typeof($v)=$v) |> esc
end

julia> @typed a=0.0
0.0

julia> a=1
1

julia> a # type-stable
1.0

[bertschi]

@clopen

I would like to avoid have to call f by f(x,av,bv) and just call it by f(x), as only x it has a meaning in my problem (a and b can change, but are used as parameter).

Well, if changing a or b changes what your function computes they are morally inputs to the function. In mathematics you may write f_{\theta}(x) or f(x | \theta) to distinguish parameters from arguments, but it’s still a function of x and \theta. Nevertheless in Julia you have several options to identify different types of arguments, .e.g.,

• use a naming convention (or custom types) for parameters: f(x, θ::Params) = ..

• use currying, i.e., f(a, b) = x -> ...

• or construct partial applications when needed, i.e., f_fixedθ = Base.Fix2(f, θ)

Adding to what has been suggested already, it is often better to not work in the global scope, i.e., define small functions and keep your parameters explicit. Here is a quick example:

# CodeA.jl

function f(x::Int, θ::NamedTuple)
    a, b = θ
    return @. b * x + a
end

function default_params()
    return (a = 2, b = 3)
end

# CodeB.jl

# include("CodeA.jl")

function main()
    θ_default = default_params()
    f_default = Base.Fix2(f, θ_default)
    compute(f_default, [1, 2, 3])

    θ_special = (a = [2, 3], b = 5)
    f_special = Base.Fix2(f, θ_special)
    compute(f_special, [1, 2, 3])
end

function compute(fx, xs)
    for x in xs
        @show fx(x)
    end
end

main()  # Fire it up

No global constants and it’s always clear what parameters and arguments are and where they come from.

[gdalle]

What I usually do in these cases is a fourth option called a callable struct:

struct MyFunction{A,B}
    a::A
    b::B
end

(mf::MyFunction)(x) = f(x, mf.a, mf.b)

mf = MyFunction(a, b)
mf(x)

[Benny]

Other than the way the method is organized, the callable struct more or less does what Fix2 does, and both approaches must reinstantiate to incorporate new values, as they do not reference global variables. Fortunately, instantiation by throwing existing instances together is very fast, and the method accesses the fields as fast as separate arguments, faster than dereferencing global variables. Neither approach benefits from having global variables in the first place, which is why the examples incorporate changes into the instances’ state instead (you wouldn’t want to incorporate default_params into methods and alter it like global variables because that would trigger compiled code invalidation).

I would say that a typed global or const Ref global (for some reason the former takes more dereferencing, but the latter requires you to manually dereference in the method) is preferable when you have arbitrarily many functions relying on the same external state, and it’d be very unwieldy to instantiate callable structs for all of them. That said, I’d still much prefer separate arguments with some defaults I leave alone.

[uniment]

@bertschi

function f(x::Int, θ::NamedTuple)
    a, b = θ
    return @. b * x + a
end

Note that you should write (; a, b) = θ to destructure by field names instead of iteration, to avoid sensitivity to ordering.

[nsajko]

@clopen

In codeB.jl, the values of the a and b variables from codeA.jl are modified, and I call the function f within codeB.jl.

The main question here is why? You’ve probably heard of the fact that using global variables is bad practice, in programming in general.

If you took a step back and gave us more context, you could perhaps get suggestions for a better design for your code.

[clopen]

Thanks for your answers and questions.

The main reason to organizing the code in this way is because I am working passing the code from Fortran to Julia, so I would like (at least at first) to kept the code similar to the Fortran version, which is organized in this way.

Finally, following the advice of @uniment, I declared the type of all ‘parameter’ variables in its definition, and I create a module called globals in CodeA.jl with all these parameter variables (as a and b variables of my example). Then I could remove a and b from the input and just call f(x), and f recognize the value of globals.a and globals.b. The number of allocations is almost 0 and the runtime has improved almost x10.

I get almost the same time following the recommendation of mapping f(x) → f(x,a,b) or similars.

Now I am considering to follow others recommendations, as to use a (mutable) structure or the parameter package to define the variables instead of the module.