These are chat archives for symengine/symengine

10th
Nov 2015
Isuru Fernando
@isuruf
Nov 10 2015 04:18
@harakas, I like the design. Can you add another method to Dict called add_dict which would take in a Python dictionary and add it to Dict? Also rename Dict to DictBasic. With these two changes, we can merge your code into master.
Ralf Stephan
@rwst
Nov 10 2015 07:37
@certik @Sumith1896 @shivamvats Is #597 PR is a superset of the #578 PR? Is there active work on it?
Ralf Stephan
@rwst
Nov 10 2015 08:31
I made a branch which merges and compiles with master. Just give a hint and I make a PR: https://github.com/rwst/symengine/tree/polys_hash_set
harakas
@harakas
Nov 10 2015 08:47
@certik Wouldn't unordered_map be even faster for substitutions?
Isuru Fernando
@isuruf
Nov 10 2015 09:00
@harakas, unordered_map should be faster as we are doing only searches in the dictionary and we don't need an order.
harakas
@harakas
Nov 10 2015 10:38

Is this intended:

x, y, z, s0, s1 = symengine.symbols('x y z s0 s1')
expr = x * (y + z)
subs_dict = {y + z: s0, x * s0 : s1};
print expr.subs(subs_dict)
print expr.subs(subs_dict).subs(subs_dict)

results in

s0*x
s1

I'd like to get s1 straight away.

Isuru Fernando
@isuruf
Nov 10 2015 10:41
That's intentional.
(x+y).subs({x:y, y:z}) would give y+z
harakas
@harakas
Nov 10 2015 10:42
Um what's the connection?
Ah. I see.
And now I don't. I'm looking at Mul:subs(), and doing one final self-replace at the end would not break this behaviour, as the initial replace is done at the start already with a return.
Isuru Fernando
@isuruf
Nov 10 2015 10:53
It's easy to break things.
For eg:
Assume we do the same Add::subs()
Earlier example would give 2*z
Another eg: I would expect (x - y + z).subs({x : y, z: x}) to be x not y.
harakas
@harakas
Nov 10 2015 11:00
I'm afraid you misunderstood my question. My question is why does the first subs() not result in s1. It's easy to achieve without breaking any of your examples.
I'm not suggesting doing double-substitutions or triple or whatever. The problem is I think within the ::subs() functions themselves -- after substitutuing their arguments, they don't test the final result for substitution.
Isuru Fernando
@isuruf
Nov 10 2015 11:13
I think doing that is double substitution. In your case y+z : s0 is substituted first and s0*x is never found and you'd want to check for s0*x after subs. In general case, it may not be the intended behaviour. sin(sin(sin(sin(a)))).subs({sin(a):a})would result in a, but it should be the programmer who should handle those cases, otherwise another person who wants just one sin replaced, can't do it.
harakas
@harakas
Nov 10 2015 11:15
I see now. That's a better example. You want subs() to be something that maps parts of equation into something fixed. In my application I desire a recursive substituter though. I guess I can just call it repeatedly until nothing changes.
Ondřej Čertík
@certik
Nov 10 2015 14:39
@rwst regarding #597 and #578, none of the students have time to work on it right now. I wouldn't make ourselves depend on it. We can take useful work from the PRs and open new ones.
I see, you did already. I used your branch to open a PR in #658. We'll see if it passes tests and so that we can review it.
Ondřej Čertík
@certik
Nov 10 2015 14:45
@harakas unordered_map would probably be faster, though presumably the substitution dictionary is small, and std::map is typically pretty fast for small dictionaries (e.g. can be faster than unordered_map).
Ralf Stephan
@rwst
Nov 10 2015 14:45
@certik In my branch "just" the two tests are failing with Piranha=yes, but not with Piranha=no.
Ondřej Čertík
@certik
Nov 10 2015 15:38
@rwst do you know why they are failing?
Ralf Stephan
@rwst
Nov 10 2015 15:44

