I try to choose a parameter vector which limits the number of complex constraints when it is easy to do so, but in general I do not see a practical way for transforming the parameters of a problem with a generic non linear constraint function G from R^m to R^p and H from R^m to R^q (m=dimension of parameter vector, p = number of equality scalar constraints, q = number of inequality scalar constraints) such that the constraints G(x) = 0 and H(x) >= 0 can be replaced by no constraints or box constraints. This would require “inverting” G, H.
Not “in general” if you want to automate it and/or the boundaries themselves are implicit, but for statistical problems, you can find a way most of the time with a bit of effort. See eg the “Transformations of Constrained Variables” in the Stan manual. I have a package,
which implements some of these, but it will soon be replaced another which has a much simpler interface and relies on some Julia v0.7 features.
2 Likes
I am working on some projects related to quantum information and machine learning. We recently implemented a framework for quantum algorithm design (though, it will still takes months to make it mature). Yao
To talk about this, I’ll talk about the problem we met in our practical work, which could be hard for other languages (compilers). And why we choose Julia.
Our first language used in research was C++ (back in 2015), however, C++ is not quite friendly for new students and most physicisist since nobody want to care about what is not their expertise. Besides, there are some old libraries in the group written in Fortran 95.
Before we know Julia, some people use MATLAB or python with numpy. MATLAB is fast indeed, but it is harder to write and has a quite long start time (though nobody will close it until everything is finished). Besides, it is commercial which means it won’t be friendly to students (like PhDs). We found that our university only bought the basic feature for us, which is totally not enough when we want to use GPUs or clusters.
On the Python side, the dynamic feature increases our speed of implementation, however, it will be much slower if you do not deal with your implementation details carefully. Yes, we do have pypy and numba, however, I can still not use high level abstractions in python, since this will reduce the performance. Usually, if you serious about performance, people will choose to use python only as a And before pybind11 was implemented, for large project like PyTorch, people have to manipulate the C API themselves.
Then it becomes extremely hard when we try to extend PyTorch with complex numbers. Because you cannot tweak python part, I tried many times, and it seems to hard to maintain such extensions. #755
I didn’t actually expect to much on packages for such a new language. Therefore, it is the language itself attract me to use.
-
we like strong typed language, this helps development a lot. Especially when you do not actually know how others implement a feature.
-
the performance, we do not expect have something really as fast as Fortran/C, but we need a reasonable balance between development efficiency and the performance. Most scientific computing requires such feature.
-
Builtin language sugar for multi-dimensional arrays (tensors). In quantum physics, tensors are the basic mathematical structure that we will study, the lack of language sugar in numpy make Julia a little bit more attractive.
-
Multiple dispatch based OO, though this is completely different from other languages. But I do think this is very useful and intuitive for physicists and mathematicians (maybe also others). I cannot talk too much about this here. I collected some blogs and posts and disclosure comments in this repo Tutorial For Physicists.
-
meta programming, Julia’s lisp like meta programming allows us to do things that is hard to accomplish in other languages and can be very useful. In C++, there is a trick called expression template for lazy evaluation. However, libraries using this trick like mshadow, and Eigen will reports a long and un-readable error message if an error occurs. Comparing to Julia, we can make a lot of things happen in compile time freely with macros. I believe the TensorOperations package is an excellent example. I also used this feature to construct Hamiltonians myself. It save the performance without loosing your abstractions. I compared tensor operations in python and Julia under the issue of TensorOperations, it is amazing fast comparing to numpy’s einsum and tensor dot. Besides, we will be able to do code generation without writing a customized compiler like ATen. And I believe in the future, things like Tensor Comprehension will also appears in Julia (with less code of course.)
-
native GPU support, the reason that we wanted to extend pytorch is all about its GPU. And we found this would be much easier with CuArrays/CUDAnative/etc. and all of them has complex number support!
-
friendly interface to C. This is not Julia only, but it is a plus indeed. I think I cannot expect to live without writing C for performance enhancement, etc. at the moment, though I believe after years Julia will be good enough to do most things without C.
In the end, just to clarify a little bit, I’m talking about other languages is bad. I cannot write too much about the disadvantage of Julia comparing to python, etc. And it is not related to the topic. But I believe Julia did solve some key problem in our projects!
9 Likes
That’s how I felt originally, but shortly after starting with Julia (3 years ago), I learned that I could fairly easily get identical or better performance in Julia than C (and I was known for performance optimizations in C/C++). I haven’t had to program in C since I’ve learned Julia (and that, even though I’ve been writing very low-level code).
4 Likes
Me too. We achieve much better performance in pure Julia comparing to some other C++ backends. (in the Yao framework) But it is just safe to say we still need C 
1 Like
On my bucket list is learning Rust, because if I do have to use something other than Julia, I’d like to use something a bit more modern 
(I’ve been using C full time since 1981), and I’d never use it alone, just as an adjunct to Julia (you can ccall into it, just like you can also ccall into Fortran)
5 Likes
Some State of the Art Packages in Julia v1.0 - Stochastic Lifestyle and posted to Hacker News
11 Likes
What is the best (*) pure Julia package alternative to Stan?
Klara, DynamicHMC, Mamba, Lora, Turing, other…?
(*) I mean fast, able to converge with complex models, dealing easily with discrete parameters…
1 Like
Turing. DynamicHMC is a little bit faster, but Turing.jl is a big interface to handle all algorithms and wrapped DynamicHMC . DynamicHMC is (~5 times) faster than Turing's HMC implementation · Issue #26 · TuringLang/AdvancedHMC.jl · GitHub . Turing has a great Stan-like API and so it can be thought of as an analogue, with the added bonus of having many different Julia-based tools for doing the actual sampling.
2 Likes
What about Klara, Mamba or Lora?
Is there anything like rstanarm or brms for Julia? (automatically convert regression models specified by a simple equation to MCMC).
Do Turing or DynamicHMC let you work easily with discrete parameters?
I’ve read Stan now can use the GPU.
Is any of the Julia’s bayesian packages able to do it too?
Turing does discrete parameters. For continuing this question though, I think you might want this to be its own thread.
I will add MIDI.jl and MusicManipulations.jl to this list, as they are high “class” or high “level” libraries that allows one to programmatically manipulate and analyze music. They are much better than their python counterparts (only MIDI) has a counter part in python), because they are made for musicians, not electrical engineers 
3 Likes
Do you know if with these it’s possible/simple to play the midi sounds? (on Linux)
Naming isolated packages does the ecosystem disregard because things julia-magically work together and therefore packages do not only add to the ecosystem but multiply.
Prime example ForwardDiff which makes everything else nicer. You model a trebuchet with RigidBodyDynamics on one day and find the optimal proportions on the next day with ForwardDiff. I still owe you to do that actually, so in the meantime I just have the example of a simple Kalman filter which suddenly can be fitted to data because we have automatic differentiation and did I say that it is the fasted filter in the world too because of StaticArrays?
9 Likes
It ought to be possible to make a package that calls an external soundfont player to play MIDI files.
2 Likes
I would not be so drastic… Pyomo is indeed (in functionality) very similar to JuMP, and it has too lots of extensions…
1 Like
The problem with things like pyomo is that they don’t actually work with anything else, and quickly become useless in cases where e.g. you have problems that are too large for their tools to support. Last I checked, pyomo still did not have decent linear algebra support (think about how crazy it is not to be able to do linear algebra for these types of problems). Let’s suppose I’m wrong for a minute and it does have that. To do it efficiently, it’s probably going to have its own custom array types. They would have had to do whatever work is required to get them to work efficiently with numpy. This cannot be done in Python, because loops in Python are too slow, so perhaps they used Cython or are just wrapping some C++ library. Maybe they have something that only works in PyPy. Ok, well now I want to do it with sparse matrices. Will this work? Nope. Maybe I can plug their variables into a numpy sparse matrix, but remember, numpy is written in C, so if you plug some python variables into one of their matrices and try to do stuff with it it’ll be too slow. Now suppose I have my own custom type of sparse array (something which is already in itself a feat in Python but trivial in Julia), and I want to do (multi-threaded) symbolic sparse multiplication with that… well, that experience is going to be pretty damn painful in Python. In Julia, it’s a couple dozen lines of code, and you have to jump through exactly 0 hoops to get it to work with JuMP.
The problem with Python packages is always the same: nothing fast enough to be useful in any but the simplest use cases can actually be written in Python. You have to write that in C, C++, or some other performant language. In practice this usually means relying on some library (e.g. numpy, scipy, pandas) which has done this for you. What if I want to add functionality to that, or god forbid, just get two of those things to work together in a way that was not already envisioned by their designers? In those cases your options are really bad, and those cases seem to happen to me too often for me to have any tolerance for Python.
12 Likes
But what Pyomo does isn’t just creating an ampl .nl model and then feeding the solver engine with it? (at least for NLP problems)
If so the model creation time is often negligible…
You would think this, but that fails to be true for large problems (e.g. O(10^6) rows and columns). So you can either have fun finding a more direct way of getting the thing into a standard form, or you have to deal with severe performance limitations of using something like pyomo to do it. Either way, you run into my anecdotes above.
Ah, ok then… I don’t have experience with such huge models, my thoughts were that for such a huge models then also the resolution phase, those managed by the solving engine, would have been huge² at least, and so the “model creation” phase, in relative terms, would have been even more negligible (as huge/huge² = 1/huge).