From: Tim K. <tim...@ce...> - 2010-06-24 13:59:39
|
Dear all, In my application, I am doing quite a lot of assignments between NumericVector instances (which always are PetscVector instances in my case), and since I am lazily most of the time working in opt mode, I didn't realize that I was sometimes assigning vectors of different type, that is in particular v = w where v is a parallel vector and w is a ghosted vector. I now ran my application in devel mode for some other reasons and noted that this violates an assert. Well, on the other hand, there seems to me no principal problem in doing such an assignment since v requires less information than w provides. Therefore, I would suggest to explicitely allow assignments of this type here. In particular, in opt mode, this seemed to work quite well (effectively, it just does VecCopy() in this case). If you guys agree, I would change PetscVector::operator= accordingly after asking in the petsc-users list whether this is an allowed operation or, if not, what else operation would be appropriate at this point (unless Jed should just comment on this right away, in which case I wouldn't have to ask...). Best Regards, Tim -- Dr. Tim Kroeger tim...@me... Phone +49-421-218-7710 tim...@ce... Fax +49-421-218-4236 www.mevis.fraunhofer.de/~tim Fraunhofer MEVIS, Institute for Medical Image Computing Universitaetsallee 29, 28359 Bremen, Germany |
From: Jed B. <je...@59...> - 2010-06-24 14:24:10
|
On Thu, 24 Jun 2010 15:59:22 +0200 (CEST), Tim Kroeger <tim...@ce...> wrote: > Dear all, > > In my application, I am doing quite a lot of assignments between > NumericVector instances (which always are PetscVector instances in my > case), and since I am lazily most of the time working in opt mode, I > didn't realize that I was sometimes assigning vectors of different > type, that is in particular > > v = w > > where v is a parallel vector and w is a ghosted vector. This is value-assignment, right? There is no behavioral difference between a "parallel vector" and a "ghosted vector" except that you can only VecGhostGetLocalForm and VecGhostUpdate the latter. You can certainly VecCopy between them. Whether your operator= should support these mismatched types is entirely up to you. I could imagine imposing an artificial restriction just to maintain consistent types so that you don't accidentally end up calling VecGhostGetLocalForm on an plain vector (with work to figure out why it's a plain vector and not VecGhost). It is fine to use VecGhost for all vectors, the cost of holding the (unused) local part in all the Krylov vectors is pretty small unless you have very small subdomains (e.g. 8^3), very bad partitions, or large overlap. The smallest possible storage requirements are with VecGhost used only for vectors holding the state (and residual); Krylov vectors, vectors held by the time integrator, checkpoints for adjoints, etc., can be held in plain vectors and copied into a VecGhost if you later decide you need the local form. The general operator= seems to make sense in this context. Jed |
From: Roy S. <roy...@ic...> - 2010-06-24 17:17:03
|
On Thu, 24 Jun 2010, Jed Brown wrote: > On Thu, 24 Jun 2010 15:59:22 +0200 (CEST), Tim Kroeger <tim...@ce...> wrote: > >> In my application, I am doing quite a lot of assignments between >> NumericVector instances (which always are PetscVector instances in my >> case), and since I am lazily most of the time working in opt mode, I >> didn't realize that I was sometimes assigning vectors of different >> type, that is in particular >> >> v = w >> >> where v is a parallel vector and w is a ghosted vector. > > This is value-assignment, right? Right. > There is no behavioral difference between a "parallel vector" and a > "ghosted vector" except that you can only VecGhostGetLocalForm and > VecGhostUpdate the latter. You can certainly VecCopy between them. Thanks for the update. > Whether your operator= should support these mismatched types is entirely > up to you. I suspect there are only two reasons why we didn't support mismatched types: when that assertion was written we only had serial and parallel vectors, a more extreme mismatch we didn't know for sure that you could VecCopy between mismatched types. These reasons no longer apply. > I could imagine imposing an artificial restriction just to > maintain consistent types so that you don't accidentally end up calling > VecGhostGetLocalForm on an plain vector (with work to figure out why > it's a plain vector and not VecGhost). As long as we don't copy the type over with operator=, this shouldn't be a concern for PetscVector. I'd say that we should at least enable "parallelvec = ghostedvec", if everyone agrees. "parallelvec = serialvec" and "ghostedvec = serialvec" would also be natural to enable. I'd just as soon eventually support every operator= combo, but the "morelocaldatavec = lesslocaldatavec" type operations would require us to add some code in addition to the plain VecCopy before they'd have the semantics I'd want. --- Roy |
From: Jed B. <je...@59...> - 2010-06-24 17:31:16
|
On Thu, 24 Jun 2010 12:16:53 -0500 (CDT), Roy Stogner <roy...@ic...> wrote: > I suspect there are only two reasons why we didn't support mismatched > types: > > when that assertion was written we only had serial and parallel > vectors, a more extreme mismatch Yes, these are different in a deep way, they have different communicators and likely different local sizes. > As long as we don't copy the type over with operator=, this shouldn't > be a concern for PetscVector. No, just recognize that LocalFunction(xghost,fghost); yplain = xghost; LocalFunction(yplain,fghost); would fail in LocalFunction instead of operator=. I think this fine, but it is a slightly later error message. > I'd say that we should at least enable > "parallelvec = ghostedvec", if everyone agrees. And ghostedvec=parallelvec. > "parallelvec = serialvec" and "ghostedvec = serialvec" Would the semantics of the latter be VecGetLocalForm(g,&lf); VecCopy(s,lf); VecRestoreLocalForm(g,&lf); or something else? Would it include a VecGhostUpdateBegin/End? If so, would that use ADD_VALUES or INSERT_VALUES? Jed |
From: Tim K. <tim...@ce...> - 2010-06-25 07:13:02
Attachments:
patch
|
On Thu, 24 Jun 2010, Roy Stogner wrote: > I suspect there are only two reasons why we didn't support mismatched > types: > > when that assertion was written we only had serial and parallel > vectors, a more extreme mismatch > > we didn't know for sure that you could VecCopy between mismatched > types. > > These reasons no longer apply. I've meanwhile implemented "parallel=ghosted", see attached patch. It solves the problem for me. (BTW, it's amazing what mistakes one can do even in such a trivial change; I had to correct my patch two times before it did what I wanted...) I'll commit that within the next days if nobody objects. > I'd say that we should at least enable "parallelvec = ghostedvec", > if everyone agrees. "parallelvec = serialvec" and "ghostedvec = > serialvec" would also be natural to enable. Conceptually yes, but as Jed pointed out, it's more difficult to implement. Hence I'd like to leave that to somebody else, as long as I don't need that for my application. (-: > I'd just as soon > eventually support every operator= combo, but the "morelocaldatavec > = lesslocaldatavec" type operations would require us to add some > code in addition to the plain VecCopy before they'd have the > semantics I'd want. Well, aren't "ghosted=parallel" and "serial=parallel" just the same as NumericVector::localize()? If that is the case, only "parallel=serial" is missing, which should, altough it's not just VecCopy(), still be quite easy. Then, the remaining directions can easily be composed of existing ones. Best Regards, Tim -- Dr. Tim Kroeger tim...@me... Phone +49-421-218-7710 tim...@ce... Fax +49-421-218-4236 www.mevis.fraunhofer.de/~tim Fraunhofer MEVIS, Institute for Medical Image Computing Universitaetsallee 29, 28359 Bremen, Germany |
From: Roy S. <roy...@ic...> - 2010-06-25 16:09:58
|
On Fri, 25 Jun 2010, Tim Kroeger wrote: > I've meanwhile implemented "parallel=ghosted", see attached patch. It solves > the problem for me. (BTW, it's amazing what mistakes one can do even in such > a trivial change; I had to correct my patch two times before it did what I > wanted...) I'll commit that within the next days if nobody objects. Looks good to me. >> I'd say that we should at least enable "parallelvec = ghostedvec", if >> everyone agrees. "parallelvec = serialvec" and "ghostedvec = serialvec" >> would also be natural to enable. > > Conceptually yes, but as Jed pointed out, it's more difficult to implement. And he's quite right; the proper definition of "natural" here is "doesn't need any more code", for which "doesn't need any more communication" is necessary yet not sufficient. > Hence I'd like to leave that to somebody else, as long as I don't need that > for my application. (-: Certainly. > Then, the remaining directions can easily be composed of existing ones. Yeah, there's no functionality here that's really *missing*, just some that's not as convenient as it should be. But it'd be easy to augment the API incorrectly if we didn't have an application really testing it, so we'll wait until someone else needs it. --- Roy |
From: Tim K. <tim...@ce...> - 2010-06-28 09:07:00
|
On Fri, 25 Jun 2010, Roy Stogner wrote: > On Fri, 25 Jun 2010, Tim Kroeger wrote: > >> I've meanwhile implemented "parallel=ghosted", see attached patch. It >> solves the problem for me. (BTW, it's amazing what mistakes one can do >> even in such a trivial change; I had to correct my patch two times before >> it did what I wanted...) I'll commit that within the next days if nobody >> objects. > > Looks good to me. Okay, it's in now. Tim -- Dr. Tim Kroeger tim...@me... Phone +49-421-218-7710 tim...@ce... Fax +49-421-218-4236 www.mevis.fraunhofer.de/~tim Fraunhofer MEVIS, Institute for Medical Image Computing Universitaetsallee 29, 28359 Bremen, Germany |