@certik the code in printer.cpp:

        std::vector<long long> out(n_vars);

with n_vars=4 creates a std::vector of length 0, capacity 0 which is probably not intended.

no wrong size 0 is intended.
Ralf Stephan
@rwst
Nov 10 2015 15:51
ah, decode takes the argument const int_type &n but is given int nvars=0 which subsequently tests !=0
what nonsense
Ondřej Čertík
@certik
Nov 10 2015 15:56
I think @Sumith1896 didn't have time to finish the code, so it is left with bugs like this.
This code tries to implement Piranha style polynomials ourselves, using only a few Piranha's primitives. It's almost as fast, but I think @Sumith1896 said it is still a tiny bit slower (something like 13ms vs 15ms if I remember from one particular benchmark).
So I think a cleaner approach is to simply use Piranha directly, and just hook it up into SymEngine.
Another reason for just using Piranha directly is that the I like the following design: having the low level polynomial code not use reference counting or any other SymEngine machinery, and just implement polynomials as fast as it can. Then the higher level code just hooks it up to SymEngine, so that you can use it as part of expressions. If some benchmarks are slower than the low level code, then more stuff needs to be moved to the low level code. Piranha can act as this low level code.
The current approach in #658 on the other hand is implementing Polynomial as a SymEngine object, so using the design above, this must be the top level wrapper. But there is no low level code. Everything is done through this Polynomial. So we are probably losing performance and it's not clear how much. If we use Piranha, then we can simply rerun our benchmarks against pure Piranha benchmarks and ensure that we didn't lose any performance due to the SymEngine wrapper.
Sumith Kulal
@Sumith1896
Nov 10 2015 16:03
Yes, there is 13ms vs 15ms problem
Ralf Stephan
@rwst
Nov 10 2015 16:03
@certik so will Piranha be a requirement for SymEngine?
Ondřej Čertík
@certik
Nov 10 2015 16:03
It should be an optional requirement. It is already.
It can't be a hard requirement, since it doesn't compile using MSVC.
Could Sage use Piranha?
Ralf Stephan
@rwst
Nov 10 2015 16:05
@certik if it's optional and I implement expansions directly on it expansions will not be available without Piranha
Ondřej Čertík
@certik
Nov 10 2015 16:05
That's a drawback indeed.
@rwst do you have a solution?
Ralf Stephan
@rwst
Nov 10 2015 16:06
blame MS they are so rich but can't provide C++11
Ondřej Čertík
@certik
Nov 10 2015 16:09
The only drawbacks of using Piranha (as @bluescarni and I discussed many times :) are:
  • it doesn't compile using MSVC, but hopefully in a few years it will
  • it uses a copylefted license (GPL), so if we make it hard dependency, then symengine itself would have to be GPL, since it won't compile without a GPL library, and it is my understanding that you can only use GPL libraries in non GPL programs if the GPL libraries are not essential part of the program, i.e. they are only optional. If @bluescarni changes it to LGPL at least, then we could depend on it, just like we already depend on GMP anyway.
