Isn’t this what you meant before by overwriting randn(2,2)?
No, I mean to use a macro to rewrite it.
You can’t do that if the randn is inside a function
The type is usually easily inferable from parameters
Not in this particular case, and not in others.
but no we’ll not accept any change to base or pacakages that uses global state like this.
Meaning it would get rejected at package registration?
You apply that macro in that function.
It is not the case only in the cases where you’ve claimed it doesn’t matter. You either don’t care about the type/performance what so ever, in which case you don’t need anything, or you want to construct something with fixed type, in which case convertion is done automatically.
It means that if it somehow passed review because someone didn’t look closely, it’ll be strongly recommanded against.
No attempt will be made to fix the breakage caused by it in any other package or Base and PR’s trying to do that will be rejected. People seeing crashes or other issues when using it will not get any support apart from being told that the package is invalid and asked to stop using it.
I’m not sure if I’m qualified to join this discussion, and sorry if I miss the point, but I kind of share the same feeling that it would be easy if I could just write 0.0 instead of zero(T) and it be a Float32 when everything else in my function is also defined as Float32’s, or if I’m clearly calculating things up to that precision.
I also really don’t like the idea of a global state at all, but something that changes the default behaviour inside one function, depending on a type parameter/ parametric variant that is supplied to it, would be something that I think wouldn’t be against the contained/functional mindset of julia, and would allow for simpler to write/(understand) code.
Aside from all this, I feel like it would also clear up some easy to make mistakes/bugs, especially for beginners.
But please correct me if I’m wrong.
Right. That’s exactly the @parse_precision macro @ChrisRackauckas proposed above does. It can even be generaized a little bit to handle a few other constructs too. As long as you still use local information to supply the type, it’ll be perfectly fine.
True,
Maybe it would be possible to also change 0.0 to Float32(0.0) etc, and any other operation that results by default in floats by using a macro, but at what point does it then get harder rather than easier. I mean it’s true that all the tools to accomplish the desired behaviour are already available, it’s just something I wouldnt mind to see and use.
OK, let’s sum up. Sorry for the long post, I’ll try to be complete.
Julia has a default precision, which is Float64. This shows up in the parsing of x = 4.2, or in x = randn(2,2). This means that any code (script or function) that uses these features will work in Float64. This behavior is standard, but by no means fundamental: there’s no particular reason to impose Float64 as default rather than Float32 after all, other than this is what’s commonly done (e.g. in numpy, matlab).
It is sometimes useful to work in a different precision. For instance, I might want to test the numerical stability of my algorithm, and test it in Float32 to see if it breaks or not. Or I might run into a particularly weird case that requires extended precision. Or I might want to experiment with doing computations in Float32. In other languages this is thoroughly annoying to do, in julia this is trivial (and efficient) as long as 1) one remembers to write generic code for all functions 2) one seeds the algorithms with the correct element type.
I’m not saying that the situation is bad. I actually don’t know any other languages where it’s as easy. However, it is annoying (at least I find it annoying, and I don’t think I’m the only one) to write completely generic code when all you want 99% of the time is to operate on arrays of Float64. It is also annoying to go hunt for constants and take care to set their types. You can ease that with macros, which is very good for a single script or function, but is annoying as well when you have many functions, each calling randn(2,2) or zeros(2,2) or whatever. Don’t get me wrong, I like generic code as much as anybody, and I know that it is bad practice (non-generic code) to call randn() without giving a type explicitly. I’m just saying that it’s the path of least resistance, and so the path everybody who’s not a base/package developer is going to take unless they have very good reasons not to.
Given all that, it would be nice if there was a way to change the default precision. (not essential, but nice). Now we come down to earth and we ask “is it implementable?” Well, it appears so:
function my_fun(n)
sum(randn(n)) #just an example but actually do something more complicated
end
display(my_fun(2))
import Base.randn
Base.randn(dims::Integer...) = randn(Float32, dims...)
display(my_fun(2))
The generated code for my_fun is as fast as if there was no overloading. So one could imagine a package (let’s not even talk about base) that would overload the randn(), zeros() etc. functions (let’s forget about the constants for the moment). Obviously such a package would be of use only as a diagnostic/experiment tool, and not to be used by serious packages.
Con: it creates a global state. I agree global states make it harder to reason about code, and can lead to subtle bugs, but they’re still extremely convenient. Usages in julia that I can find (by typing set TAB in the REPL), that affect the behavior of code: set_zero_subnormals, setprecision, setrounding. I’ve not understood the status of https://github.com/JuliaLang/julia/pull/23205, but if it’s merged, then add that to the list. The number of threads is also a global state.
@yuyichao, you feel very strongly that this would be a bad idea, even as a registered external package (which, as far as I can tell, doesn’t imply any sort of acknowledgement from the julia devs about the quality of the package, and whose bugs are the responsability of the package developer, not julia devs). You take a similar position in the RNG PR. Is this because 1) you don’t believe in the use case 2) you don’t want another global state 3) there is something fundamentally wrong about redefining functions in Base, and it will make julia misbehave in some way? 1) and 2) are duly noted, but, with all due respect, I’m not sure why that should be grounds for rejection of a package. You seem to imply 3) in some of your comments (on this and on the RNG PR) but I’m still not sure why (and others don’t seem sure either). Could you clarify?
It’s type-piracy and so re-defining Base.randn on standard types could break other people’s packages just by a user doing using Package. That’s not a good safe behavior to promote.
Most of the time in any real situation the types are inferred by the arguments. I would argue that instead of my_fun getting n, you should pass in the array. In which case, that’s sum which will already use the eltype to define zero(T) and all of that. Really, the only special cases are array constructors where sometimes you just want to say “create an NxN array”, which is a statement that isn’t complete so we obviously have to fill in some defaults (1-based indexed, Float64, etc) with the ability to change the defaults. Of course you can come up with other edge cases, but in those you can obviously just have a default argument which is for passing internal types to use. That’s really not that bad…
But if you think about it, with differential equations the user has to give you an array to start, so you can use that type. In optimization the user has to give an initial array, so you can use that type. For linear solving the user gives you a matrix and an array, so you use those types. In machine learning the user gives you arrays of parameters and an initial condition, so you use that type. In plotting the user gives you an array of values, so you use that type. When writing special function approximations the user gives you an x to approximate at, so you can use that type. In almost every case it’s very clear, without any user intervention, where to get the type. An in the few other cases, you just add one argument when you want to switch from the default. And doing it this way is fully explicit, no global state, and compiles all of the abstractions away. Is this really such a big deal?
As I stated earlier, I only see it as a big deal for creating the arrays to give to the library functions, since for example with machine learning I want to make everything Float32. In that case, a macro which creates a local scope that parses everything to a different type is both explicit and syntactically simple, so I see that as a good solution (where you have to hardcode some array constructor edge cases).
I thoroughly believe that engineering a whole non-local system for changing defaults for types that everyone is supposed to plug into is far too much over-engineering that would only make things more complicated, and the elegance that would be gained for writing small scripts would be lost when trying to read a script since you will not always know the “precision setting” it’s supposed to be run at.
Really? If I redefine randn(2,2) to mean randn(Float32,2,2) instead of randn(Float64,2,2), what will it break?
About your other points: you speak as a package developer, in which case I 100% agree, not only is there no need to change the default precision, but it should be avoided at all cost and the precision should be taken as input from the user. I speak as a user, in which case you do need to set the precision somewhere (or it will be set to Float64 for you). In such a situation, I will typically write a few tens of functions, maybe put them into a module if I’m feeling fancy, then have a script with some setup code calling those functions. In this situation I don’t care about the “user” (which is me) calling with fancy types, because I know I do everything in my setup code in Float64. And then I figure it would be nice to try Float32 or BigFloat, and I have to rewrite my functions to be generic (which I don’t want to because this is code that I will throw away in two days when I figure out my method is crap anyway) This is the use case in which it would be very useful to set a default precision. I think I read a blog post not so long ago urging the community to target users and not package developers, now where could that have been 
Nope. That’s not a default precision, it’s a type. Precision may be a property of the type. It’s perfectly fine to change this default property of a type (and we have that for BigFloat) but it’s absolutely not to change the type.
Nope. There’s very easy convertions and also integer / math constants.
The type is one of the most fundamental thing in julia so it’s as fundamental as that.
You don’t have to write generic code. In fact, a lot of the code I write are not and I have no interest in making them generic. I also won’t ever run them on random input types. You just can’t write generic code while pretending you don’t want to write it.
And we have very good tools to catch those.
Nope, it’s the path that’s the hardest to maintain and is almost guarantee to crash from time to time.
Nope, this will break a lot of code and the overwite is not immediately visible in a function.
Well, afaict this official forum gets a lot of questions about package. If we bang question about the proposed package here as well as any official channel, I’ll have no problem with it.
Nop. If you actually read that issue, you’ll notice that I do support being able to change the global RNG. So 1) I do believe the use case there (but not this one) 2) I think that global state is fine (because RNG is very special, and it doesn’t apply to this one) 3) and there’s everything fundamentally wrong with redefining functions in base and it’ll make everything misbehave especially when you are chaning the return type!!! And just to make it clear, I’m only objecting to the implementatin in the RNG PR, not it’s idea and the main difference is that pseudo RNG fundamentally requires global state and swapping RNG does not and must not change the return type.
I’ve made it very clear in the thread already. You are destroying package precompiliation, system image, have an API with very strange behavior (not immediately visible) and in your case, you’ll break everyone’s code.
Everyone that expect what the function to do what they are documented to do.
It sounded like a fun way to practice macro foo, so here is a gist that transforms
@parsetypes (Int8,Float16) 1.0+1
@parsefunction rand Float32 rand(3,3)
to
Float16(1.0)+Int8(1)
rand(T,3,3)
I’m not going to take it any farther, since I don’t actually need it.
@yuyichao I’m not going to engage further, this is getting nowhere. But there is a default precision for floats in julia (or if you prefer, a default floating point type), and there is no fundamental reason why it should be double precision. You might not believe anyone would want to change that, but I do. It is still (after 30 replies on this topic) unclear to me why it would be impossible to change the default behavior. There must be something deeply wrong with me both for wanting an obviously stupid feature and not understanding why it obviously can’t be done, and let’s leave the discussion at that.
I have been reading up on Cassette.jl, and it looks like the perfect way to test what can and cannot be done. I’ll pick this up when there is a usable version and try and make a prototype implementation. Then we can discuss exactly how stupid this is on a concrete code. What fun we will have.
You make it sound like a historical accident, but it isn’t. For scientific computation, the trade-off is usually between economizing on CPU and memory use, and worrying about conditioning & catastrophic cancellation. Float64 turns out to be a sweet spot for a lot of numerical work. That’s why many (most?) languages for scientific computing use it as the default.
I think that everyone is sympathetic to the fact that sometimes you want to use other types with more/less precision (or range, for integers). There is no argument about this. It’s just that yet another global state turns out to be a bad solution.
You may disagree about this, but I frequently find that making more-or-less “global” arguments (ie something which every single function in my code use as is 99% of the time) explicit leads to improved code. The additional bonus is that it meshes well with the compilation model of Julia. In an MCMC library I am working on, I made the RNG an argument everywhere, and found it really convenient.
You make it sound like a historical accident, but it isn’t. For scientific computation, the trade-off is usually between economizing on CPU and memory use, and worrying about conditioning & catastrophic cancellation. Float64 turns out to be a sweet spot for a lot of numerical work. That’s why many (most?) languages for scientific computing use it as the default.
That strongly depends on what you’re doing, but in many applications, if your algorithm is ill-conditioned it will blow up in Float64 anyway, and if it is well-conditioned, you lose maybe 4 digits to error accumulation, which still brings you to 4 usable digits. I’m not an expert in this, but it seems to me that 32 bits was the “default” at one point: C “float” mean 32 bits, for instance. I imagine that’s because double were not widely implemented in hardware. Then hardware got fast at doing doubles, and double precision became the default for scientific computing. It seems to me that the primary reason why it is so is because it’s so inconvenient to change precision in many languages, and if you have to pick one, you pick the safest that’s reasonably fast on the architecture you target. I don’t believe it is the sweet spot for all applications: it just happens to be the most sane default to pick if you have to pick one. As a result, everybody is locked into doing everything in double precision. Julia is a good opportunity to change that and allow users to experiment with what they need for their particular use case. It is already in the current form, it would just be simpler with a global set_precision.
Again, I’m not arguing with the fact that passing arguments around is a much cleaner solution. I’m just saying that most users (not package developers) are not in the business of writing clean code, and will take the path of least resistance (which is this case means writing non-generic code, which makes it harder to change precision later on). It is not a moral failing of users. Having a way to set the default precision, like basically anything that changes the global state (default RNG, treatment of subnormals, bigfloat precision…), is of use in user code. It does not absolve anyone from writing nice and generic library code, which should not depend on the default precision anyway.
Yes and it’s not what I prefer, it is what it really is and it makes a fundamental difference as for what you can do and what you can’t.
I don’t have particular opinion on that.
As I said, the macro that does local change to your code is perfectly fine. However, because it’s a type and not a precision and so you must not change the return types of other functions.
I wouldn’t call it “deeply wrong with you” but I believe the piece you are missing is that you don’t think types are fundamental. It is. Everything in julia relies one it.
Again, why is this exactly? Isn’t this what defaultRNG() does in https://github.com/JuliaLang/julia/pull/23205 ? Or is that PR invalid? (I’m still unclear on what the conclusion of that discussion is) As far as I understand, every function calling into a function that is redefined will get recompiled, which is fine by me. Of course it will confuse the expectations of Base and library code which are that zeros(2) will produce a Float64 array, but as the untyped array-producing functions are more user- than library-facing, I’m not sure that will be a big problem - I’ll see when I get around to implementing it.
I might be very wrong about the behavior of redefining functions, I just haven’t seen any evidence to the contrary; if I am please explain. Please understand how frustrating it is to be told something does not work and I should not do it under penalty of mysterious horrors, when all appearances, simple code, and an open PR that seems functional shows that it works fine.
An illustrative example of why it could be harmful to redefine functions:
julia> Base.similar(A::Matrix, ::Type{T}, dims::Dims{2}) where T = similar(A, Int)
julia> rand(2,2)*rand(2,2)
ERROR: InexactError()
Stacktrace:
[1] convert(::Type{Int64}, ::Float64) at ./float.jl:679
[2] matmul2x2!(::Array{Int64,2}, ::Char, ::Char, ::Array{Float64,2}, ::Array{Float64,2}) at ./linalg/matmul.jl:679
[3] generic_matmatmul!(::Array{Int64,2}, ::Char, ::Char, ::Array{Float64,2}, ::Array{Float64,2}) at ./linalg/matmul.jl:478
[4] *(::Array{Float64,2}, ::Array{Float64,2}) at ./linalg/matmul.jl:146
This I get : it is not redefining functions in itself which is problematic, it is that the redefined function will be, well, redefined, and cause bad interactions with other functions calling that function (what I call “confuse the expectations” above). My hope is that library code does not call array-building functions without specifying a type (because that is non-generic programming), and therefore this kind of things will not happen. The point I’m confused about is whether it is legitimate at all to redefine functions (with the understanding that, of course, it will change behavior)
Sorry for hopping in, but what is wrong with using the functions’ overloads to have ::Type{T} as their first argument?
Above, you mentioned about C and some historical reasons of why float corresponds to single precision. Yet, when you write in C/C++ 1.0 it is double. float x = 1.0 uses implicit conversion from double to float there. A similary behaviour can be obtained, if you are not willing to change the way you write your code when calling randn, as simply as:
function myfunc(::Type{T}, n, m) where T<:AbstractFloat
A::Matrix{T} = randn(n,m) # or, A::AbstractMatrix{T} if you are going to use sprandn(n,m,0.1)
return A
end
myfunc(n, m) = myfunc(Float64, n, m)
function otherfunc(::Type{T}, n, m, o) where T<:AbstractFloat
result::T = n+m+o
return result
end
otherfunc(n, m, o) = otherfunc(Float32, n, m, o)
otherfunc(Float32, 5, 6.0f0, BigFloat(10))
where you declare your local variable to be of some type, which in turn calls the corresponding convert(T, value) function implicitly. This is more or less what happens when you write float x = 1.0. Obviously, this has also been pointed out for types by @yuyichao: