-
Notifications
You must be signed in to change notification settings - Fork 47
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Add upper_sqrt #524
Add upper_sqrt #524
Conversation
I didn't use an input type of !real, as the input could be any type. The output is complex though. I forgot the transpose the name to sqrt_upper, but corrected that afterwards. I don't know if my choice of Generic Types is the best choice. |
Thank you for transposing the name, that's great. My approach there is it's hierarchical, and conceptually it goes "with" the other kind of You said the inputs could be real, but you've only allowed complex inputs, and not even all of those (no Could you undo the Please could you also adjust the new test so it uses |
I'm sorry, I don't understand where I did that. In the tests I provided I test the (complex) sqrt_upper of real negative numbers. If I add !real to the input parameter, the test fails. So where do you suggest I use $AF? By the way, the input could be integer (from the users point of view), not only floating or complex. I guess it is somewhat promoted before csqrt is called. I admit that I don't know all that takes place. |
I guess you are right. I used the English custom of using the adjective before the noun. In Spanish we usually put the noun before the adjective, so raizcuadrada_superior (sqrt_upper) would be natural using both points of view. |
Ha. My editor removes trailing spaces automatically when I open any text file. I changed my Changes to your previous Changes. I wonder, why it is important? |
I'll do, but I have to understand first some more about how the types are assigned, so I can predict the correct type for the expected result when writing the tests, especially in the case of real inputs and complex outputs. |
Mainly because the Changes file is for the maintainer to update, which here is me. But in a very practical way, if I fix a different bug (and note it in Changes) while you're still updating this, suddenly there'd be a merge conflict because of competing changes to that file. That would create pointless extra work. |
If an ndarray of non-supported type is given, a supported one will be picked to convert it to, in this case the last one in the list, which is why, if you look, the You might gather that the changes I'm asking for would be vastly less time and effort for me to make myself; the idea of asking you to make them is so you'd be able to make better PP operations going forward :-) It's also pretty clear the developer docs could be even better. |
First of all, I hope that my previous comment explains how types work in practice a bit better. Please ask if not, since it's important that people other than me understand it :-) One technique that's perfectly reasonable (according to me) in making tests, is to write the test with a reasonable input, have somewhat of a guess at the "expected" value for A great technique for how to do that hardcoding can be seen in |
I still have to polish my GH manners, sorry :) |
I changed the GenericTypes to all the complex types. I understand that input may be of any type but it gets converted to one of the complex types before the sqrt is taken, which I think is correct. Otherwise, it would fail for negative real numbers. I'm not sure if I should change the default. Now it is cfloatMaybe I should change the name to csqrt_upper, to make it explicit that the result would be complex even if the input is real. I used now the Test::PDL is_pdl function. There is a strange result though: long is promoted to cfloat and float is promoted to cfloat, which I think is all right, but double and ldouble are also promoted to cfloat, thus loosing precision. I know I could change the default complex to cldouble, but that may be wasteful. Is there some way to promote automatically floats to cfloats, doubles to cdoubles and ldoubles to cldoubles? |
I understand and appreciate your effort :). I looked at r2C, but haven't understood yet how it manages to convert float to cfloat, double to cdouble and ldouble to cldouble. |
Ha ha, no need! That's going to be a per-project convention, and that's the one for PDL, is all. |
It will not, because you are calling the
The name-change is actually a great idea, please do that. You do need to follow my suggestion and allow all floating-point types with
Noted, which is great, thank you. You may not have noticed but you are defining
Please don't test against a calculation from the main result (
I already spelled out how type conversions work. You are ignoring my suggestion to allow all floats (including real ones), so when you give it one of those, you are giving it an unsupported type. It picks the last-listed one, as already previously spelled out. A further problem: |
We have been dephased: I had changed csqrt to sqrt. When I use csqrt and $AF, the imaginary part of the square roots of real numbers are silently removed, so that sqrt_upper(-1) becomes complex 0, a mistaken result. I woud expect i(). I have attempted an alternative solution, which is to use PMCode to explicitly convert real inputs to complex inputs where necessary using r2C. Would this be too wastful? This seem to produce the correct results and the correct types. Integer types become cdouble, float becomes cfloat, double becomes cdouble, the complex types yield result with their same precission. The only surprise was ldouble which becomes cdouble, not cldouble. It turns out that r2C: ldouble->cdouble. |
I changed the name, but haven't used $AF as explained in my previous response |
That had been done. |
Yes, I noticed but I don't understand. There are NAN's when taking sqrt of real negative numbers, but there shouldn't be when those are converted to complex. |
Cool! But it's not been pushed yet, so I can't see it. |
I think you need to use both
Generally speaking, PDL_IF_GENTYPE_REAL(complex,) $GENERIC() in = $i(), tmp = csqrt(in);
This was a bug in |
It seems that is what I needed, a test to check the input type (I had seen it before, and forgot about it). I removed my call to r2C, introduced something similar to your suggestion and rewrote the tests. It seems to work (here, at least). I also realized there is no csqrt function that could do, for example csqrt(-1) and return i. So I added csqrt_right, which is simply a call to csqrt after dealing with the types. It yields as usual a positive real part, so the result is in the right half of the complex name. For consistency, I renamed csqrt_upper to csqrt_up, which is also shorter. Maybe csqrt_right may be simply renamed to csqrt. Thanks a lot for your patience. |
By the way, 'make' after a small change used to take a long time in old versions of PDL. Now it is much faster as it only makes what is needed given the change, I guess. However, after renaming a subroutine in Ops.pd and running 'make', the renamed routine is not found. I guess the '.so' libraries are not remade. If I run 'perl Makefila.PL' and 'make' again, it takes a little longer, but now the new routines are correctly found. Is this the expected behavior? |
Another doubt: If I admit all input types ($A instead of $AF) to permit integer inputs, then long and longlong get converted to cfloat instead of cdouble. Which conversion is preferable? |
lib/PDL/Ops.pd
Outdated
tmp=csqrt(tmp); | ||
//if(creal(tmp)<0){ | ||
// tmp = -tmp; | ||
//} |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Please remove this commented code. https://en.cppreference.com/w/c/numeric/complex/csqrt indicates that csqrt
is already guaranteed by the standard to be in the right-hand sight, so in fact this operation adds no value, except that it guarantees returning a complex output which sqrt
wouldn't for negative real inputs. I think in that case, this should just be called csqrt
, what do you think?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Maybe I'm not understanding. Before adding the proposed csqrt_right I tried calling csqrt from pdl code and it failed:
pdl> p csqrt(-1)
Undefined subroutine &main::csqrt called
I think it would be useful to have a csqrt in pdl that would answer i, the DWIM result. That is why I proposed csqrt_right, in analogy to csqrt_up. On the other hand sqrt(-1) should complain or return Nan, as there is no real square root of -1. The careful programer could write instead sqrt(cdouble(-1)) and get the desired result, but for large ndarrays I guess it would have the same problem you commented above: it would require converting and storing the whole ndarray before acting with the complex function. Explicit conversion is as wasteful as implicit conversion, right?
The comments were to show that the idea behind this routine is the same as that behind csqrt_upper but with a different branch cut. They are commented though as the C csqrt does that internally. I'll remove then though.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Any conversion has a cost, but if it's what's intended, then it's fine. If it's not what's intended, that's a problem. My thinking for the name is simply that while as you say the expectation (thanks to the C-like nature of PDL) is that sqrt
of real negative number will give NaN
, this would give a new option. I think that it's in keeping with the C-like nature to just give it the same name as the C function, so csqrt
. I wouldn't call its behaviour "internal", since the standard specifies it must do that.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I changed the name from csqrt_right to the standard csqrt, for the case the user wants the usual complex square root of any input (real, interpreted as complex with zero imaginary part, or complex). If I have understood all the conversation, my $out=csqrt($real_in);
may be better than my $out=sqrt(cfloat($real_in));
or than my $cfloat_in=cfloat($real_in); my $out=sqrt($cfloat_in);
(or any other complex conversion) if there is no need for keeping the probably large converted input. So I believe csqrt
may be useful although there is a generic sqrt
.
lib/PDL/Ops.pd
Outdated
'infinitesimally' below the negative real axis. | ||
EOF | ||
Code => q{ | ||
PDL_IF_GENTYPE_REAL(complex,$GENERIC()) tmp=$i(); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Please change this to what I said:
PDL_IF_GENTYPE_REAL(complex,) $GENERIC() tmp=$i();
I didn't know that complex
by itself was a valid type, but because of how that macro operates (it's either the first arg, or the second arg), apparently it is but you'd be throwing away the precision if the input were long double. Did you try the thing I put in the comment, by literally copy-pasting it? If so, did it work?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It doesn't work:
$ make
...
lib/PDL/Ops-pp-csqrt_up.c: In function ‘pdl_csqrt_up_readdata’:
lib/PDL/Ops-pp-csqrt_up.c:97:59: error: expected ‘=’, ‘,’, ‘;’, ‘asm’ or ‘__attribute__’ before ‘tmp’
97 | PDL_IF_GENTYPE_REAL(complex,) PDL_Float tmp=(i_datap)[0];
| ^~~
lib/PDL/Ops-pp-csqrt_up.c:97:59: error: ‘tmp’ undeclared (first use in this function); did you mean ‘tms’?
97 | PDL_IF_GENTYPE_REAL(complex,) PDL_Float tmp=(i_datap)[0];
| ^~~
...
I hadn't tried until now copy&pasting your line because I didn't understood it. I thought the idea was to either put 'complex tmp' or '$GENERIC() tmp', but not 'complex $GENERIC() tmp', which I guess is what your macro call does. Somehow, the way I wrote it seems to produce the expected complex precission, according to the 'is_pdl' tests.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ironically here we're falling victim to PDL's hiding of C types. What we need here is complex float
where the real type is float
. So all we need to do to stick with that idiom and get what we need is (untested, but should work):
$TFDEGCH(PDL_CFloat,PDL_CDouble,PDL_CLDouble,$GENERIC(),$GENERIC(),$GENERIC()) tmp=$i();
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It didn't work with $GENERIC()'s, but it did with PDL_CFloat,...
I didn't know the $T.... macros.
Curiosity: How do you manage actual commas ',' within the code you want to expand? They seem to separate the alterenatives.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The only way to get commas inside a $T
or similar macro would be to enclose it in brackets, which normally would introduce other difficulties. Even that is a relatively new addition by me, and before that it was strictly comma-separated. It's not really for anything where commas would occur inside the brackets.
The How would you feel about making some notes (i.e. PR-ing) in e.g. PDL::DeveloperGuide or maybe PDL::Core::Dev (I recently added some notes at the top there about this new scheme of repo layout)? I have the disadvantage I don't truly know what's surprising or needs documenting, because I know it fairly well. |
This is due to the By the way, I asked you to rebase past the changes I made on |
It would be great, but I'll ask you to remind me in 1.5 weeks. (I have to do some intensive teaching) |
I really don't know which is better.
I did 'git rebase -i master', and then C-c C-c, but maybe I made a mistake. Then I couldn't push to origin/upper_sqrt. So I pulled from origin/upper_sqrt and failed. I 'set git config pull.rebase true' and then pulled and pushed succesfully. I'm not sure how to proceed to undo the damage and remake it correctly. My knowledge of git is very fragmented. |
Fair enough, let's leave it as is then.
There isn't any damage :-) You just need to do |
I tried to do an interactive rebase, reordered the commits and force pushed. Hope it is in the correct order now. Then I had to do it again, as I found another commit to master :) I hope not to forget what I'm learning :) |
Found another of your commits was out of order in my branch (that related to #517). |
Thank you for the rebase and force-push! The way that's normally done is you'd rebase against current A further opportunity when doing a
Make sense? If not, please ask :-) |
I don't understand the first suggestin. I combined my revisions as you suggested, but I can't push them: Pushing to github.com:wlmb/pdl.git So what I thought I should do is to incorporate your latter changes before mine and push again. But you say I can rebase against current master and delete any commits that aren't mine. Including the three previous and the more recent ones. Is that correct? I'm not sure what is the correct way of doing that: you mentioned 'git rebase -i master'. I guess I should update master first, but then, what should I do in the interactive rebase session? What should I do before pushing? |
The request is:
That would be the most important step, that removes the slightly absurd situation where this PR shows you trying to merge my code, which is already on The reduction in the number of commits would be a further, separate process. You would do the Edited to add: after reducing the number of commits, |
Thank you! I can tell that you left out the "synchronise your I cherry-picked the I've now force-pushed that to your branch, which is permitted when you select (or don't change) the "allow maintainer to push to your branch" option. As soon as that's happy, I'll merge this. Thank you very much for yet another contribution to PDL! |
On Thu, Jan 16, 2025 at 01:51:23PM -0800, mohawk2 wrote:
Thank you! I can tell that you left out the "synchronise your `master` from the upstream" step, because there were several commits different.
Hi,
I recall synchronizing it. The most recent commits in my master were when I rebased
c2cc6e4 master origin/master add CPPFLAGS building "pdl" executable to help Debian packaging
cf594d3 use OPTIMIZE building "pdl" executable
64b52cf add "update_data_from" method as streamlined get_dataref/upd_data
b0619ed r2C and i2C now handle long double both real and complex - #524
84d0abc capture current wrong upward propagate_badflag - #517
and currently, my master is behind by 4 commits only,
Add tests of r2C 
wlmb authored and mohawk2 committed 1 hour ago
89d5b51
Add csqrt_up and tests
wlmb authored and mohawk2 committed 1 hour ago
176dcba
Added csqrt and tests
wlmb authored and mohawk2 committed 1 hour ago
2bc153c
credit for new csqrt - PDLPorters#524
mohawk2 committed 1 hour ago
45b3795
and there is no other pending commit.
Sorry for all the mistakes I made on the way. There is much to learn about git.
I cherry-picked the `r2C` tests onto `master` immediately, ...
I've now force-pushed that to your branch, which is permitted when you select (or don't change) the "allow maintainer to push to your branch" option. As soon as that's happy, I'll merge this.
I didn't understand that last comment. Didn't you merge it already? Do
I have to do something else?
Thank you very much for yet another contribution to PDL!
Thank you for all your guidance and patience.
Best regards,
Luis
…--
o
W. Luis Mochán, | tel:(52)(777)329-1734 /<(*)
Instituto de Ciencias Físicas, UNAM | fax:(52)(777)317-5388 `>/ /\
Av. Universidad s/n CP 62210 | (*)/\/ \
Cuernavaca, Morelos, México | ***@***.*** /\_/\__/
GPG: 791EB9EB, C949 3F81 6D9B 1191 9A16 C2DF 5F0A C52B 791E B9EB
|
I don't know what to tell you - at the time of my writing the above, I pulled from your repository, and your first commit was after roughly 4 commits before the then-current
Git is definitely very powerful, but takes some getting used to. It took me months to get to a basic level, and many more months (a couple of years, really) before I got as comfortable as I am with it now. I'll bet the maths behind your photonics work didn't take five minutes to learn!
No, shortly after I wrote that I merged your PR. All is well!
Thank you for sticking with this. I know it's not easy. When you have capacity, it will be really valuable if you write on one of the open documentation issues (say #433) what you now know you didn't know, so we can at least hint at that in PDL::DeveloperGuide. |
On Thu, Jan 16, 2025 at 09:15:58PM -0800, mohawk2 wrote:
> > Thank you! I can tell that you left out the "synchronise your `master` from the upstream" step, because there were several commits different.
>
> Hi, I recall synchronizing it.
I don't know what to tell you - at the time of my writing the above, I pulled from your repository, and your first commit was after roughly 4 commits before the then-current `master`. It's not the end of the world.
I think I understand what happened. I tried to follow instructions
without understanding them. You said:
Thank you for the rebase and force-push! The way that's normally done is you'd rebase against current master, which means you can delete any commits that aren't yours. That would mean *you can delete those 3 from me*.
I found that strange, but that is what I did: I synchronized my fork,
pulled from my master, rebased my branch, and then *I explicitly
deleted* your last commits before force pushing.
Regards,
Luis
…--
o
W. Luis Mochán, | tel:(52)(777)329-1734 /<(*)
Instituto de Ciencias Físicas, UNAM | fax:(52)(777)317-5388 `>/ /\
Av. Universidad s/n CP 62210 | (*)/\/ \
Cuernavaca, Morelos, México | ***@***.*** /\_/\__/
GPG: 791EB9EB, C949 3F81 6D9B 1191 9A16 C2DF 5F0A C52B 791E B9EB
|
The "deleting my commits" part was fine! I suspect what didn't happen correctly was either the sync-ing your fork, or pulling from your |
Added upper_sqrt to take the complex square root of an arbitrary number and choose the branch that lies in the upper half complex plane.