Sumith Kulal
@Sumith1896
Nov 10 2015 16:10
@rwst Really sorry for that. Let me know if you need to figure out something in that codebase, I'll try to help.
I'll have a look at that PR parallely too
Ondřej Čertík
@certik
Nov 10 2015 16:11
Piranha is also very slow and needs a lot of memory to compile, but that can be fixed by simply including it in a .cpp file only, so only that cpp file will be slow to compile, but the rest of symengine will be fast. So that's not an issue I think.
Ralf Stephan
@rwst
Nov 10 2015 16:12
@Sumith1896 if ka::decode complains when given an empty vector which size would you give the vector?
Francesco Biscani
@bluescarni
Nov 10 2015 16:13
The switch to LGPL will be done, I haven't gotten around yet to change the headers and the license files in the repo :) I was planning to use exactly the same licensing as GMP
dual LGPL+GPL
so hopefully it's zero overhead licensing-wise - at least for all those projects that are using GMP already
Ondřej Čertík
@certik
Nov 10 2015 16:14
@bluescarni thank you. That would indeed completely fix the licensing issue.
Ralf Stephan
@rwst
Nov 10 2015 16:15
@Sumith1896 see the gdb backtrace in symengine/symengine#658
Ondřej Čertík
@certik
Nov 10 2015 16:15
@rwst do you see a technical reason why Piranha couldn't be included in Sage?
Ralf Stephan
@rwst
Nov 10 2015 16:16
@certik no, yet another library which needs an interface maintained
Ondřej Čertík
@certik
Nov 10 2015 16:17
@rwst if it's used in SymEngine only, then you don't need to maintain the interface, only to ensure it builds (it's header only, so at least no linking issues)
And actually no build issues either, as those become symengine build issues, which we take care of.
Ralf Stephan
@rwst
Nov 10 2015 16:18
@certik I'm not thinking in Sage terms these days
Ondřej Čertík
@certik
Nov 10 2015 16:18
Why not?
Ralf Stephan
@rwst
Nov 10 2015 16:19
no life in that
well, you asked 8-S
Ondřej Čertík
@certik
Nov 10 2015 16:21
I personally don't use Sage much, but I know it has a large community, and so I want to provide as solution that can technically be used by both Sage and SymPy, and thus unify the communities.
Sumith Kulal
@Sumith1896
Nov 10 2015 16:21
ka is Piranha's kronecker array.
Hmm I didn't come across this
Sorry I'm on phone, the replies are messed up :)
@bluescarni can help with that?
Francesco Biscani
@bluescarni
Nov 10 2015 16:23
ka is short for kronecker_array, the class is documented here: http://bluescarni.github.io/piranha/doxygen/classpiranha_1_1kronecker__array.html Let me check a second
So the decode function can work with an empty vector, but then it requires an empty integer in input
If the input is not zero then something sounds wrong
Ralf Stephan
@rwst
Nov 10 2015 16:25
but it gets a 4
Francesco Biscani
@bluescarni
Nov 10 2015 16:26
right... so the size of the vector that is passed in input to that function, how is that established?
the size should be equal to the number of variables in the polynomial
Ralf Stephan
@rwst
Nov 10 2015 16:27
it is created via std::vector<long long> out(n_vars); with n_vars =0.
but n=it->first which is 4
Francesco Biscani
@bluescarni
Nov 10 2015 16:28
that sounds incorrect... where is the n_vars coming from? If we are operating on the ring of polynomials with zero variables, we are essentially dealing with scalar values aren't we?
harakas
@harakas
Nov 10 2015 16:29

More substitution stuff:

>>> x,y,z = symbols('x y z')
>>> (1+2*x*y).subs(2*x*y, z)
1 + 2*x*y
>>> (1+2*x*y).subs(x*y, z)
1 + 2*z
>>>

So I can't substitute keys that contain numbers. Is this a bug or a feature?

Isuru Fernando
@isuruf
Nov 10 2015 16:30
@harakas, it's a bug. #557
Francesco Biscani
@bluescarni
Nov 10 2015 16:30
just to give a bit of background: encode() takes a vector of integers and packs it into a single integer, decode does the opposite. decode infers the number of integers to unpack from the size of the vector you pass in as input
Ralf Stephan
@rwst
Nov 10 2015 16:31
@bluescarni I have no idea any further of what I wrote. Would calling decode(0,0) in case of a scalar be correct?
Francesco Biscani
@bluescarni
Nov 10 2015 16:32
the first parameter of decode should be a vector, into which it will write the result of unpacking the second parameter
harakas
@harakas
Nov 10 2015 16:32
What about (x*y*z).subs(x*y,z)? In sympy that goes to z**2, while in symengine nothing happens currently.
Francesco Biscani
@bluescarni
Nov 10 2015 16:33
what I am trying to say is that if the second parameter is a nonzero value, then the size of the vector cannot possibly be zero
Ralf Stephan
@rwst
Nov 10 2015 16:33
@bluescarni yes I meant that decode(vector of size 0, 0)
Francesco Biscani
@bluescarni
Nov 10 2015 16:33
ah yes that should work I think
it won't do anything, but it will not throw either :)
Ralf Stephan
@rwst
Nov 10 2015 16:34
@bluescarni well then it need not be called at all
Ondřej Čertík
@certik
Nov 10 2015 16:36
@rwst the solution regarding Piranha that I have in mind is that we should just use it, make it optional (which technically is pretty easy) and get the features that we want working. If people want to use them, they have to install Piranha. See how things go. If lots of people would like to have those features without Piranha, perhaps on MSVC, then we can provide our own implementation of the polynomials (perhaps simpler and a bit slower, but working without Piranha). I expect most users will install symengine via Conda, Hashdist, Sage or apt-get, and for all these platforms we can easily build it with Piranha. So it's not a big deal I think, and it can still be optional, since I see big value in being able to build and test symengine with MSVC (as it ensures good portability of the C++ code, together with gcc and clang).
Francesco Biscani
@bluescarni
Nov 10 2015 16:36
So basically you are trying to unpack a value of 4, right. 4 in binary is 00000100. This number is then chopped up in digits groups according to the size of the input vector. So if size is 2 then the two component of the vectors will be, after decoding [0000,0100]
but if the size is 4 then the components of the vector will be: [00,00,01,00]
this is roughly how kronecker codification would work on a 8-bit architecture
it's just fancy bit packing basically
you need the vector size in order to know how to interpret a sequence of digits
Ralf Stephan
@rwst
Nov 10 2015 16:38
@certik give me a hint please on how to do it optionally, you just mean an empty function ifndef piranha?
Ondřej Čertík
@certik
Nov 10 2015 16:40
Ralf Stephan
@rwst
Nov 10 2015 16:40
@bluescarni actually n_vars should be 1, so the bug is elsewhere not in that call
but thanks
Ondřej Čertík
@certik
Nov 10 2015 16:41
@rwst and you can raise an exception if you need to call Piranha from some function that has to be compiled even if Piranha is not available.
Francesco Biscani
@bluescarni
Nov 10 2015 16:41
@certik I agree on the value of compiling on multiple platforms - even though MSVC historically chokes on many standard constructs and forces you to modify standard code rather often. I hope they will either improve it or just switch to clang as their backend eventually
@rwst no problem
Ondřej Čertík
@certik
Nov 10 2015 16:41
@rwst And you can use the unique_ptr with a forward declaration of a class idiom to hide the actual Piranha implementation into a .cpp file
@bluescarni there is also an Intel compiler, and I found that all the changes that I needed to do for Intel were also needed for MSVC.
And for symengine, it hasn't been a big deal to support all the compilers, since we use C++11, but not very fancy constructs.
Francesco Biscani
@bluescarni
Nov 10 2015 16:43
@certik weirdly enough, the Intel compiler looked actually better than MSVC last time I tried :)
Ralf Stephan
@rwst
Nov 10 2015 16:43
enough for today, good evening gents
Ondřej Čertík
@certik
Nov 10 2015 16:43
@bluescarni I think it is a bit better. :)
Francesco Biscani
@bluescarni
Nov 10 2015 16:44
@rwst o/
Ondřej Čertík
@certik
Nov 10 2015 16:44
Good evening @rwst. It's morning here.
Isuru Fernando
@isuruf
Nov 10 2015 16:45
@harakas, We haven't implemented (x*y*z).subs(x*y,z) to be z**2. I thought sympy only does this in replace, but now I'm confused what the difference is in sympy's subs and replace
Ondřej Čertík
@certik
Nov 10 2015 17:06
@isuruf there are use cases for both behaviors ---- we want to have some kind of subs, that is as fast as possible, but might not be the most powerful. We also want to have a very powerful subs that does all that @harakas mentioned, but it will be slower. Users will then choose which one to use for the given situtation.